Advanced Search

Gabriel samaila

35 Reputation

4 Badges

6 years, 141 days
Contact Gabriel samaila

MaplePrimes Activity


These are replies submitted by Gabriel samaila

Re- Use separate dsolve for each Nu...

Gabriel samaila 35
June 11 2021
0 0

@Carl Love 

thank you so much for your response, is now running perfectely. However i notice that there is a problem of convergency. I think that is why the author ishak used shooting method before solving it. I dont know if you have idea of shooting method?

In his paper, the last figure is round (complete circle) and in the present code, if  i extent eta to any thing more than ten, the solution will not converge but in Ishark work, the solution converger. Or do you think of other ways? Please help because i want to extended his work.

new.mw

 
 

Download new.mw

 

Re- Use separate dsolve for each Nu...

Gabriel samaila 35
June 04 2021
0 0

@Carl Love 

Thank you for attending to my problem

attach is the code you requested for, possibly i did not implement it as instructed. please see the code for more details on what i want to do (written in upper case with red colour)

new.mw

 

Re- Use separate dsolve for each Nu...

Gabriel samaila 35
June 01 2021
0 0

@Carl Love 

Thank you for attending to my problem

attach is the code you requested for, possibly i did not implement it as instructed. What actually i am trying to do is to plot f'' against S(episolon) as in Ishak paper. i have already ploted figure 1 (theta against eta) but i could not plot f'' against episolon

ishak.mw          

Link to paper abstract
[article pdf uploaded without copyright permission, removed by moderator]

Re- Use separate dsolve for each Nu...

Gabriel samaila 35
May 24 2021
0 0

@Carl Love 

thank you so much for your response, you alway respond to my request i really thank you. However, i could not run it using maple 16, i may have to look for the recent version and try it.

ExportMatrix...

Gabriel samaila 35
September 06 2020
0 0

@Carl Love 

Thank you for your response. However  I tried to do as you suggested but i couldnt get it right. 

Attach is your code you ever sent me to solve a particular problem. i have modify it to solve Reddy et al paper which i want to extend. Could you please show me how to get the numeric data of those plots?.

thank you in antipation 

reddy.mw
Download reddy.mw

Exporting maple plot into MATLAB...

Gabriel samaila 35
September 04 2020
0 0

@tomleslie 

Thanks for your response I will try and another way I think is to use print instead of plot in maple then I get the .dat then i can import it into MATLAB. 

Is that possible?

homotopy...

Gabriel samaila 35
May 26 2019
0 0

@tomleslie 

you always help me sort every sort of problem on maple, i'm very graful for all your post. i just wish i'm smart like you.

i really appreciate your effort.

homotopy...

Gabriel samaila 35
May 25 2019
0 0

@tomleslie 

thank you so much for your response is very helpful,  i think  the code is not supporting Maple 2016 version, see the error massage. Also is it possible to write boundary condition for u0,w0,f0 and theta0 different from ui,wi,fi and theta i for i=1..N?

thank you
 

> 

  restart;
  PDEtools[declare](f(x),theta(x),u(x),w(x), prime=x):
  N:= 4;

f(x)*`will now be displayed as`*f

  

theta(x)*`will now be displayed as`*theta

  

u(x)*`will now be displayed as`*u

  

w(x)*`will now be displayed as`*w

  

`derivatives with respect to`*x*`of functions of one variable will now be displayed with '`

  

4

 (1)
> 

  F:=      x -> local i;
                add
                ( p^i*f[i](x),
                  i=0..N
                );
  Theta:=  x -> local i;
                add
                ( p^i*theta[i](x),
                  i = 0..N
                );
  U:=      x -> local i;
                add
                ( p^i*u[i](x),
                  i = 0..N
                );
  W:=      x -> local i;
                add
                ( p^i*w[i](x),
                  i = 0 .. N
                );

Error, reserved word `local` unexpected

  
> 

  HPMEq := collect
           ( (1-p)*diff(F(x), x$2)+p*(diff(F(x), x$2)-k1*diff(F(x), x)-k2*F(x)),
              p
           );
  HPMEr := collect
           ( (1-p)*diff(Theta(x),x$2)+p*(diff(Theta(x),x$2)-k11*diff(Theta(x),x)+k12*diff(U(x),x)^2+k13*diff(W(x),x)^2+k14*Theta(x)),
              p
           );
  HPMEs := collect
           ( (1-p)*diff(U(x),x$2)+p*(diff(U(x),x$2)-R*diff(U(x),x)-A-k8*W(x)-k7*U(x)+k5*Theta(x)+k6*F(x)),
              p
           );
  HPMEt := collect
           ( (1-p)*diff(W(x),x$2)+p*(diff(W(x),x$2)-R*diff(W(x),x)+k9*U(x)-k10*W(x)),
              p
           );

