## 11749 Reputation

15 years, 295 days

## Problematic start...

Since

eval(rhs(eqn13),{f(z)=1,e(z)=1,z=0});
evaluates to

-6210.378246*I

i.e. an imaginary number, you may want to rethink the mathematics of the problem. Is there a wrong sign under the squareroot in eqn13?

Preben Alsholm

## Left hand side not a name...

You are trying to assign to a sequence of sums.

In the help page for assign it says that you can assign to names or functions, where 'function' means an unevalated  function call, like f(8).

Preben Alsholm

## Do the circles differently...

If I understand your intention correctly, then one solution is to draw only parts of the circles. Find find the relevant angles and use a parametric representation, as I have done below.

I have not changed much of your original code.

restart;
with(plots):
S1:=unapply(u*FresnelS(u)+((1/Pi)*cos((1/2)*Pi*u^2))-(1/Pi),u);
C1:=unapply(u*FresnelC(u)-((1/Pi)*sin((1/2)*Pi*u^2)),u);
R1:=1;
u1:=evalf(sqrt(0.4));
E1:=evalf(Pi*R1*S1(u1)+R1);
D1:=evalf(Pi*R1*C1(u1)+1);
a:=evalf(Pi*R1*(u1));
p1:=plot([1+a*FresnelC(u),a*FresnelS(u)-E1,u=0..u1]):
p2:=plot(-1.064877386,x=0..1):
p3:=plot(1.064877386,x=0..1):
p4:=plot([1+a*FresnelC(u),-a*FresnelS(u)+E1,u=0..u1]):
p5:=plot(-1.064877386,x=-1..0):
p6:=plot(1.064877386,x=-1..0):
p7:=plot([-1+(-a*FresnelC(u)),a*FresnelS(u)-E1,u=0..u1]):
p8:=plot([-1+(-a*FresnelC(u)),-a*FresnelS(u)+E1,u=0..u1]):
X:=1+a*FresnelC(u1);
t1:=fsolve(D1+R1*cos(t)=X,t);
C1:=plot([D1+R1*cos(t),R1*sin(t),t=-t1..t1]):
C2:=plot([-D1-R1*cos(t),R1*sin(t),t=-t1..t1]):
display(p7,p1,p2,p3,p4,p5,p6,p8,C1,C2,scaling=constrained);

Preben Alsholm

## Yet another way...

You already got some excellent solutions.

Here is another.

eq:=diff(y(x),x,x)+y(x)^2*x^2=x^2;
a:=[-0.6,-0.4,2.4,3.4];
sol:=y1->subs(dsolve({eq,y(0)=0,D(y)(0)=y1},numeric,output=listprocedure),y(x)):
#sol(y1) is the numerical procedure for evaluating y(x) (with D(y)(0)=y1), so that
sol(-0.6)(0.12345);
#will return the value of y(0.12345)
#and
plot([seq(sol(a[i]),i=1..4)],0..2);
# will plot the 4 different solutions in one plot.

Preben Alsholm

## Homework...

Do you expect us to do your homework?

(Excuse me if I was mistaken).

Preben Alsholm

## Same index in two different sums...

Two sums inside each other are using the same summation index. This leads to questions of interpretation.

If you change the expression to

f := Sum( exp(-beta*(2*x2-1)-alpha*(2*x2-1))
/ Sum( Sum( exp(-beta*(4*x0*x1-2*x0-2*x1+1)-alpha*(2*x0-1)), x1=0..1 ),
x0=0..1 ),
x2=0..1 );
where I have replaced x0 in three places by x2 the interpretation is clear and there is no problem.

Preben Alsholm

## Use subs...

subs(solutions,P);

will give you -6.960403739*10^(-13).

Preben Alsholm

## Use dsolve/numeric with infinity = 20...

You could use dsolve/numeric with infinity = 20 (or whatever):

m := 1;
de := diff(y(x), x, x, x)+((m+1)*(1/2))*y(x)*(diff(y(x), x, x))+m*(1-(diff(y(x), x))^2) = 0;
L:=dsolve({de,y(0) = 0, D(y)(0) = 0, D(y)(20) = 1},numeric);
plots:-odeplot(L,[x,y(x)],x=0..8);

The graphs on your photo are not graphs of y(x), but you can use odeplot to plot functions of y(x) also.

Preben Alsholm

## Missing 'and' and commas ...

You need 'and' between the boolean expressions. Also a few commas in the print statements were forgotten.

Here is the corrected version:

