Carl Love

Carl Love

28055 Reputation

25 Badges

12 years, 358 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

Define eta thus:

eta:= t-> piecewise(abs(H(t)/omega) <= 1, arccos(-H(t)/omega), 0);

If decimal numbers can be expressed as simple fractions, it's usually best to express them that way. This often leads to simpler answers.

Sol:= dsolve({diff(x(t), t$2) + x(t) = cos(4*t/5)/2, x(0)=0, D(x)(0)=0});

plot(eval(x(t), Sol), t= 0..20*Pi);

 

Don't use display. Just use F(1). Also, your code will be more robust if you replace animate with plots:-animate.

Here is your procedure, syntax corrected:

LU:= proc(A::Matrix)
uses LA= LinearAlgebra;
local
     i, j, k, n:= LA:-RowDimension(A),
     U:= Matrix(LA:-Dimensions(A)),
     L:= Matrix(LA:-Dimensions(A))
;
     if A[1,1]=0 then error "A[1,1] can't be 0" end if;
     U[1,1]:= A[1,1];
     U[1, 2..n]:= A[1, 2..n];
     L[2..n, 1]:= A[2..n, 1] /~ U[1,1];
     for i from 2 to n-1 do
          U[i,i]:= A[i,i] - add(L[i,k]*U[k,i], k= 1..i-1);
          if U[i,i]=0 then error "0 on diagonal" end if;
          for j from i+1 to n do
               U[i,j]:= A[i,j] - add(L[i,k]*U[k,i], k= 1..i-1);
               L[i,j]:= (A[j,i] - add(L[j,k]*U[k,i], k= 1..i-1))/U[i,i]
          end do
     end do;
     U[n,n]:= A[n,n] - add(L[n,k]*U[k,n], k= 1..n-1);
     if U[n,n]=0 then error "0 on diagonal" end if;
     L,U
end proc:   

Some errors that you had:

  1. You had an extra end do; that should've been obvious.
  2. You had n as both a parameter and a local.
  3. You tried to use square brackets [ ] for algebraic grouping; only round parentheses ( ) can be used for this.
  4. You tried to return multiple values in two separate statements: L; U; That needs to be L,U.
  5. Matrix should be with an uppercase M. The lowercase matrix refers to an older structure that shouldn't be used anymore.

I also introduced the error statement, whose use is better than just printing the word "error".

I didn't change the logic of your program at all. I claim neither that it produces nor that it doesn't produce a correct LU factorization. I haven't tested that at all.

Here's your program---debugged, cleaned up, and modernized:



nevilee:= proc(X::list, Y::list, x)  
local i, j, n:= nops(X), Q:= Matrix((n,n), Y, scan= columns);    
     for i from 2 to n do  
          for j to i-1 do    
               Q[i,j+1]:= ((x-X[i-j])*Q[i,j]-(x-X[i])*Q[i-1,j])/(X[i]-X[i-j])
          od  
     od;
     < <X> | Q >  
end proc:


nevilee([-2, -1, 0, 1, 2], [1/4, 1/2, 1, 2, 4], .5);

Matrix(5, 6, {(1, 1) = -2, (1, 2) = 1/4, (1, 3) = 0, (1, 4) = 0, (1, 5) = 0, (1, 6) = 0, (2, 1) = -1, (2, 2) = 1/2, (2, 3) = .8750000000, (2, 4) = 0, (2, 5) = 0, (2, 6) = 0, (3, 1) = 0, (3, 2) = 1, (3, 3) = 1.250000000, (3, 4) = 1.343750000, (3, 5) = 0, (3, 6) = 0, (4, 1) = 1, (4, 2) = 2, (4, 3) = 1.5, (4, 4) = 1.437500000, (4, 5) = 1.421875000, (4, 6) = 0, (5, 1) = 2, (5, 2) = 4, (5, 3) = 1.0, (5, 4) = 1.375000000, (5, 5) = 1.406250000, (5, 6) = 1.412109375})

(1)

 



Download NevilleInterp.mw

otat:= (a::`=`, b::list(name=algebraic), c)-> map2(eval, isolate(a,y), b):

otat(x-y = 1, [x=10, x=11]);

     [y = 9, y = 10]

I've translated the algorithm that you posted into Maple. However, the algorithm does not seem to be correct. Here's the translation:

`mod/SFF`:= proc(f::polynom, p::posint)
local
     i:= 1,
     R:= 1,
     x:= indets(f, name)[],
     c, w, g, y, z
;
     if x=NULL then return f mod p end if;
     g:= diff(f,x);
     c:= Gcd(f,g) mod p;
     w:= Quo(f,c,x) mod p;
     while w <> 1 do
          y:= Gcd(w,c) mod p;
          z:= Quo(w,y,x) mod p;
          R:= R*z^i;
          i:= i+1;
          w:= y;
          c:= Quo(c,y,x) mod p
     end do;
     R * (c mod p)
