Preben Alsholm

13578 Reputation

22 Badges

19 years, 172 days

MaplePrimes Activity


These are replies submitted by Preben Alsholm

Since you are able to use pointplot right after restart, i.e. you don't need plots:-pointplot, you must have a maple.ini file or some start up code.

So what is in that?

@mmcdara I agree.  In fact why should Maple do anything to 1h and 5 min. I myself like that better than 65 min.

@nm Very nice. I would, however, add one requirement:
 

if not has(eq,diff) or not T::`+` or not has(eq,yx) then return eq end if;

Your test case number 2 isn't an ode in y(x) so ought not be handled as such in my view.

The same can be written:

if not ( has(eq,diff) and T::`+` and has(eq,yx) ) then return eq end if;

@nm Number 3 and number 9 work as expected since the very first line in q:

if not has(eq,diff) then return eq end if;

makes those return unevaluated.
Number 7, however, certainly shows a weakness that should be handled.

@nm Thanks for making these points.
The first corrects a silly mistake on my part.
The second point is more tricky.
Here is an attempt:

p:=proc(eq::equation,{known::function:=NULL}) local d,fu,res;
  if not has(eq,diff) then return eq end if;
  d:=indets(indets(eq,specfunc(diff)),function) minus {known};
  fu:=op(remove(has,d,diff));
  if numelems({fu})>1 and known=NULL then 
    error "%1 has derivatives of more than one function. Use a second argument: known=function",eq end if;
  res:=selectremove(has,(lhs-rhs)(eq),fu);
  res[1]=-res[2]  
end proc:
########################
pde:=diff(u(x, t), t, t) + 3 + 2*diff(u(x, t), t) + 4*t + x^2 + x^3/3 + diff(u(x, t), t, x, x) + diff(u(x, t), x, x, x, x) = x*t^2;
ode:=y(x)+diff(y(x),x$2)+x=sin(x);
ode1:=y(x)+diff(y(x),x$2)+diff(f(x),x)=sin(x);
########
p(ode);
p(pde);
p(ode1);
p(ode1,known=f(x));

A much easier solution is to require a second argument, the function to solve for.

q:=proc(eq::equation,yx::function) local res;
  if not has(eq,diff) then return eq end if;
  res:=selectremove(has,(lhs-rhs)(eq),yx);
  res[1]=-res[2]  
end proc:
q(ode,y(x));
q(ode1,y(x));
q(pde,u(x,t));


 

For the fun of it I changed the ode's right hand side to 1/sqrt(Pi):

ode1:=diff(y(x),x)+x*y(x)=1/sqrt(Pi);
dsolve(ode1); # No Pi
ic:=y(0)=0;
Student:-ODEs:-ODESteps([ode1,ic])

That works fine.

@Carl Love Thanks Carl. I will use that too from now on.

@acer Thank you acer.

@acer I use that kind of thing in my maple.ini file , but in a version using sscanf called to my attention by Joe Riel a few years ago.
https://www.mapleprimes.com/maplesoftblog/209424-The-Maple-Command-Line-Interface-Now-In-Color

I need the digit after the period.
So I would use:
 

restart;
str:=convert(kernelopts(version),string);
### Also parse: 
N:=parse(str[7 .. StringTools:-Search(".",str)+1]); # i.e. +1 instead of -1

This makes it possible in my maple.ini file to distinguish between versions as in
if N> 2023.2 then .... end if.

So I have been using this:
 

restart;
VERSION:=kernelopts('version');
   N:=op(sscanf(VERSION, "Maple %f[0-9]"));  
   forget(evalf); #Otherwise the remember table will have an entry like 2019.2 = (2019.2, 10)
# I don't quite see the need for forget anymore.

 

@nm Try this:

restart;

sol:=y(x) = (exp(RootOf(-sin(x)*tanh(1/2*_Z+1/2*c__1)^2+sin(x)+exp(_Z)))+sin(x))/sin(x);
ode:=diff(y(x),x)-cot(x)*(y(x)^(1/2)-y(x)) = 0;

sol1,sol2:=allvalues(sol);

