## 355 Reputation

11 years, 304 days

## Syntax issues...

There were quite a few syntax issues.

ode.mw

## piecewise...

The best is to use piecewise, see for details maple help

restart:with(plots):
p1 := proc (x, y) piecewise(x^2+y^2 <= 1 ,x*y-y^2,0) end proc;
plot3d(p1(x,y), x = -1 .. 1, y = -1 .. 1);

To avoid the error you are getting, just use single quotes

p1 := proc (x, y) if x^2+y^2 <= 1 then x*y-y^2 else 0 end if end proc;

plot3d('p1(x,y)', x = -1 .. 1, y = -1 .. 1);

## printlevel:=3...

As mentioned by Carl Love, your question is unclear. But if you want to see the ouput of the two loops then you need to use printlevel:=3 before the loops.

## Wrong analytical solution...

The so-called analytical solution you proposed is not satisfying the ode,

`Eq:=diff(w(x), x\$6)+diff(w(x), x\$4)+diff(w(x), x\$2)-w(x) = -90;`

`sol1:=w(x)=C1*sinh(x)+C2*cosh(x)+C3*sin(x)+C4*cos(x)+C5*sin(x)+C6*cos(x)+90;`

`odetest(sol1, Eq);`

## odetest...

To verify a solution to an ode, you can use odetest to do so.

 > restart:
 > #p1:=1:p2:=1:p3:=1:chi:=1:omega:=2:
 > Eq:=diff(w(x), x\$6)+A*(diff(w(x), x\$4))+B*(diff(w(x), x\$2))+(-chi*omega^2+C)*w(x) = 0;
 (1)
 > sol:=dsolve(Eq):
 > simplify(expand(odetest(sol, Eq)))
 (2)
 > H:=subs(36*A*B+108*chi*omega^2-108*C-8*A^3+12*sqrt(-12*A^3*chi*omega^2+81*chi^2*omega^4+54*A*B*chi*omega^2+12*A^3*C-3*A^2*B^2-162*C*chi*omega^2-54*A*B*C+12*B^3+81*C^2) = E, (4*I)*A^2*sqrt(3)-I*sqrt(3)*E^(2/3)-(12*I)*B*sqrt(3)-4*A^2-4*A*E^(1/3)-E^(2/3)+12*B = F, (4*I)*A^2*sqrt(3)-I*sqrt(3)*E^(2/3)-(12*I)*B*sqrt(3)+4*A^2+4*A*E^(1/3)+E^(2/3)-12*B = G, 4*A^2-2*A*E^(1/3)+E^(2/3)-12*B = S, rhs(sol)):
 > simplify(subs(E=36*A*B+108*chi*omega^2-108*C-8*A^3+12*sqrt(-12*A^3*chi*omega^2+81*chi^2*omega^4+54*A*B*chi*omega^2+12*A^3*C-3*A^2*B^2-162*C*chi*omega^2-54*A*B*C+12*B^3+81*C^2),subs(S=4*A^2-2*A*E^(1/3)+E^(2/3)-12*B ,subs(G=(4*I)*A^2*sqrt(3)-I*sqrt(3)*E^(2/3)-(12*I)*B*sqrt(3)+4*A^2+4*A*E^(1/3)+E^(2/3)-12*B,subs(E=36*A*B+108*chi*omega^2-108*C-8*A^3+12*sqrt(-12*A^3*chi*omega^2+81*chi^2*omega^4+54*A*B*chi*omega^2+12*A^3*C-3*A^2*B^2-162*C*chi*omega^2-54*A*B*C+12*B^3+81*C^2),subs(F=(4*I)*A^2*sqrt(3)-I*sqrt(3)*E^(2/3)-(12*I)*B*sqrt(3)-4*A^2-4*A*E^(1/3)-E^(2/3)+12*B,subs(E=36*A*B+108*chi*omega^2-108*C-8*A^3+12*sqrt(-12*A^3*chi*omega^2+81*chi^2*omega^4+54*A*B*chi*omega^2+12*A^3*C-3*A^2*B^2-162*C*chi*omega^2-54*A*B*C+12*B^3+81*C^2),H))))))):
 > sol1:=w(x)=%:
 > odetest(sol1, Eq)
 (3)
 >

