Carl Love

Carl Love

28095 Reputation

25 Badges

13 years, 100 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

The simplest solution is to simply not use gamma as a variable. The next simplest solution is to begin your code with the line

local gamma;

This latter solution will only work in relatively recent versions of Maple.

If either of these solutions is unacceptable, there are still other alternatives.

Zero-point and negative-point users have never appeared on the Users' list in the three years that I've been here. That's why I give zero-point users a vote up every chance I get (except for spammers, of course). Since only Questions, Posts, and Answers can be voted up---and Mcloone has never made any of those---I never had a chance to give him a vote up. Since only Questions, Posts, and Answers get indexed, you won't find him by searching for words from his Comments and Replies.

The following worksheet shows that there are two solutions for N which produce equally minimal residuals and which both lead to n = 14.


restart:

a:= (.0625*N+.0215)/2:

b:= a - .043:

c:= b/(N/2-1)+.0215:

SolN:= [solve(b/c=n, N)]:

nops(SolN);

2

(1)

N:= SolN[1]:

res:= abs(.0645*0.25*(833-(2*N-1)) + n*.043 - sum(sqrt(a^2 - (.043+i*c)^2), i= 1..n)):

plot([seq([k,eval(res, n= k)], k= 1..30)]);

 

eval(res, n= 14);

0.2320950e-1

(2)

eval(N, n= 14);

39.38892063

(3)

N:= SolN[2]:

plot([seq([k,eval(res, n= k)], k= 1..30)]);

 

eval(res, n= 14);

0.2320950e-1

(4)

eval(SolN, n= 14);

[39.38892063, 1.27507937]

(5)

 


Download SolveForN.mw

Yes, it's trivial. For example,

(a,b,c,d):= 'LinearAlgebra:-RandomVector(3)' $ 4:
m:= LinearAlgebra:-LinearSolve(<b|c|d>, a):
#Verify solution:
a = m[1]*b + m[2]*c + m[3]*d;

You just stumbled upon the real assumption, and you assumed that it had something to do with the simplification of the result to 0. But the real assumption is irrelevant. If you simply use simplify(a-c), you will also get 0. In Maple, the expression a - b (for previously stored algebraic expressions a and b) will return 0 if and only if a and b are stored as exactly the same expression. For example:

a:= (x+1)^2:
b:= x^2 + 2*x + 1:

a-b;

simplify(%);

     0

There are many ways in Maple to prove that two algebraic expressions are mathematically equal. This is an oft-debated topic here. Checking that the difference simplifies to 0 is one way. It has been proven that there is no general algorithm for it.

To avoid this error (which I call invalid indirect reference to a procedure parameter), I prefer to use procedures rather than expressions to as great an extent as is feasible (or reasonable). Specifically, in this case, I'd make W1, W2, and W all procedures.

W1:= (A,omega,k)-> A*cos(omega*t-k*x);
W2:= (A,omega,k)-> W1(A,omega,-k);
W:= W1+W2;

SW:= (A,omega,k)->
     plots:-animate(
          plot,
          [[W1,W2,W](A,omega,k), x= -4..4, y= -4..4,
           color= [red,green,blue], scaling= constrained],
          t= 0..5,
          frames= 10
):

SW(2,2*Pi,5);

Note that the sum of two procedures, W1+W2, can be used as a procedure and that a list of procedures, [W1,W2,W], can also be used as a procedure: In both cases, they can be applied to a sequence of arguments.

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);

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