Preben Alsholm

13743 Reputation

22 Badges

20 years, 342 days

MaplePrimes Activity


These are replies submitted by Preben Alsholm

By saying "I have written a demonstration worksheet to show this problem", I presume that you intended to include that worksheet in your question, but I don't see any link to it.

@adel-00 If you just want to plot the quantity 0.25*z(t)^2+Re(y(t))^2 then you can just use res from before:

plots:-odeplot(res,[t,0.25*z(t)^2+Re(y(t))^2],0..2);
# but if you want to get actual numbers you can use output=listprocedure, which would be OK to use in any case.
#Notice that this time I have closed the assignments to resL, Y, Z by a colon. If you replace colon with semicolon you will know why (a huge output: it is usually hidden. It is surely not intended. I shall report it as a minor bug).

resL:=dsolve(dsys union {x(0)=0,y(0)=0,z(0)=-1/2},numeric,output=listprocedure):
#You can use odeplot as before:
plots:-odeplot(resL,[t,0.25*z(t)^2+Re(y(t))^2],0..2);
Y,Z:=op(subs(resL,[y(t),z(t)])):
Y(0.12345);
Z(0.12345);
#You could if you wish define the quantity as W:
W:=t->0.25*Z(t)^2+Re(Y(t))^2;
W(.987654);
#And you can plot W (same plot as above):
plot(W(t),t=0..2);




@RAfossey Here is a stepwise solution by laplace transform. As you are doing I'm assigning to a,b,c,d immediately although that doesn't seem to be what the problem formulation called for.

restart;

with(inttrans):
a:=1;b:=1;c:=1;d:=1;
sys:=diff(x(t),t,t)=a*x(t)+b*y(t),diff(y(t),t,t)=c*x(t)+d*y(t);
ics:=x(0)=1,D(x)(0)=0,y(0)=0,D(y)(0)=1;
#Applying laplace elementwise (notice the tilde ~ ) to the system:
laplace~([sys],t,s);
eval(%,{ics});
#Solving for the laplace transforms:
L:=solve(%,{laplace(x(t),t,s),laplace(y(t),t,s)});
#Elementwise inverse transformation
sol:=invlaplace~(L,s,t);
#Extracting the solutions:
X,Y:=op(subs(sol,[x(t),y(t)]));
#Plots:
plot([X,Y],t=0..2);
#This is a parametric plot:
plot([X,Y,t=0..2]);


@wenny You should try executing the worksheet again. I did not get that error message.
If you get that error once again then try this:

indets({g1,g2,g3},name);

you should be seeing the unassigned names and hopefully what you see are just the variables.

@ I'm apalled to learn that some math professor out there cares a hoot about how you indicate the end of proof as long as you do it!

@jonlg Maple simply cannot do the integral even when T1, T2, U0, L0 are given concrete values.

I dropped the case t<=0 because I thought (perhaps wrongly) that you had no interest in t<0.

piecewise works like this: if the first condition is satisfied then the following expression is returned. If not then the next condition is evaluated, etc. If none of these is satisfied then 0 is returned. That 'otherwise' case is specified explicitly in my example as the last argument; when given explicitly it can be anything.

Could you give us a (simple) example of what you have in mind, thereby answering questions like 'Does the parameter have to be determined also or can it chosen arbitrarily' ?

@Konstantin@ A common (and recommended) method is to include the integrals in the system.

Suppose as a simple example that you want to find the integral int(z(x)*q(x),x=0..a), where q is a known function.
Then you could include the ode diff(w(x),x)=z(x)*q(x) with w(0)=0 in the system:

sys:={diff(z(x), x$2)+z(x) = -2, diff(w(x),x)=z(x)*q(x), z(0)=0, z(1)=1, w(0)=0}.

Then the matrix returned with the option=Array(...) will include w-values, and of course w(a) will be your integral.

Did you actually try this:

dsolve({eq,y(0)=1,D(y)(0)=0},y(x),type='series',order=10);

