Question: How do I change variables

I want to change my variables. The coding below shows the influence of the Nb and Nt on the local Nusselt number (-theta') against the convection Bi.

But now I want it to be  the influence of the Bi on the local Nusselt number (-theta') against the Nb and Nt.

Can anyone help me?

NULL

with(plots)

DE1 := diff(f(eta), `$`(eta, 3))+f(eta)*(diff(f(eta), `$`(eta, 2)))-(diff(f(eta), eta))^2 = 0

diff(diff(diff(f(eta), eta), eta), eta)+f(eta)*(diff(diff(f(eta), eta), eta))-(diff(f(eta), eta))^2 = 0

(1)

DE2 := diff(theta(eta), `$`(eta, 2))+Pr*f(eta)*(diff(theta(eta), eta))+Pr*Nb*(diff(theta(eta), eta))*(diff(phi(eta), eta))+Pr*Nt*(diff(theta(eta), eta))^2 = 0

diff(diff(theta(eta), eta), eta)+Pr*f(eta)*(diff(theta(eta), eta))+Pr*Nb*(diff(theta(eta), eta))*(diff(phi(eta), eta))+Pr*Nt*(diff(theta(eta), eta))^2 = 0

(2)

DE3 := diff(phi(eta), `$`(eta, 2))+Le*(diff(phi(eta), eta))+Nt*(diff(theta(eta), `$`(eta, 2)))/Nb

diff(diff(phi(eta), eta), eta)+Le*(diff(phi(eta), eta))+Nt*(diff(diff(theta(eta), eta), eta))/Nb

(3)

BC1 := f(0) = 0, (D(f))(0) = 1, (D(f))(10) = 0

f(0) = 0, (D(f))(0) = 1, (D(f))(10) = 0

(4)

BC2 := theta(10) = 0, (D(theta))(0) = -Bi*(1-theta(0))

theta(10) = 0, (D(theta))(0) = -Bi*(1-theta(0))

(5)

BC3 := phi(0) = 1, phi(10) = 0

phi(0) = 1, phi(10) = 0

(6)

``

  getRes:= proc(a, x);
                local sol;
                if   type(x, numeric)
                then sol:= dsolve
                           ( eval
                             ( {BC1, BC2, BC3, DE1, DE2, DE3},
                               [Nt = a, Nb = a, Pr = 2, Le = 5, Bi = x]
                             ),
                             numeric,
                             output = listprocedure,
                             abserr = 0.0001,
                             maxmesh = 1024,
                             initmesh = 512
                           );
                else return 'procname(a,x)';
                fi;
                return -eval(diff(theta(eta), eta), sol)(0)
           end proc:

  cl:=[red, green, blue, yellow, black]:
  L:=[0.2, 0.4, 0.6,0.8,1.0]:
  display
  ( [ seq
      ( plot
        ( getRes( L[k], x),
          x=0..5,color = colorList[k], legend = ["Nb=Nt" = L[k]], legendstyle = [location = top],
          color=cl[k]
        ),
        k=1..5
      )
    ]
  )
               

 

 

 

Download odeProb2.mw

Please Wait...