@tomleslie 

Thanks for the suggestions,

if we change the values of H=0.1, 0.2,0.3  we would get three graphs one after the other, i am attaching the maple codes,  and sample graphs also .

instead of applying individual values of H is it possible to apply all H values in some procedure such that if we run the programme all the graphs should be display one after the other.

my_codes.mw

@acer  thanks , for your prompt reply,  in the above graph i have taken t=Pi/3, if we take t:=0,Pi/4,Pi/3; then how we will  plot for different values of t  

  sample graph

@tomleslie  thanks for your prompt reply , 

can we get the above result like this form

MM1:=(1/(2*lambda))*[C[o]^2*exp(2*lambda*z)+D[o]^2*exp(-2*lambda*z)]+(1/(2))*[C[o]*exp(lambda*z)+D[o]*exp(-lambda*z)]-(r^4/(2))*[C[o]*exp(lambda*z)+D[o]*exp(-lambda*z)]+r^2*(int(A[1](z), z));
 

@acer  thanks a lot can we further apply the same codes to collect C[o] and  C1 terms for  P3:=int(u[1](r,z),r); also 

@acer  I am unable to combine and collect  , please see the codes 

restart:
with(DETools):
H:=proc(ee)
     local temp;
     temp:=sort([op(indets(ee,And(polynom(integer,r),
                                  satisfies(u->degree(u,r)>0))))],
                (a,b)->degree(a,r)>degree(b,r) and length(a)>length(b));
     subs(map(u->u=freeze(u),temp),ee);
end proc:

P[o](z):=C[o]*exp(lambda*z)+D[o]*exp(-lambda*z):
u[o](r,z):=(1-r^2)*diff(P[o](z),z):
v[o](r,z):=(2*r-r^3)*diff(P[o](z),z,z):
     
A1:=P[o](z)*u[o](r,z):
A2:=lambda^2*P[o](z)*(1+r^2)+v[o](r,z):
v[1](r,z):=(gamma1^2/16)*diff(P[o](z),z,z)*(2*r-r^3):
u[1](r,z):=G*lambda*(C1^2*exp(2*lambda*z)-D1^2*exp(-2*lambda*z))*(((1-r^4)/16)-((1-r^6)/72))+(lambda/2)*(C1*exp(lambda*z)-D1*exp(-lambda*z))*(1-r^4)-A1*(1-r^2)/4:
gamma1:=4/lambda:
P1:=gamma1^2*diff(u[o](r,z),z,z)+diff(u[1](r,z),r,r)+(1/r)*diff(u[1](r,z),r)-G*(gamma1^2*u[o](r,z)*diff(u[o](r,z),z)+v[1](r,z)*diff(u[o](r,z),r)):
P2:=int(P1,z):
P[1](z):=simplify(P2):
ans2:=thaw(factor(combine(expand((P[1](z)))))):
 

@acer  thanks a lot 

Carl Love 13652  thanks for excellent results and graphs . can you help to solve this article by applying FEM variational principle only .

Finite Element Analysis of MHD Flow of Micropolar Fluid over a Shrinking Sheet with a Convective Surface Boundary Condition

@Carl Love 13647 Thanks a lot for your excellent post 

please attach maple codes file

@Carl Love  small change in there in Galerkin FEM and Variational FEM as discussed in Gupta paper , the matlab codes in not provided by them. Some one else given matlab code for non linear Differential equation solving by Galerkin FEM ,  if  some body help to convert matlab code to maple then there is a possibility to write the maple codes for system of Non-linear coupled differential equation using variational principle (FEM)as given by Gupta, et al. 

@Rouben Rostamian  

we want to solve the above differential equation using Galerkin Finite Element Method with maple codes and compare the result with exact solution, once if any one can able to write the maple codes using FEM then it can be helpful to solve non-linear coupled differential equation also.

@Carl Love Thanks for responding.methodology.pdf