local det;
det:=b^2-4*a*c;
if a=0 and b=0 and c=0 then print( "All x is solution")
elif a=0 and b=0 and c<>0 then print ("No solution")
elif a=0 and b<>0 and c<>0 then print("One solution:",x=-c/b)
elif det>0 then print("Two real solutions:",x1=(-b+sqrt(det))/(2*a),x2=(-b-sqrt(det))/(2*a))
elif (det<0) then print("Two complex solutions:",x1=(-b+sqrt(det))/(2*a),x2=(-b-sqrt(det))/(2*a))
else print("Two identical solutions:",x=-b/(2*a))
end if;
end proc;

Be aware that your procedure is designed for numerical input only. Thus quad(r,s,t); doesn't work unless r, s, and t are assigned numerical values.

Preben Alsholm

## Random hermitian matrices...

Here is one way,

with(LinearAlgebra):
G:=RandomTools:-Generate(complex(distribution(Normal(0,.5))), makeproc=true);
Gd:=RandomTools:-Generate(distribution(Normal(0,.5)), makeproc=true);
A:=RandomMatrix(3,3,generator=G,outputoptions=[shape=hermitian]);
B:=RandomMatrix(3,3,generator=Gd,outputoptions=[shape=diagonal]);
C:=A+B;
Eigenvalues(C);

Preben Alsholm

## Looks like an exercise I...

Looks like an exercise I could have given my first semester students.
Maple sees xysin as a variable name, xysin(x-y) is an application of the function xysin to x-y.
Maple doesn't know any function by that name, so cannot plot it.

You need two multiplication signs.

Why grid = [2,2]?

Preben Alsholm

## The integral is a problem...

It seems that the basic problem is the presence of an integral. Besides, in the ODE you have D(phi(t))(t).

I have corrected that, so here is what I have,

restart;
dtOM:= (phi,a)-> -alpha*M^(alpha+4)/phi^(alpha+1) + 1/(3*Pi^2)*int(2*g^2*phi^2*k^4*(a*sqrt*(g^2*phi^2+k^2)/T+1)*exp(-a*sqrt*(g^2*phi^2+k^2)/T)*(g^2*phi^2+k^2)^(-3/2), k= 0..infinity):
RhoDE:= (phi,a)->M^(4(alpha+1))/phi^alpha+D(phi)/2*a^2:
RhoDM:=a-> C/a^3:
RhoR:=a->Pi^2/15*(T/a)^4:
alpha:=1: M:=1: g:=1: T:=1: C:=1:
RhoDE(phi,a)(t);
RhoDM(a(t));
RhoR(a(t));
dsys:={D(a)(t)=sqrt(RhoDE(phi,a)(t)+ RhoDM(a(t)) + RhoR(a(t)))*a(t),
D(D(phi))(t) + 3*a(t)*D(phi)(t)/D(a)(t)*dtOM(phi(t),a(t))=0,
a(1000)=2,
phi(1000)=3,
D(phi)(1000)=4
};
dsolve(dsys,[a(t),phi(t)], numeric);

resulting in a warning and an error message:

Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters)
Error, (in DEtools/convertsys) unable to convert systems with integrals containing the dependent variables

The warning may have to do with the variable of integration. The error message is self-explanatory.

Preben Alsholm

## Powers ARE too high...

By looking at the procedure 'residue' in line 8 you will see the loop

for i from 0 to 5 do

If you substitute for 5 a higher number, you get an answer.

In experimenting with this be aware that some things are remembered, so a restart can be necessary.

Try the following, where 'residue2' is simply 'residue' with 5 replaced by 6.

restart;
residue((x^7+1)^4/x^30, x=0);
showstat(residue);
residue2:=subs(5=6,eval(residue));
showstat(residue2);
debug(residue2);
residue2((x^7+1)^4/x^30, x=0);
restart;
residue2:=subs(5=6,eval(residue));
residue2((x^7+1)^4/x^30, x=0);

You may even make yourself a new procedure that takes the top loop variable as an optional third argument:
I have used 2 forget commands to avoid a restart and to make behavior consistent.

restart;
Residue:=proc(f, a::(name = anything),n::posint:=5)
forget(residue); forget(series);
if n=5 then residue(f,a) else subs(5=n,eval(residue))(f,a) end if
end proc;
Residue((x^7+1)^4/x^30, x=0);
Residue((x^7+1)^4/x^30, x=0,6);
Residue((x^7+1)^4/x^30, x=0);

Preben Alsholm

## Shifting from using hardware floats to s...

It sounds like with Digits <=15 computations are done with hardware floats, whereas if Digits > 15 Maple has to use software floats.

See e.g. the help page for evalhf.

Preben Alsholm

 First 146 147 148 Page 148 of 148
﻿