Preben Alsholm

MaplePrimes Activity


These are replies submitted by Preben Alsholm

@Met28 You cannot solve the problem via Maple. It is mathematics. The improper integral is simply not convergent. I illustrated that by my Maple statements.

Using the usual definition of convergence in Calculus, you need to show that the improper integrals over the 3 intervals
0..10.69,  10.69..b, b..infinity (where b is chosen so that 10.69 < b < infinity)
are all 3 convergent.

Thus to prove divergence you only need to show that one of the 3 is divergent. As it happens, they are divergent all 3.

Divergence of the first can be proven like this:
res:=int((x*0.187549975)/(x^2-10.69^2), x=0..a) assuming a>0,a<10.69;
limit(res,a=10.69,left);


@Wackeraf I executed all the statements and then examined the odes only. I didn't check the physics, i.e.  whether your derivation of the 4 odes La1, La2, La3, La4 is correct or not. I took them as your code gave them.
That you cannot solve for the highest order derivatives can sometimes be solved by differentiating one of the equations, i.e. using non-algebraic means. It is not obvious to me that that is a way out here.

###
Notice, that in A actually two of the rows are the same:

A[1,..];
A[4,..];

Did you try computing the integral before attempting minimization?

int(0.1e-3+.5*t+0.2e-2*t^2-b*t-a, t = 0 .. 300);

Output:

@acer I tried changing initial condition in your version as I did in my comment above. It doesn't work.

Just checking existing conditions doesn't work either:

ONE_SOLUTION(initial);

The initial condition for D(h) can be included in the parameters, however:

ODE_SOLUTION := dsolve({ODE, h(0) = 0,D(h)(0)=y1}, numeric, method=rosenbrock_dae,parameters = [y1,A, B, C, E, F]):
ODE_SOLUTION('parameters'=[Pi/2,0,0,1,1,1]):






@ There is ambiguity in determining the initial condition even in the case you consider.
The initial condition for u, i.e. u(0)=u1, satisfies cos(u1)=0. Thus there is a choice to be made.
After having found the solution as given, we can do

sol(initial=[0,0,-Pi/2]);

and when plotting we find that we have another solution.


@Carl Love I couldn't agree more!

@Alejandro Jakubi Thank you very much, I shall try to contact Gaston Gonnet.

@Markiyan Hirnyk Yes, but you are using two calls to eval.
You cannot do
eval(f(a*x), x = a, a = 9);

@acer Using method=nonlinearsimplex worked in Maple 15 on my simple example, which is of the same kind as the one described by joha.

@acer I noticed that joha uses Maple 14. On a simple example I tested, this device doesn't work in Maple 15 (and thus most likely not in Maple 14). It does work in Maple 18 (and in Maple 16).

I think it might be necessary to use the operator form with a user supplied gradient as shown by you in 2010 (July 14).

@Benrath plot3d(x*y,x=y..1,y=0..1);

@nm I think you are being rather unreasonable.
When I do
?pi
it brings up a help page, which starts out saying
"For information on 3.14159..., see Pi."

I think that that is as helpful as it gets!

And that
?foo
brings up a help page about Units,length, (foot)
just shows that the help system is trying to be helpful.

In your worksheet you don't attempt any solution of the system (sys), but you are finding equilibria by using solve. I don't get any error at all, but the output:



i.e. infinitely many equilibria, here parametrized by v, x, and z.

If you want to solve sys, you will undoubtably have to do it numerically. For that the constants m and q have to be made concrete and initial (or boundary) values must be known.

@sharena2 Why do you think results will get better using shooting?

However, if you insist, I suggest shooting from eta=L, which will only require two parameters.
For comparison I continue from the code I gave above. Thus 'sol' is the output from dsolve/bvp.

#Shooting from the right with the 2 parameters fL and vL:
ICSL:={F(L) = 0, H(L) = 2, f(L) = fL, g(L) =-fL, u(L) = 0,v(L)=vL};
resL:=dsolve(ODE union ICSL,numeric,parameters=[fL,vL],output=listprocedure);
#We shall need these two:
ff,uu:=op(subs(resL,[f,u](eta)));
#Getting parameter values from sol (no shooting):
pp:=subs(sol(L),[f(eta),v(eta)]);
#Quick graphical test:
resL(parameters=pp):
plots:-display(Array([[seq(plots:-odeplot(resL,A[i],0..L),i=1..3)],[seq(plots:-odeplot(resL,A[i],0..L),i=4..6)]]));

#There are 2 parameters  (fL and vL) and 2 requirements: f(0) = 0 and u(0) = 1
pL:=proc(fL,vL) option remember;
  resL(parameters=[fL,vL]);
  [ff(0),uu(0)-1]
end proc:
p1:=proc(fL,vL) pL(_passed)[1] end proc:
p2:=proc(fL,vL) pL(_passed)[2] end proc:
#Testing pL on the parameters from sol (pp):
pL(op(pp)); #Should be close to zero
ppL:=fsolve([p1,p2],pp);
pp;
ppL-pp;
#With no so good a guess:
fsolve([p1,p2],[.4,0]);
resL(parameters=%):
plots:-display(Array([[seq(plots:-odeplot(resL,A[i],0..L),i=1..3)],[seq(plots:-odeplot(sol,A[i],0..L),i=4..6)]]));


You normally would do something like (example)

if not type(aaa,numeric) then error "%1 should be of type numeric",aaa end if;

First 133 134 135 136 137 138 139 Last Page 135 of 231