Carl Love

## 27177 Reputation

11 years, 338 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

## coeff(b, x, 0)...

In this case, you need to circumvent the automatic simplification of x^0 to 1. You do this with

coeff(b, x, 0).

## plots:-surfdata...

We make Points1 into an (M+1) x (N+1) x 3 listlistlist:

A:= [seq](Points1[(M+1)*k+1..(M+1)*(k+1)], k= 0..N);

plots:-surfdata(A);

And you can add the same options that you used for the pointplot3d command.

## map Re...

If they can be considered the same thing and also be considered the same as 1, then it is trivial:

Re ~ (sol);

Then proceed with whatever else you were doing with the repetitions.

## Easier to remove repetitions than to fin...

If your goal is simply to remove the repetitions from a list L, and you don't care about the order, then it is trivial: convert to a set:

{L[]};

If you want to keep it a list, but you still don't care about the order, it is still trivial: convert to a set, then convert back to a list:

[{L[]}[]];

When order doesn't matter, removing repetitions, whether or not they exist, is much easier than detecting or finding them: The former can be done entirely in Maple's kernel.

## A list is not an rtable...

An rtable is a Matrix, Vector, or Array. A list is not an rtable. FindRepetitions returns a list (so converting it to a list doesn't do anything). Why don't you use the numelems command that acer recommended? If your Maple version is less than 15, use nops for a list.

## Another three simple ways...

If the implicit-form solution can be solved for the independent variable, as it can in this case, then there are several more very easy ways to get it.  Here's three: (Sorry, explicit uploading of worksheets seems broken right now, so I am uploading this in plaintext (prettyprint=1).)

`restart;ode:= x - y(x)/diff(y(x), x) = y(x);                            y(x)                               x - -------- = y(x)                           d                                       --- y(x)                                  dx            (1)dsolve(ode);                       /        /     x    \      \             y(x) = exp|LambertW|- --------| + _C1|                       \        \  exp(_C1)/      /isolate(subs(y(x)= y, %), x);                     x = -y (-_C1 + ln(y))(2) The name of the conversion below is y_x regardless of the names actually used in the ODE.convert(ode, y_x);                       d              x(y)                      --- x(y) = -1 + ----                       dy              y  dsolve(%);                    x(y) = (-ln(y) + _C1) y(3) The order of  the substitutions below is significant.subs(diff(y(x),x)= 1/diff(X(y),y), y(x)= y, x= x(y), X= x,  ode);                             / d      \                        x(y) - y |--- x(y)| = y                             \ dy     /    dsolve(%);                    x(y) = (-ln(y) + _C1) y`

## y*conjugate(y) = abs(y)^2...

For any complex y, y*conjugate(y) = |y|^2 (although I would've guessed that you already knew this identity). The identity can be applied automatically via

simplify(z);

The resulting equation can be solved for abs(y) via

solve(%, abs(y));

indicating that there is a circle of solutions in the complex plane.

## while loop and Norm...

How about this? With the data in your worksheet, this took exactly 200 iterations. The resulting plot is a straight line (as well as can be seen with the naked eye). Note the the default for Norm is the infinity norm, which is equivalent to checking the absolute values of all of the entries. In other words, the convergence criterion is that the distance for both coordinates of all n-1 points is less than epsilon.

ErrX:= Vector(n, [0,1]):
ErrY:= Vector(n, [0,1]):
rho:= 1.2:  #relaxation factor
epsilon:= 10.^(-5):
while LinearAlgebra:-Norm(< ErrX, ErrY >) >= epsilon do
for k from 2 to n do
Xk:= X[k] - rho*GX(X, Y, k)/GXX(X, Y, k);
ErrX[k]:= Xk - X[k];
X[k] := Xk;
Yk:= Y[k] - rho*GY(X, Y, k)/GYY(X, Y, k);
ErrY[k]:= Yk - Y[k];
Y[k]:= Yk
end do
end do;

## Matrix filled by recursive proc...

restart;
m:= proc(r::nonnegint, c::nonnegint)
option remember;
if r*c=0 then 0 else c*thisproc(r-1,c)+a*thisproc(r-1,c-1) end if
end proc:
m(0,0):= 1:

M:= Matrix(5,5, m);

Matrix indices start at 1, not 0, so the (trivial) 0th row and column are not shown. If you want indices that start somewhere other thn 1, you need to use an Array. But then you won't get the neat tabular layout in the GUI.

Your procedure myp can be reduced to

myp:= (a,b,c)-> zip(`+`, a, b)*c;

Is it possible to use something like
#map2(myp,[\$1..4],{[x,y,z,w],[x2,y2,z2,w2],[x3,y3,z3,w3]},[4,10,0.5]);
to achieve this?

seq(map2(myp, [\$1..4], {[x,y,z,w],[x2,y2,z2,w2],[x3,y3,z3,w3]}, s), s= [4,10,0.5]);

(I don't know why the line above breaks bad. It appears as one line in the editor, and I've tried several times to fix it.)

Is there a document I can read apart from the ?map to look for some general tips to 'simutanous' computing both symbolically and numerically in Maple?

See ?operators,elementwise and ?curry.

And is there any documentation on using multicore CPUs?

## rtable_num_elems(..., All)...

When applied to an rtable (Matrix, Vector, Array), nops is not the number of elements; rather, it is the number of operands in the internal structure. You can see the internal structure of A via op(A). For the number of elements, use rtable_num_elems(A, 'All').

## Some plots...

I'd like to show you the relevant commands, but I also want you to do your homework. So I'll compromise by posting the following without comments, most output elided, and without an attached file.

restart;
X:= t-> 3*cos(t):  Y:= t-> 2*sin(t):  t0:= 0.5:
kappa:= eval(VectorCalculus:-Curvature(<X(t),Y(t)>), t= t0):
dist:= (X1,X2)-> sqrt((X1[1]-X2[1])^2+(X1[2]-X2[2])^2):
slope:= (D(Y)/D(X))(t0):
x0:= X(t0):  y0:= Y(t0):
eq1:= (Cy-y0)/(Cx-x0) = -1/slope:
eq2:= dist([Cx,Cy],[x0,y0]) = 1/kappa:
fsolve({eq1,eq2}, {Cy,Cx});
{Cx = 4.13904333495107, Cy = 2.19319067267216}
fsolve({eq1,eq2}, {Cx,Cy}, avoid= %);
{Cx = 1.1264520363911743, Cy = -0.27548851825535093}
assign(%);
tan_eq:= slope*(x-x0)+y0:
Ellipse:= plot([X,Y,0..2*Pi], color= blue, scaling= constrained, thickness= 2):
Circle:= plot([Cx+cos(t)/kappa, Cy+sin(t)/kappa, t= 0..2*Pi], color= red):
Center:= plot([[Cx,Cy]], style= point, symbolsize= 25, symbol= cross, color= red):
Tangent:= plot(tan_eq, x= 0..4, color= green, thickness= 2):
plots:-display([Ellipse, Circle, Center, Tangent]);

## Missing a factor of V...

There are some problems with your terminology. You say gg is a vector. It doesn't make sense to raise a vector to a power. Perhaps g is a vector and gg represents the vector's dot product with itself, which would be a scalar. What you call the "gradient of psi" is the derivative of V with respect to psi. You can get this as

diff(V, psi);

But that will not express the derivative in terms of V.

You expression GVpsi is off by a factor of V. In other words, the numerator of the coefficient in front of the expression should be V instead of 1. It is difficult to get Maple to produce this expression as the derivative of V. However, you can verify it by subtracting diff(V, psi) and simplifying, obtaining 0. The expression can be derived "by hand" by a process called "logarithmic differentiation", which you can look up in a calculus textbook, which is often used to take the derivative of complicated expressions of the form f(x)^g(x).

 > restart;
 > V1:= (1-alpha+alpha*g^(1-1/psi))^(1/(1-1/psi));

 > dV:= diff(V1, psi);

 > dV2:= V*(alpha*V^(1/psi-1)*g^(1-1/psi)*ln(g)-ln(V))/psi^2/(1-1/psi);

 > eval(dV2, V= V1);

 > simplify(% - dV, symbolic);

 >

## Should be parentheses only...

All algebraic grouping in Maple, no matter how deeply nested, is done with ordinary round parentheses (), not square brackets [] or curly braces {}.

## A much simpler way: evalindets...

Trying to separate assigned and unassigned names can get very confusing, especially if you try to generalize it by putting the names in a set and then iterating over that. Here's a much simpler method to set certain SpecialNames to in an expression ex if they are not assigned.

SpecialNames:= {h||(251..256)}:
evalindets(ex, name, x-> `if`(x in SpecialNames, 0, x));

You can make up weird situations where this won't work, like if one name in the set is assigned to another name in the set. But you're not going to do that weird thing, right?

Let me know how that works for you.

 First 367 368 369 370 371 372 373 Last Page 369 of 389
﻿