end proc:

The usage is like this:

f:= x^6 + randpoly(x):

SFF(f) mod 3;

To compare with Maple's prepackaged algorithm, use

Sqrfree(f) mod 3;

Your worksheet prompts the user for input. That is a very, very crude way (IMO) to use/write a Maple program. IMO, the only acceptable way for a Maple program to accept input is through a procedure's parameters.

Your wa and wb should be w*a and w*b.

When setting the values of the parameters ab, and w, you need to use := rather than =.

On the right side of eqx, you use "naked" x; this should be x(t).

On the right side of eqz, you use y; I'm guessing that this should be z(t).

The squaring operation can also be done on a list of lists using elementwise operators:

(`^`~)~(L,2);

For editing code in a worksheet, try Code Edit Regions. You can find them on the Insert menu. 

I see that you're using Maple 18. There are two nuisances with Code Edit Regions in Maple 18. I think that the first of these has been corrected in Maple 2015 and that the second hasn't. The first nuisance is that when you execute the code in the Region, if there's an error, it doesn't tell you where the error is nor does it place the cursor near the error. The second nuisance is that like most worksheet-embedded components, it's very hard to get your mouse cursor precisely on the control points that let you resize the Region, which you'll very likely need to do.

I don't know about Volterra equations, but your integro-differential equation can be easily converted to a third-order ODE-IVP which can be easily solved symbolically. So, I say, "Why bother with Volterra equations?" To do the solving entirely in Maple:


restart:

(intE, d0):= (diff(u(x), x) = 1+x-x^2+int((x-t)*u(t), t= 0..x), u(0)=3):

d1:= intE;

diff(u(x), x) = 1+x-x^2+int((x-t)*u(t), t = 0 .. x)

(1)

d2:= diff(d1, x);

diff(diff(u(x), x), x) = 1-2*x+int(u(t), t = 0 .. x)

(2)

d3:= diff(d2, x);

diff(diff(diff(u(x), x), x), x) = -2+u(x)

(3)

sys:= %dsolve({d3, seq((D@@k)(u)(0)=eval(rhs(d||k), x= 0), k= 0..2)});

%dsolve({diff(diff(diff(u(x), x), x), x) = -2+u(x), u(0) = 3, (D(u))(0) = 1, ((D@@2)(u))(0) = 1})

(4)

Sol:= value(sys);

u(x) = 2+exp(x)

(5)

#Verify solution:
eval(intE, u= unapply(rhs(Sol),x));

exp(x) = exp(x)

(6)


Download IntegroDiffEqn.mw

Here's another way to do it without unapply, using subs. This method is extremely powerful for three reasons:

  1. Whereas unapply can only be applied to a single expression, this subs method can be applied to any procedure whatsoever.
  2. This method can be used after the fact: i.e., after f, g, and h have already been defined.
  3. It is very fast, based on the built-in command subs.

However, my personal opinion is that imprudent use of this method can make code difficult to read. So, for this particular example, I'd use Kitonum's second method (my first choice) or unapply.

The subs method is shown in the line offset by (**) below.

f:= 2*x:  g:= x^2:
h:= x-> piecewise(x <= 5, _f, _g):
(**) h1:= subs([_f= f, _g= g], eval(h)); (**)
h1~([1,2,8]);

I show two other unrelated things in the above code fragment, one related to map and the other to piecewise. The former is that map can often be replaced by the elementwise operator ~. See ?elementwise. The latter is that inequalities compounded with and are usually unnecessary (and wasteful) as piecewise conditions. The conditions are processed left to right; the fact that the nth condition is being evaluated implies that conditions 1, ..., n-1 are false. I very often see even experienced Maple programmers using these unnecessary compound conditions in piecewise

 

Since both the commands D and diff already provide this functionality, I don't know why you'd want to create the procedure MultiDiff. But, assuming that you have some reason for doing it, ...

MultiDiff:= proc(f, n::posint, x::name) (D@@n)(f)(x) end proc;

The use command has several subleties; I find it very difficult to use (pun intended). Use fully qualified procedure names instead. That is, replace references to Color with ColorTools:-Color.

You've made a common Maple error that has nothing to do with either map or piecewise. I don't know if this error has an official name, so I'll call it an attempted indirect reference to a procedure parameter. (Perhaps that'll become its official name.)

Your h is a procedure (also called a function) whose parameter is x. Your expressions f and g depend on global x. It's not the same x. When you placed f and g within h's body, you were hoping that the indirect reference to x would be to the parameter x. It doesn't work that way. The command unapply is used to resolve these indirect references so that all references to x within an expression refer to the parameter x. Unfortunately, the help page ?unapply doesn't clearly state that.

First 242 243 244 245 246 247 248 Last Page 244 of 395