(-k1*(diff(F(x), x))-k2*F(x))*p+diff(diff(F(x), x), x)

  

(-k11*(diff(Theta(x), x))+k12*(diff(U(x), x))^2+k13*(diff(W(x), x))^2+k14*Theta(x))*p+diff(diff(Theta(x), x), x)

  

(-R*(diff(U(x), x))-A-k8*W(x)-k7*U(x)+k5*Theta(x)+k6*F(x))*p+diff(diff(U(x), x), x)

  

(-R*(diff(W(x), x))+k9*U(x)-k10*W(x))*p+diff(diff(W(x), x), x)

 (2)
> 

#
# renamed equ[1] to equ1
#
  for i from 0 to N do
      equ1[i] := coeff(HPMEq, p, i) = 0;
  end do;
#
# renamed equa[1] to equ2
#
  for i from 0 to N do
      equ2[i] := coeff(HPMEr, p, i) = 0;
  end do;
#
# renamed equat[1] to equ3
#
  for i from 0 to N do
      equ3[i] := coeff(HPMEs, p, i) = 0;
  end do;
#
# renamed equati[1] to equ4
#
  for i from 0 to N do
      equ4[i] := coeff(HPMEt, p, i) = 0;
  end do

diff(diff(F(x), x), x) = 0

  

-k1*(diff(F(x), x))-k2*F(x) = 0

  

0 = 0

  

0 = 0

  

0 = 0

  

diff(diff(Theta(x), x), x) = 0

  

-k11*(diff(Theta(x), x))+k12*(diff(U(x), x))^2+k13*(diff(W(x), x))^2+k14*Theta(x) = 0

  

0 = 0

  

0 = 0

  

0 = 0

  

diff(diff(U(x), x), x) = 0

  

-R*(diff(U(x), x))-A-k8*W(x)-k7*U(x)+k5*Theta(x)+k6*F(x) = 0

  

0 = 0

  

0 = 0

  

0 = 0

  

Warning, inserted missing semicolon at end of statement

  

diff(diff(W(x), x), x) = 0

  

-R*(diff(W(x), x))+k9*U(x)-k10*W(x) = 0

  

0 = 0

  

0 = 0

  

0 = 0

 (3)
> 

#
# Sequentially solve the system(s)
#
  for i from 0 to N do
      cons:= f[i](-1) = 1, f[i](1) = 1,
             theta[i](-1)=1, theta[i](1)=1,
             u[i](-1)=1, u[i](1)=1,
             w[i](-1)=1, w[i](1)=1;
      assign(dsolve( [ equ1[i], equ2[i], equ3[i], equ4[i], cons ]));
  od:

Error, (in dsolve) found the following equations not depending on the unknowns of the input system: {f[0](-1) = 1, f[0](1) = 1, u[0](-1) = 1, u[0](1) = 1, w[0](-1) = 1, w[0](1) = 1, theta[0](-1) = 1, theta[0](1) = 1}

  
> 

#
# What happens next??
#
# OP can access any of the "individual" functions,
# for example theta[2](x) or w[3](x), just by entering
# these expressions, as in
#
  collect(theta[2](x), x);
  collect(w[3](x), x);
#
# Or one can get the "sum" terms for these functions
# just by entering F(x), Theta(x), U(x) or W(x). These
# expresssions are a bit lengthy
#
# For example, U(x) is shown below
#
  collect(Theta(x), x);

theta[2](x)

  

w[3](x)

  

Theta(x)

 (4)
> 

 


 

Download HOM2.mw

Table...

Gabriel samaila 35
March 07 2019
0 0

@tomleslie 

i said may be because I check it every seems to be same but I’m surprised is not coming up.

table...

Gabriel samaila 35
March 07 2019
0 0

@tomleslie 

thank you so much for all your comment, they are all helpful, i have attach the source of the governing equation probalily i did not input it correctely. The table is not include in the pdf but i found the table in journal book compile by the author. i have tried to attach the journal book but is very heavy it couldnt attach.Sc.docx


 

> 

  restart;
  kernelopts(version);

`Maple 16.00, IBM INTEL NT, Feb 23 2012, Build ID 731223`

 (1)
> 