It works for me in Maple 18, Maple 15, and Maple 12. I'm pretty confident that it works in all versions from 12 and up. Maybe it also works even in some earlier releases. I don't have access to those on this machine.

As for series:

series(cos(x),x=0,10);

@Markiyan Hirnyk The system DE1, DE2 is clearly NOT autonomous, since rhs(DE2) depends on t explicitly through cos(t).
When you tried using DEtools[autonomous] you had a wrong independent variable at the end.

DEtools[autonomous]({DE1, DE2}, {y, z}, t); #t not x

@Konstantin@ I should have written 1 and 3, not 1 and 2. Notice that I go on quoting those points.

The puzzling part about 1 is that read for me works on .m files in the sense that something is read. 'read' doesn't produce an output, it just reads stuff. In this case this means that S1 is assigned to the list you see, when you write S1; after having executed the read statement.
I thought that your original problem was that indeed something is read, but that something (S1) doesn't work as expected. It is difficult to read your original question any other way.

@wenny You should use the numerical solver, fsolve instead of solve. The integrals are inert integrals (Int instead of int). When fsolve gets to work it will evaluate the integrals numerically.
You will have more of a chance of success if you provide starting values or ranges.

@Konstantin@ I use Maple 18.02 on Windows 7. Which Maple version are you using?

I just had occasion to repeat the recipe I gave and with the same result as mentioned earlier.

I'm puzzled by your points 1 and 2:

that you cannot read the file you saved as an .m file
and
that when you save with any other extension everything works. (!!?)

I tried with 2 different extensions: .mpl and .txt. It didn't make the read-in value of S1 any better. But this time the reason was another, I got the error you report in your point 3. So do you really mean that everything works for you if you save with a different extension than .m ?

I tried your ode with initial values instead of boundary values. There a similar remedy I attempted with the procedure `dsolve/numeric/SC/IVPdcopy` didn't succeed. But see the remark about the range option at the bottom.

If you try

F1(0.12345);

you get the error:

Error, (in F1) argument #2, must be the target rtable

That may be because in F1 there is an Array named _dat whose second element contains pointers to places in the memory. These pointers are unique to the session. (I hope this is a somewhat accurate description).
Thus when saved in one session and read into another they will not point to anything (relevant).

showstat(F1);

The first element of _dat (out of 4) is a (large) procedure. The remaining were in my session:
2 = (Array(1..3, [4466537074,4466537218,4466537314])), 3 = [x, z(x), diff(z(x),x)], 4 = []

pointto(4466537074); #Error

If you try the same in the session where you save S1 you will see that the similar pointto command outputs a procedure. The 3 numbers refer to the procedures for x, z(x), diff(f(x),x) respectively.

If you take another look at the boundary value problem you will see that in F1 there is an Array called data quite like _dat. Its second element is an Array of three pointers (addresses).

That the copying device described above works and the pointers are apparently not needed after that is interesting.
#### The range argument in the ivp-problem
Trying with the range argument in the ivp problem showed that indeed a similar device worked here. Executing and copying in a fresh worksheet the procedure `dsolve/numeric/SC/IVPval` and assigning the pasted version to `dsolve/numeric/SC/IVPval` in our "read"-session was successful.

CONCLUSION: To avoid the trouble use output=Array(...) if you need to save data between sessions.

@wenny Certainly it is possible that there is a solution in some other range than {gam0=0..230, gam1=0.1..10}.

So do you have any guesses as to where a solution (if any exists) should lie?
Should gam0 and/or gam1 be positive?

The reason I tried this:

evalf(eval([g1,g2],{gam0=0.01,gam1=3.62}));

was simply to check that actual numbers came out.
If I couldn't find some gam0 and gam1 for which numbers were produced there would be no reason to try fsolve.

Does the problem come from some application (in physics or engineering e.g.)?

First 130 131 132 133 134 135 136 Last Page 132 of 231