methodology for FEM by variational principle, please see the methodology in the paper they have used 2X2 matrix in that, if maple code can be written for this problem than it can be written for system of non linear differential coupled equation also as given by Gupta et.al paper . I am also trying to understanding how maple is implemented for FEM and how it can be coded and error can be plotted with graphs.

@Carl Love  thanks for  responding the comment i further trying to understand the equations  Eqs 23-31. further i dug  one more article of gupta el. al there they have given some  methodology . As far as my knowledge i think  they have used matlab to solve those equations .  As   ( Bagh Ali 5 ) particulary want procedure to solve FEM using maple by  variational formula . so it is very helpful that  a procedure can be written on maple so that comparision with analytical , FEM  and error also can be found with the help of graph and numerically.

FEM_paper2.pdfFEM_paper2.pdf

@Carl Love  sorry for reminding again and again , i am also enthusiastic about how to solve non- linear differential equations by  using FEM( maple codes). what ever codes you have applied to solve the non-linear differential equation i have applied the same codes to some other problem still i am getting error . I am attaching the complete paper which i have try to solve , but particularly  ( Bagh Ali 5 ) in his comment ask to solve the non-linear differential  equation  using FEM by maple.

I am also want to interested to see how the nonlinear ODE can be solved using  FEM in maple.

here the maple codes and paper attaching

restart:
Digits:= 15:

ODEs:= [
   #Eq 1/2/3:
   (1+K)*diff(u(y),y$2) + K*diff(N(y),y) + theta(y)-M=0,
   #Eq 10/19/25:
   (1+K/2)*diff(N(y),y$2) -K*(2*N(y)+diff(u(y),y$1))  =  0,
   #Eq 11/20/26:
   diff(theta(y),y$2)=  0
]:                                                            BCs:= [
   #Eq 4:
   u(0) = 0,theta(0) = R, N(0) =0, 
   #Eq 5:
   u(1) = 0,theta(1) = R, N(1) = 0 
   
];                                                               param_names:= [K,M, R];                                       <ODEs[]>;

     BCs:=[uApplyFunction(0)=0,?ApplyFunction(0)=R,NApplyFunction(0)=0,u

       ApplyFunction(1)=0,?ApplyFunction(1)=R,NApplyFunction(1)=0]
                            param_names:=[K,M,R]
                           Matrix(%id = 163720632)
Solve:= 
   subs(
      _P= param_names,
      proc({
         K::realcons:= 5,    
         R::realcons:= 0.5,    
         M::realcons:= 5
      })
         userinfo(1, Solve, param_names=~ _P);
         dsolve(
            eval(
               convert(ODEs, set) union convert(BCs, set),
               param_names=~ _P
            ),
            numeric
         )
      end proc
   )
;
proc({K::realcons := 5, M::realcons := 5, R::realcons := 0.5})  ...  end;
infolevel[Solve]:= 1:

P:= M:
vals:= [0, 5, 10]:
sols:= [seq(Solve(P= v), v= vals)]:
colors:= [red, green, blue, black,orange]:
for F in [u,N](y) do
   print(plots:-display(
      [seq(
         plots:-odeplot(sols[k], [y,F], color= colors[k], legend= [P= vals[k]]),
         k= 1..nops(vals)
      )],
      labeldirections= [horizontal,vertical]
   ))
od:
Solve: [K = 5 M = 0 R = .5]
Solve: [K = 5 M = 5 R = .5]
Solve: [K = 5 M = 10 R = .5]

sols(1):
P:= K:
vals:= [0, 5, 10]:
sols:= [seq(Solve(P= v), v= vals)]:
colors:= [red, green, blue, black,orange]:
for F in [u,N](y) do
   print(plots:-display(
      [seq(
         plots:-odeplot(sols[k], [y,F], color= colors[k], legend= [P= vals[k]]),
         k= 1..nops(vals)
      )],
      labeldirections= [horizontal,vertical]
   ))
od:
paper.pdf