imparter

155 Reputation

8 Badges

12 years, 199 days

MaplePrimes Activity


These are questions asked by imparter

How to display  3 figures (matrixplot(A),matrixplot(B),matrixplot(C)) in a row  

plots:-matrixplot(A):

plots:-matrixplot(B):

plots:-matrixplot(C):

thanks in advance

Dear maple user ,while solving the system of pdes i am getting the errors and unable to plot the curves. please help me to  rectify the errors and  plot the curves.

Thanks in advance

 

restart;
with(PDETools):
G:=1300:beta:=0.0075:epsilon:=4.9*(10)^(-6):Pa:=100:alpha:=24:Pv:=97:
sys1:= {G*(beta^2*u(r,z)*diff(u(r,z),z)+v(r,z)*diff(u(r,z),r)) = -diff(p(r,z),z) +beta^2*diff(u(r,z), z$2)+diff(u(r,z),r$2)+(1/r)*diff(u(r,z),r),G*(beta^2*u(r,z)*diff(v(r,z),z)+v(r,z)*diff(v(r,z),r)) = -diff(p(r,z),r) +beta^2*diff(v(r,z), z$2)+diff(v(r,z),r$2)+(1/r)*diff(v(r,z),r)-(v(r,z)/r^2),beta^2*diff(u(r,z),z)+diff(v(r,z),r)+(v(r,z)/r)=0  }:

IBC:=(D[1](u))(0, z) = 0,v(0,z)=0,phi*(D[1](u))(1, z)+u(1,z) = 0,v(1,z)=epsilon*(p(1,z)+(Pa/alpha)-1),p(r,-1)=0,p(1,z)=(Pv-Pa)/alpha:
sol := pdsolve(sys1, IBC, numeric):
r:=0:

p1 := sol:-plot(u(r, z), phi = 0, numpoints = 100, z = -1 .. 1, color = ["Blue"], legend = ["u(r,z)"]):

p2 := sol:-plot(u(r, z), phi = 0.15, numpoints = 100, z = -1. 1, color = ["red"], legend = ["u(r,z)"]):

p3 := sol:-plot(u(r, z), phi = 0.4, numpoints = 100, z = -1 .. 1, color = ["green"], legend = ["u(r,z)"]):

plots:-display({p1, p2,p3});

 

Any one help me  to remove the error.

I want to plot the curve for different values of alpha.

here is my codes.

thanks in advance 

 

restart:
with(linalg):with(plots): 
ge[1]:=diff(u[1](x,t),t)=alpha*diff((u[2](x,t)-1)*diff(u[1](x,t),x),x)+(16*x*t-2*t-16*(u[2](x,t)-1))*(u[1](x,t)-1)+10*x*exp(-4*x):
ge[2]:=diff(u[2](x,t),t)=diff(u[2](x,t),x$2)+alpha*diff(u[1](x,t),x)+4*(u[1](x,t)-1)+x^2-2*t-10*t*exp(-4*x): 
bc1[1]:=u[1](x,t)-1: 
bc1[2]:=u[2](x,t)-1: 
bc2[1]:=3*u[1](x,t)+diff(u[1](x,t),x)-3:
bc2[2]:=5*diff(u[2](x,t),x)-evalf(exp(4))*(u[1](x,t)-1):
IC[1]:=u[1](x,0)=1: 
IC[2]:=u[2](x,0)=1: 
NN:=2: 
N:=2:
L:=1:
for i to NN do  
dydxf[i]:=1/2*(-u[2,i](t)-3*u[0,i](t)+4*u[1,i](t))/h: 
dydxb[i]:=1/2*(u[N-1,i](t)+3*u[N+1,i](t)-4*u[N,i](t))/h:
dydx[i]:=1/2/h*(u[m+1,i](t)-u[m-1,i](t)); 
d2ydx2[i]:=1/h^2*(u[m-1,i](t)-2*u[m,i](t)+u[m+1,i](t)):od:
 for i to NN do bc1[i]:=subs(diff(u[1](x,t),x)=dydxf[1],
diff(u[2](x,t),x)=dydxf[2],u[1](x,t) 
=u[0,1](t),u[2](x,t)=u[0,2](t),x=0,bc1[i]):od: 
for i to NN do bc2[i]:=subs(diff(u[1](x,t),x)=dydxb[1],
diff(u[2](x,t),x)=dydxb[2],u[1](x,t) 
=u[N+1,1](t),u[2](x,t)=u[N+1,2](t),x=L,bc2[i]):od:
for i to NN do eq[0,i]:=bc1[i];eq[N+1,i]:=bc2[i]:od: 
for i from 1 to N do eq[i,1]:=diff(u[i,1](t),t)= subs(diff(u[1](x,t),x$2) =
subs(m=i,d2ydx2[1]), 
diff(u[2](x,t),x$2) = subs(m=i,d2ydx2[2]),diff(u[1](x,t),x) =
subs(m=i,dydx[1]),diff(u[2](x,t),x) = subs(m=i,dydx[2]),u[1](x,t)=u[i,1](t), 
u[2](x,t)=u[i,2](t),x=i*h,rhs(ge[1])):od:

for i from 1 to N do eq[i,2]:=diff(u[i,2](t),t)= subs(diff(u[1](x,t),x$2) =
subs(m=i,d2ydx2[1]), 
diff(u[2](x,t),x$2) = subs(m=i,d2ydx2[2]),diff(u[1](x,t),x) =
subs(m=i,dydx[1]),diff(u[2](x,t),x) = subs(m=i,dydx[2]),u[1](x,t)=u[i,1](t),
u[2](x,t)=u[i,2](t),x=i*h,rhs(ge[2])):od: 

