Question: Error, (in dsolve/numeric/bvp/convertsys) unable to convert to an explicit first-order system

restart:

Digits:= trunc(evalhf(Digits)); #generally a very efficient setting

15

(1)

Setup of BVP system:

#ordinary differential equations:
ODEs:= [
   #Eq 1:
   A1*(diff(f(x), x, x, x))/(A2*phi)-(diff(f(x), x))^2-M^2*(f(x))+f(x)*(diff(f(x), x, x)),

   #Eq 2:
   A4*Pr*phi*(diff(Theta(x), x, x))/A3+f(x)*(diff(Theta(x), x))+Q*Theta(x)
   
   #All these ODEs are implicitly equated to 0.
]:

<ODEs[]>; #Display the ODEs.

Vector(2, {(1) = A1*(diff(diff(diff(f(x), x), x), x))/(A2*phi)-(diff(f(x), x))^2-M^2*f(x)+f(x)*(diff(diff(f(x), x), x)), (2) = A4*Pr*phi*(diff(diff(Theta(x), x), x))/A3+f(x)*(diff(Theta(x), x))+Q*Theta(x)})

(2)

Params := Record(fw = .2, M = .5, Q = .5, Pr = 6.2, phi = 0.5e-1, rf = 997.1, kf = .613, cpf = 4179, btf = 0.3e-4, p1 = 0.1e-1, p2 = 0.5e-1, p3 = 0.5e-1, rs1 = 5100, ks1 = 3007.4, cps1 = 410, bs1 = 0.2e-3, rs2 = 2200, ks2 = 5000, cps2 = 790, bs2 = 0.5e-3, rs3 = 3970, ks3 = 40, cps3 = 765, bs3 = 0.4e-3, A1 = B1*p1+B2*p2+B3*p3, B1 = 1+2.5*phi+6.2*phi^2, B2 = 1+13.5*phi+904.4*phi^2, B3 = 1+37.1*phi+612.6*phi^2, B4 = (ks1+2*kf-2*phi*(kf-ks1))/(ks1+2*kf+phi*(kf-ks1)), B5 = (ks2+3.9*kf-3.9*phi*(kf-ks2))/(ks2+3.9*kf+phi*(kf-ks2)), B6 = (ks3+4.7*kf-4.7*phi*(kf-ks3))/(ks3+4.7*kf+phi*(kf-ks3)), A2 = 1-p1-p2-p3+p1*rs1/rf+p2*rs2/rf+p3*rs3/rf, A3 = B4*p1+B5*p2+B6*p3, A4 = 1-p1-p2-p3+p1*rs1*cps1/(rf*cpf)+p2*rs2*cps2/(rf*cpf)+p3*rs3*cps3/(rf*cpf))

Record(fw = .2, M = .5, Q = .5, Pr = 6.2, phi = 0.5e-1, rf = 997.1, kf = .613, cpf = 4179, btf = 0.3e-4, p1 = 0.1e-1, p2 = 0.5e-1, p3 = 0.5e-1, rs1 = 5100, ks1 = 3007.4, cps1 = 410, bs1 = 0.2e-3, rs2 = 2200, ks2 = 5000, cps2 = 790, bs2 = 0.5e-3, rs3 = 3970, ks3 = 40, cps3 = 765, bs3 = 0.4e-3, A1 = B1*p1+B2*p2+B3*p3, B1 = 1+2.5*phi+6.2*phi^2, B2 = 1+13.5*phi+904.4*phi^2, B3 = 1+37.1*phi+612.6*phi^2, B4 = (ks1+2*kf-2*phi*(kf-ks1))/(ks1+2*kf+phi*(kf-ks1)), B5 = (ks2+3.9*kf-3.9*phi*(kf-ks2))/(ks2+3.9*kf+phi*(kf-ks2)), B6 = (ks3+4.7*kf-4.7*phi*(kf-ks3))/(ks3+4.7*kf+phi*(kf-ks3)), A2 = 1-p1-p2-p3+p1*rs1/rf+p2*rs2/rf+p3*rs3/rf, A3 = B4*p1+B5*p2+B6*p3, A4 = 1-p1-p2-p3+p1*rs1*cps1/(rf*cpf)+p2*rs2*cps2/(rf*cpf)+p3*rs3*cps3/(rf*cpf))

(3)

LB, UB := 0, 1; BCs := [`~`[`=`](([f(x), diff(f(x), x), Theta])(LB), [fw, 1, 1])[], `~`[`=`](([diff(f(x), x), Theta])(UB), [0, 0])[]]

[(f(x))(0) = fw, (diff(f(x), x))(0) = 1, Theta(0) = 1, (diff(f(x), x))(1) = 0, Theta(1) = 0]

(4)

NBVs := [A1*(diff(f(x), x, x))(0) = C*`*f`, -A4*(diff(Theta(x), x))(0) = `Nu*`]; Nu := `Nu*`; Cf := `C*__f`; x0 := Array([LB])

NULL

