Preben Alsholm

13598 Reputation

22 Badges

19 years, 197 days

MaplePrimes Activity


These are replies submitted by Preben Alsholm

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

@Aung Why do you insist on turning
2*hypergeom([1/2, -k/2], [1 - k/2], csc(omega*T)^2);
into standard functions?

If you have any reason for using that expression it should be OK without turning it into standard functions.

1 2 3 4 5 6 7 Last Page 2 of 228