## No way!...

As for as I know, the dsolve has already given you the general solution in compact form.

 > restart:
 > #p1:=1:p2:=1:p3:=1:chi:=1:omega:=2:
 > Eq:=diff(w(x), x\$6)+p1*(diff(w(x), x\$4))+p2*(diff(w(x), x\$2))+(-chi*omega^2+p3)*w(x) = 0;
 (1)
 > dsolve(Eq):
 > simplify(expand(rhs(%))):
 >

## Continuation options...

You can also opt for the continuation option.

 > restart; with(plots):
 > Digits:=20:
 > rho := 1: nu := .1:
 > B := (diff(psi(u), u\$2))/cosh(2*u);
 (1)
 > eq0 := (1000-999*lambda)*((diff(B, u\$2)+2*(diff(B, u))+B*(4+rho*sqrt(3)*(diff(psi(u), u))/nu)))=0;
 (2)
 > bc0 := (D(psi))(0.174568e-1) = 0, (D(psi))(-0.174568e-1) = 0, psi(-0.174568e-1) =0.00001, psi(0.174568e-1) = 70;
 (3)
 > sol := dsolve({bc0, eq0},numeric,continuation=lambda);
 (4)
 > odeplot(sol, [u, psi(u)], u = -0.174568e-1 .. 0.174568e-1)
 > odeplot(sol, [u, diff(psi(u),u)], u = -0.174568e-1 .. 0.174568e-1)
 >

## Maple 2016...

In Maple 2016, I have executed your sheet without any problem.

## contourplot(,options)...

You can use the many options available  with the 'contourplot' to improve the quaility.

Here I used the `contours=s` and `grid`

contourplot(-(1/2)*y^2-(1/2)*x^2-(1-.3)/sqrt((x+.3)^2+y^2)+((-1)*.3)/sqrt((x-1+.3)^2+y^2), x = -1.5 .. 1.5, y = -1.5 .. 1.5,grid=[60,60],contours=120);

You can get the 2d look from the 3d just by using orientation option

contourplot3d(-(1/2)*y^2-(1/2)*x^2-(1-.3)/sqrt((x+.3)^2+y^2)+((-1)*.3)/sqrt((x-1+.3)^2+y^2), x = -1.5 .. 1.5, y = -1.5 .. 1.5, view = -2 .. -1.3, axes = boxed,style=patchnogrid, grid=[160,160], orientation=[90,360], lightmodel=light4, shading=zhue,  transparency=0.3);

## Guide Maple...

As mentioned by @ Kitonum, there is no built-in Maple command to find Fourier series of a function.

But you can tell Maple to help you avoid all the integration and summation.

To find a Fourier series of a function, first we need to find the Euler-Fourier formulas (a_m and b_n) and then sum it.  This can be done by writting three procedures like this

 > restart:
 > a := proc(f,k)    int(f(x)*cos(k*Pi*x/L),x=-L..L)/L; end;
 (1)
 > b := proc(f,k)    int(f(x)*sin(k*Pi*x/L),x=-L..L)/L; end;
 (2)
 > fss := proc(f,n)       a(f,0)/2+sum(a(f,k)*cos(k*Pi*x/L)+b(f,k)*sin(k*Pi*x/L),k=1..n);    end;
 (3)
 > f := x ->x-floor((x+2)/4)*4;L:=2:
 (4)
 > a(f,k);
 (5)
 > b(f,k);
 (6)
 > fss(f,2);
 (7)
 > plot({f(x),fss(f,50)},x=-12..12,numpoints=1000,color=[red,green], discont,linestyle=[1,3]);
 >

This not my own idea. I got it somewhere on the internet.

## plot3d or implicit plot...

restart:with(plots):

x:=2*m*p/(m^2+n^2+p^2);

y:=2*n*p/(m^2+n^2+p^2);

z:=(p^2-(m^2+n^2))/(m^2+n^2+p^2);

p:=2:

plot3d([x,y,z],n=-1..1,m=1..1,grid=[30,30],scaling=constrained,axes=boxed);