sol1:=simplify(expand(sol1));
sol2:=simplify(expand(sol2));
sol1_0:=eval(sol1,c__1=0);
sol2_0:=eval(sol2,c__1=0);
eval(rhs(sol1_0),x=-Pi/2); # 0
eval(rhs(sol2_0),x=-Pi/2); # 4
dsolve({ode,y(-Pi/2)=0});# y(x) = 0
solDs2:=dsolve({ode,y(-Pi/2)=4});
solDs2:=simplify(value(solDs2)) assuming x<0,x>-Pi;
plot(rhs(solDs2)-rhs(sol2_0),x=-Pi..0,style=point); # basically zero
plot(rhs(solDs2),x=-Pi..0,0..15,color=red); p1:=%: 
plot(rhs~([sol1_0,sol2_0]),x=-Pi..0,0..15); p2:=%:

plots:-display(Array([p1,p2]));

Interestingly you find two solution just as dsolve does, but one of the dsolve solutions is just the trivial y(x)=0.
Both yours and dsolve's can be correct since the initial value problem with y(-Pi/2) = 0 doesn't satisfy the requirements for uniqueness.
Note added:
Your sol1_0 isn't in fact a solution to the ode in the interval I consider here.
The value of sol1_0 at x=-Pi/4 is 1 - 2*2^(1/4) + sqrt(2) which is positive.
Using dsolve with that initial condition gives a solution that is valid on the interval
-27*Pi/32 ..-5*Pi/32.
On that interval sol1_0 and the dsolve solution differ:

@jestrup Sorry, now I get it. 
By "This could be relevant when using Maple for a test"  you mean that you don't want students to be able to access the AI formula assistant during a test.
As long as the students have a copy of Maple on their computers, I don't see any way you can do that except maybe for turning off your wi-fi.

@vv I tried the following for solve and for SolveTools:-SemiAlgebraic.
 

sol3:=solve({torus1<0,x>=0,y>=0,z>=0}, [x,y,z]);
numelems(sol3); # 3
sol2:=SolveTools:-SemiAlgebraic({torus1<0,x>=0,y>=0,z>=0}, [x,y,z]);
numelems(sol2); # 3
sol3[1];
sol2[1]; #same as above
sol3[2];
sol2[2]; #same as above
sol3[3];
sol2[3]; #same as above

sol3 and sol2 appears to say the same pairwise.

I understand the output from solve in this way:
 

restart;
torus:= (x^2+y^2+z^2 + R^2-r^2)^2 - 4*R^2*(x^2+y^2):
torus1:=eval(torus,[r=1,R=4]);
sol1:=solve(torus1<0, [x,y,z]); # It should be easy for Maple
numelems(sol1);
sol1[1]; # x = 0 is a requirement for this solution
sol1[2]; # x = 0 is also a requirement for this solution
### Thus this doesn't make any sense:
eval(sol1, [x=3,y=2,z=0]);

Note: From the output of eval(torus1, [x=3,y=2,z=0]) , however,  we can conlude that solve didn't find all solutions at all.

It works for me with a second order ode. In that case the equation is not turned into a system of first order.
 

restart;
ode:=diff(y(x),x)=0;

the_output:=Student:-ODEs:-ODESteps(ode,y(x));

latex(the_output)
########### 2nd order
ode := 3*diff(y(x), x, x) + 4*diff(y(x), x) + 2*y(x) = 0;

the_output:=Student:-ODEs:-ODESteps(ode,y(x));
latex(the_output);

MaplePrimes24-04-14_latex.mw

Note: If you continue the worksheet with the third order ode:
diff(y(x), x, x, x ) + 3*diff(y(x), x, x) + 4*diff(y(x), x) + 2*y(x) = 0;
Then after the "Too many levels of recursion" error the worksheet is ruined.
Try removing all output ( Ctrl+D), then after restart you are NOT back to normal.
You must copy the input code to a new worksheet.

I tried starting from x = 1 and found that for a solution to exist at zero in the limiting sense x -> 0 right,
it is necessary that D(y)(1) = 2*y(1).

restart;
ode:=x^2*diff(y(x),x$2)-2*y(x)=0;
sol1:=rhs(dsolve({ode,y(1)=y1,D(y)(1)=dy1}));
expand(sol1);
u:=-1/3*dy1/x + 2/3*y1/x; # The last 2 terms
expand(u*x*3);
solve(%,{y1});
sol12:=rhs(dsolve({ode,y(1)=dy1/2,D(y)(1)=dy1}));
plot(eval(sol12,dy1=2),x=0..1);
plots:-animate(plot,[sol12,x=0..1],dy1=-2..2);
1 2 3 4 5 6 7 Last Page 1 of 227