Solve := module () local nbvs_rhs, Sol, Dsolve, ModuleApply, AccumData, ModuleLoad; export SavedData, Pos, Init;  nbvs_rhs := `~`[rhs](:-NBVs); Dsolve := proc (Sys, Params::(set(name = realcons))) option remember; Sol := dsolve(Sys, _rest, 'numeric'); AccumData(Params); eval(Sol) end proc; ModuleApply := subs(_Sys = {:-BCs[], :-NBVs[], :-ODEs[]}, proc ({ fw::realcons := Params:-fw, Pr::realcons := Params:-Pr, M::realcons := Params:-M, Q::realcons := Params:-Q, phi::realcons := Params:-phi }) Dsolve(_Sys, {_options}, {_rest}[]) end proc); AccumData := proc (params::(set(name = realcons))) local n, nbvs; if Sol::rtable then nbvs := seq(n = Sol[2, 1][1, Pos(n)], n = nbvs_rhs) else nbvs := `~`[`=`](nbvs_rhs, eval(nbvs_rhs, Sol(:-LB)))[] end if; SavedData[params] := Record[packed](params[], nbvs); return  end proc; ModuleLoad := eval(Init); Init := proc () Pos := proc (n::name) local p; option remember; member(n, Sol[1, 1], 'p'); p end proc; SavedData := table(); return  end proc; ModuleLoad() end module

NULL

colseq := [red, green, blue, brown]

#parameter values that remain fixed for the entire set of plots:
Pc:= phi=0.05:
 

#parameter values that remain fixed with each of the four plots::
Ps:= [
   [fw=0.2, Pr=6.2, M=0.5],
   [fw=0.2, Q=0.3, M=0.5],
   [fw=0.2, Pr=6.2, Q=0.3],
   [Q=0.3, Pr=6.2, M=0.5]
]:

#parameter value for each curve
Pv:= [
   Q=[0.2, 0.4, 0.6, 0.8],
   Pr=[0.7, 1.4, 2.1, 2.8],
   M=[0.6, 1.2, 1.8, 2.4],
   fw=[1, 2, 3, 4]
]:
      

for i to nops(Ps) do
   plots:-display(
      [seq(
         plots:-odeplot(
            Solve(lhs(Pv[i])= rhs(Pv[i])[j], Ps[i][], Pc),
            [x, Theta(x)], 'color'= colseq[j], 'legend'= [lhs(Pv[i])= rhs(Pv[i])[j]]
         ),
         j= 1..nops(rhs(Pv[i]))
      )],
      'axes'= 'boxed', 'gridlines'= false,
      'labelfont'= ['TIMES', 'BOLDOBLIQUE', 16],
      'caption'= nprintf(
         cat("\n%a = %4.2f, "$nops(Ps[i])-1, "%a = %4.2f\n\n"), (lhs,rhs)~(Ps[i])[]
      ),
      'captionfont'= ['TIMES', 16]
   )
od;

Error, (in dsolve/numeric/bvp/convertsys) unable to convert to an explicit first-order system

 

for i to nops(Ps) do plots:-display([seq(plots:-odeplot(Solve(lhs(Pv[i]) = rhs(Pv[i])[j], Ps[i][], Pc), [x, D(f(x))], 'color' = colseq[j], 'legend' = [lhs(Pv[i]) = rhs(Pv[i])[j]]), j = 1 .. nops(rhs(Pv[i])))], 'axes' = 'boxed', 'gridlines' = false, 'labelfont' = ['TIMES', 'BOLDOBLIQUE', 16], 'caption' = nprintf(cat(`$`("\n%a = %4.2f, ", nops(Ps[i])-1), "%a = %4.2f\n\n"), `~`[lhs, rhs](Ps[i])[]), 'captionfont' = ['TIMES', 16]) end do

Error, (in dsolve/numeric/bvp/convertsys) unable to convert to an explicit first-order system

 

ParamPlot2d := proc (Y::{`module`, procedure}, X::(name = range(realcons)), FP::(list(name = realcons)), { dsolveopts::(list({name, name = anything})) := [] }) plot(proc (x) options operator, arrow; Y(Solve(lhs(X) = x, FP[], 'abserr' = 0.5e-4, 'interpolant' = false, 'output' = x0, dsolveopts[])) end proc, rhs(X), 'numpoints' = 25, 'axes' = 'boxed', 'gridlines' = false, 'labelfont' = ['TIMES', 'BOLDOBLIQUE', 16], 'caption' = nprintf(cat(`$`("%a = %4.2f, ", nops(FP)-1), "%a = %4.2f"), `~`[lhs, rhs](FP)[]), 'captionfont' = ['TIMES', 16], _rest) end proc

#procedure that extracts Nusselt number from dsolve solution:
GetNu:= (Sol::Matrix)-> Sol[2,1][1, Solve:-Pos(:-Nu)]:

Q:= [0.2, 0.4, 0.6]:
plots:-display(
   [seq(
      ParamPlot2d(
         GetNu, fw= 1..4, [M= 0.5],
         'dsolveopts'= [Q= Q[k], Pr=6.2,  phi=0.05],
         'legend'= [Q= Q[k]], 'color'= colseq[k], 'labels'= [fw, Nu]
      ),
      k= 1..nops(Q)
   )]
);

Error, invalid input: ParamPlot2d expects value for keyword parameter dsolveopts to be of type list({name, name = anything}), but received [[.2, .4, .6] = .2, Pr = 6.2, phi = 0.5e-1]

 

NULL

Download surface_dinesh_paper.mw  please help me to solve the problem

Please Wait...