for i to NN do u[0,i](t):=(solve(eq[0,i],u[0,i](t))):od:

 for i to NN do u[N+1,i](t):=(solve(eq[N+1,i],u[N+1,i](t))):od:

 h:=L/(N+1): 

for i from 1 to N do eq[i,1]:=eval(eq[i,1]):od: 
 for i from 1 to N do eq[i,2]:=eval(eq[i,2]):od:

eqs:=seq(seq((eq[i,j]),i=1..N),j=1..NN): 
Y:=seq(seq(u[i,j](t),i=1..N),j=1..NN): 

 ICs:=seq(u[i,1](0)=rhs(IC[1]),i=1..N),seq(u[i,2](0)=rhs(IC[2]),i=1..N): 

sol:=dsolve({eqs,ICs},{Y},type=numeric,stiff=true,maxfun=1000000,abserr=1e-6,relerr=1e-5,output=listprocedure):

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)
for j to NN do for i to N do U[i,j]:=subs(sol,u[i,j](t)):od:od: 

for i to NN do U[0,i]:=subs(u[1,1](t)=U[1,1],u[1,2](t)=U[1,2],
u[2,1](t)=U[2,1],u[2,2](t)=U[2,2],u[0,i](t)):od:
 for i to NN do U[N+1,i]:=eval(subs(u[N,1](t)=U[N,1],u[N,2](t)=U[N,2],
u[N-1,1](t)=U[N-1,1],u[N-1,2](t)=U[N-1,2],u[N+1,i](t))):od:
tf:=1.: 
M:=30: 
T1:=[seq(tf*i/M,i=0..M)]: 
PP:=matrix(N+2,M+1): 
for i from 1 to N+2 do PP[i,1]:=evalf(subs(x=(i-1)*h,rhs(IC[1]))):od: 
for i from 1 to N+2 do for j from 2 to M+1 do
PP[i,j]:=evalf(subs(t=T1[j],U[i-1,1](t))):od:od:
 
G1:=[seq([ seq([(i-1)*h,T1[j],PP[i,j]], i=1..N+2)], j=1..M+1)]: 
t=0.02: 
pars:=[0.1,0.5,1,2,5];
clr:=[black,red,green,gold,blue];  
for m from 1 to 5 do 
G1[m]:=plot([seq(subs(alpha=pars[m],G1[i])],i=0..N+1)],thickness=3,color=clr[j]):od:  

display({seq(G1[i],i=1..5)},title="Figure ",axes=boxed,labels=[x,u]);
restart:
 

Hai

Help required to plot the graphs for system  difference schemes .

 I am attaching the codes and sample graphs and  maple query but still getting  error.

here is the codes 

restart;
restart; Digits := 1;
Pr:=0.01:E:=1:a:=0:N:=10:

`Δt`:=0.01:`Δy`:=0.01:

#Discritization scheme

for i from 1 by 1 while i<=N do;  
end:
for j from 0 by 1 while j<=N do;
end:

eq1[i, j] := (U[i, j+1]-U[i, j])/`&Delta;t` = (1/2)*Gr*(theta[i, j+1]+theta[i, j])+(1/2)*Gc*(C[i, j+1]+C[i, j])+(U[i-1, j+1]-2.*U[i, j+1]-2.*U[i, j]+U[i+1, j])/(2.*`&Delta;y`)^2-(1/2)*M*(U[i, j+1]+U[i, j]):

eq2[i, j] := (theta[i, j+1]-theta[i, j])/`&Delta;t` = (1/Pr)*(theta[i-1, j+1]-2*theta[i, j+1]+theta[i+1, j+1]+theta[i-1, j]-2*theta[i,j]+theta[i+1,j])/(2.*`&Delta;y`)^2-E*((1/`&Delta;y`)^2*(U[i+1, j]-U[i, j])^2):

eq3[i, j] := (C[i, j+1]-C[i, j])/`&Delta;t` = (1/Sc)*(C[i-1, j+1]-2*C[i, j+1]+C[i+1, j+1]+C[i-1, j]-2*C[i,j]+C[i+1,j])/(2*`&Delta;y`)^2-(K/2)*(C[i, j+1]-C[i, j]):
end do;  
Error, reserved word `end` unexpected
end do:
Error, reserved word `end` unexpected

# initial conditions
U[i, 0] := 0:
theta[i, 0]:= 0:
C[i, 0] := 0:
 

NULL;
U[0,j]:=exp(a*j*`&Delta;t`):
theta[0,j]:=(j*`&Delta;t`):
C[0,j]:=j*`&Delta;t`:

U[N,j]:=0:
theta[N,j]:=0:
C[N,j]:=0:

sys := ([seq])(seq(eq[i, j], j = 0 .. N), i = 1 .. N):
nops(sys);
vars:=indets(sys):
nn := Matrix(N+1, N+1,(i, j)-> U[i-1, j-1]):
##
p:=proc(kk) local U_res,A;
  U_res:=solve(eval(sys,k=kk),vars);
  A:=eval(nn,U_res);
  plots:-matrixplot(A)
end proc;

plots:-(U,[M],Y=0..4);

 

https://www.mapleprimes.com/questions/225117-How-To-Solve-This-Error-In-Maple

Dear maple user  any one suggest me how to solve  second order coupled differential equation using galerkin finite element method for 8 elements and 10 elements using maple codes

 

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