f:=x^2+y^2=z^2;

implicitplot(f,n=-1..1,m=-1..1);

## Input in maple...

Just type it down in maple

restart:
B[d]:=(-d*w[1]+w[3])/(((-2*w[1]+w[3])^2)-((-d*w[1]+w[3])^2)):

B[o]:=(-2*w[1]+w[3])/(((-2*w[1]+w[3])^2)-((-d*w[1]+w[3])^2)):

A[o]:=w[1]*(alpha[o]-c[o]-t[o])+2*w[2]*e[o]:

A[d]:=w[1]*(alpha[d]-c[d]-t[d])+2*w[2]*e[d]:

q[o]:=B[o]*A[o]-B[d]*A[d]:

q[d]:=B[o]*A[d]-B[d]*A[o]:

U[d]:=w[1]*(q[d]*(alpha[d]-q[d]-d*q[o]-c[d]-t[d])-C)+w[2]*(e[d]*q[d]+e[o]*q[o])+w[3]*((1/2)*(q[d]+q[o])^2);

## Something like this?...

restart:

for i to 2 by 1 do
n[i] := (1/2)*i
end do;

for j to 2 by 1 do
eq[j] := diff(f(eta), eta, eta)+(1/(4*eta^2)+(3.*n[j]+1)/(n[j]+1)*(NU/(4*a*eta)-NU*(eta/a)^((n[j]+1)/(2.0*n[j]))/(4*a*eta)))*f(eta) = 0:
end do;

for j to 2 do
dsolve(eq[j]);
end do;

## dsolve?...

As @Rouben mentioned, you need to post the ODEs you want an answer to.

Checking the maple help page for dsolve, I found an example which might help

Define a system of ODEs.

sys_ode := diff(y(t),t) = x(t), diff(x(t),t) = -x(t);

If the unknowns are not specified, all differentiated indeterminate functions in the system are treated as the unknowns of the problem.

dsolve([sys_ode]);

{x(t) = _C2 exp(-t), y(t) = -_C2 exp(-t) + _C1}

Define initial conditions.

ics := x(0)=1, y(1)=0;
ics := x(0) = 1, y(1) = 0

Solve the system of ODEs subject to the initial conditions ics.

dsolve([sys_ode, ics]);

You can then extend for a system of second order ODE's fairly easy. While trying for that if maple is unable to solve analytically your unseen problem then you can use numeric.

GoodLuck

## Might help...

 > restart:
 > ODE:=diff(f(eta),eta,eta,eta)-(f(eta)*diff(f(eta),eta,eta)+beta*(1-diff(f(eta),eta)^2))=0;
 (1)
 > BCS:=f(0)=0,D(f)(0)=0,D(f)(infinity)=-1;
 (2)
 > d1:= (f[i+1]-f[i-1])/(2*h);
 (3)
 > d2:=(f[i+1]-2*f[i]+f[i-1])/h^2;
 (4)
 > d3:= (f[i+1]-3*f[i]+3*f[i-1]-f[i-2])/h^3;

 (5)
 > ode1:=d3-f[i]*d2-beta(1-d1^2)=0;
 (6)
 > beta:=0.5;
 (7)
 > h:=0.5;
 (8)
 > ode1;
 (9)
 > eta[0]:=0;
 (10)
 > for i from 0 to 9 do eta[i+1]:=eta[i]+h; end do:

Using the BCS

 > f[0]:=0:
 > f[-1]:=f[1]:
 > eq[10]:=((f[10]-f[9])/(h))=-1:
 > for i from 0 to 9 do eq[i]:=ode1; end do;
 (11)
 > soln:=fsolve({eq[0],eq[1],eq[2],eq[3],eq[4],eq[5],eq[6],eq[7],eq[8],eq[9],eq[10]},{f[-2],f[1],f[2],f[3],f[4],f[5],f[6],f[7],f[8],f[9],f[10]});
 (12)
 > assign(soln);
 > data:=[seq([eta[0]+i*h,f[i]],i=0..10)];
 (13)
 >
 > plot([data],x=eta[0]..10);
 >