#
# Define ODEs
#
  ODEs:= {diff(phi(r), r, r)+(1)/(r)*diff(phi(r),r)-Re*Sc*v(r)*(diff(phi(r), r))+Nt*(diff(theta(r), r, r)+(diff(theta(r), r))/r)/Nb = 0,
         diff(v(r), r, r)+(diff(v(r), r))/r-(Ha^2/(1-eta)^2+1/r^2)*v(r)-Re*v(r)*(diff(v(r), r)) = 0,
         diff(theta(r), r, r)+(1)/(r)*diff(theta(r),r)+Ec*Pr*(diff(v(r), r)-v(r)/r)^2-Pr*Re*v(r)*(diff(theta(r), r))+Nb*(diff(theta(r), r))*(diff(phi(r), r))+Nt*(diff(theta(r), r))^2 = 0};

{diff(diff(phi(r), r), r)+(diff(phi(r), r))/r-Re*Sc*v(r)*(diff(phi(r), r))+Nt*(diff(diff(theta(r), r), r)+(diff(theta(r), r))/r)/Nb = 0, diff(diff(v(r), r), r)+(diff(v(r), r))/r-(Ha^2/(1-eta)^2+1/r^2)*v(r)-Re*v(r)*(diff(v(r), r)) = 0, diff(diff(theta(r), r), r)+(diff(theta(r), r))/r+Ec*Pr*(diff(v(r), r)-v(r)/r)^2-Pr*Re*v(r)*(diff(theta(r), r))+Nb*(diff(theta(r), r))*(diff(phi(r), r))+Nt*(diff(theta(r), r))^2 = 0}

 (2)
> 

#
# Define boundary conditions: note that these restrict
# the solution to the range eta < r < 1.0 (or possibly
# 1.0 < r < eta, depending on whether eta>1 or eta<1)
#
# In the parameter list provided by the OP, eta is
# hardwired to 0.5, so valid solution will be in the range
#
#  0.5 <= r <= 1.0
#
  BCs := { phi(1) = 0, phi(eta) = 1, theta(1) = 0,
           theta(eta) = 1, v(1) = 0, v(eta) = 1
         }:
> 

#
# Define values for Ha and Sc which are used
# to generarte output table
#
  HaVals:=[0, 2, 4, 6]:
  ScVals:=[0.5, 1.0, 2.0]:
#
# Solve the ODE system for all pairs of values in
# the lists HaVals, ScVals, and park all solutions
# in a Maple Table()
#
  for j from 1 by 1 to numelems(HaVals) do
      for k from 1 by 1 to numelems(ScVals) do
        #
        # Generate parameter list
        #
          pVals:= [ eta = .5,       Ha = HaVals[j],
                    Sc = ScVals[k], Nt = .1,
                    Nb = .1,        Re = 1,
                    Ec = 0.01,    Pr = 10
                  ];
        #
        # Produce ODE solution for this parameter list
        #
          sol[j,k]:=dsolve( `union`( eval(ODEs, pVals),
                                     eval(BCs, pVals)
                                   ),
                             numeric,
                             output=listprocedure
                           );
      od:
  od:
> 

#
# generate plots for theta(y), phi(r), v(r) for
# all Ha and Sc values.
#
# Note that theta(r) and v(r) don't seem to vary
# as Sc changes value (but phi(r) does)
#
  pfunc:= z-> plot
              ( [ seq
                  ( seq
                    ( eval(z(r), sol[i,j])(r),
                      i=1..numelems(HaVals)
                    ),
                    j=1..numelems(ScVals)
                  ) ],
                  r=0.5..1,
                  title=cat( convert(z, string),
                             "(r) vs r"
                           ),
                  titlefont=[times, bold, 24],
                  legend=[ seq
                           ( seq
                             ( typeset
                               ( "Ha= ", HaVals[i],
                                 ", Sc= ", ScVals[j]
                               ),
                               i=1..numelems(HaVals)
                             ),
                             j=1..numelems(ScVals)
                           )
                         ],
                  legendstyle=[font=[times, bold, 16]]
              ):
  pfunc(theta);
  pfunc(phi);
  pfunc(v);

 

 

 
> 

#
# Define function for the output table. This could
# be any expression involving  phi(y), theta(y), v(y)
# so just for illustration, let's use
#
#   phi(y)*theta(y)/v(r)5
#
  func:=phi(r):
#
# Define the value of 'y' at which the output is
# required. The boundary conditions require that
# 0.5 < r < 1.0. For illustration, let's use
#
#   r=0.75
#
  argr:=0.75:
#
# Print header lines for output table
#
  printf
  ( "%30s\n %8s%12.6f%12.6f%12.6f\n\n",
    "Sc","Ha", seq
               ( ScVals[j],
                 j=1..numelems(ScVals)
               )
  );

  seq
  ( printf
    ( " %8d%12.6f%12.6f%12.6f\n%",
      HaVals[k],
      seq
      ( eval
        ( func, sol[k,j]
        )(argr),
        j=1..numelems(ScVals)
      )
    ),
    k=1..numelems(HaVals)
  );

                            Sc
       Ha    0.500000    1.000000    2.000000

        0    0.146498    0.161720    0.192380
        2    0.207675    0.218275    0.239489
        4    0.284420    0.290105    0.301407
        6    0.325064    0.328319    0.334773

  
> 

 


 

Download new_odeProb6.mw

table...

Gabriel samaila 35
February 25 2019
0 0

@tomleslie 

eventhough is not specified at what value of phi(y) is given but i think is at phi(0). 

 

Attachment...

Gabriel samaila 35
February 20 2019
0 0

@acer 

sorry about the attachemen, i have been trying to get the attached table using same code but i couldnt get it.Sc.docx


 

> 

restart;
kernelopts(version);

`Maple 16.00, IBM INTEL NT, Feb 23 2012, Build ID 731223`

 (1)
> 

ODEs:= {diff(phi(y), y, y)-Re*Sc*v(y)*(diff(phi(y), y))+Nt*(diff(theta(y), y, y)+(diff(theta(y), y))/y)/Nb = 0,
         diff(v(y), y, y)+(diff(v(y), y))/y-(Ha^2/(1-eta)^2+1/y^2)*v(y)-Re*v(y)*(diff(v(y), y)) = 0,
         diff(theta(y), y, y)+Ec*Pr*(diff(v(y), y)-v(y)/y)^2-Pr*Re*v(y)*(diff(theta(y), y))+Nb*(diff(theta(y), y))*(diff(phi(y), y))+Nt*(diff(theta(y), y))^2 = 0};

{diff(diff(phi(y), y), y)-Re*Sc*v(y)*(diff(phi(y), y))+Nt*(diff(diff(theta(y), y), y)+(diff(theta(y), y))/y)/Nb = 0, diff(diff(v(y), y), y)+(diff(v(y), y))/y-(Ha^2/(1-eta)^2+1/y^2)*v(y)-Re*v(y)*(diff(v(y), y)) = 0, diff(diff(theta(y), y), y)+Ec*Pr*(diff(v(y), y)-v(y)/y)^2-Pr*Re*v(y)*(diff(theta(y), y))+Nb*(diff(theta(y), y))*(diff(phi(y), y))+Nt*(diff(theta(y), y))^2 = 0}

 (2)
> 

BCs := {phi(1) = 0, phi(eta) = 1, theta(1) = 0, theta(eta) = 1, v(1) = 0, v(eta) = 1};

{phi(1) = 0, phi(eta) = 1, v(1) = 0, v(eta) = 1, theta(1) = 0, theta(eta) = 1}

 (3)
> 

HaVals:=[1, 2, 5, 10]:
for j from 1 by 1 to numelems(HaVals) do
    pVals := [eta = .5, Ha = HaVals[j], Sc = .8, Nt = .1, Nb = .1, Re = 2, Ec = 0.1e-1, Pr = 10];
    sol[j]:=dsolve( `union`(eval(ODEs, pVals), eval(BCs, pVals)),
                     numeric,
                     output=listprocedure
                  );
od:
displ_opt:= titlefont=[times, bold, 20],
            legend=[seq(typeset("Ha= ", HaVals[j]),j=1..4)],
            legendstyle=[font=[times, bold, 16]],
            color=[red, green, blue, black]:
plot( [seq( eval(theta(y), sol[j])(y),j=1..4)],
      y=0.5..1,
      title="theta(y) vs y",
      displ_opt
    );
plot( [seq( eval(phi(y), sol[j])(y),j=1..4)],
      y=0.5..1,
      title="phi(y) vs y",
      displ_opt
     );
plot( [seq( eval(v(y), sol[j])(y),j=1..4)],
      y=0.5..1,
      title="v(y) vs y",
      displ_opt
   );

 

 

 
> 

 

> 

 


 

Download odeProb4.mw

Maple 2016...

Gabriel samaila 35
February 16 2019
0 0

@acer 

yes I’m using maple 2016 and code is running now

Running...

Gabriel samaila 35
February 16 2019
0 0

@tomleslie 

i said the code has been running in maple 2016 but I had to remove the command I have shown you above before it work. 

Thank you

maple 2016...

Gabriel samaila 35
February 15 2019
0 0

@tomleslie 

thank you for your recent post, i'm using maple 2016. the code can run  perfectely without the command

size=[800,600],

1 2 3 4 Page 1 of 4