## 120 Reputation

2 years, 205 days

## How to obtain multiple solution...

Maple

How to obtain the multiple solution and graph given in the paper.

Stefan Blowing and Slip Effects on Unsteady Nanofluid Transport Past a Shrinking Sheet: Multiple Solutions

https://doi.org/10.1002/htj.21470

Can anyone help to get solutions.

## How to solve 2 ODE equation with Boundar...

Maple

Equations

ODES := (diff(f(eta), `\$`(eta, 4)))/((1-phi1)^2.5*(1-phi2)^2.5*((1-phi2)*(1-phi1+phi1*rhos1/rhosf)+phi2*rhos2/rhosf))+S*(f(eta)*(diff(f(eta), `\$`(eta, 3)))-3*(diff(f(eta), `\$`(eta, 2)))-eta*(diff(f(eta), `\$`(eta, 3)))-(diff(f(eta), eta))*(diff(f(eta), `\$`(eta, 2)))) = 0,

(khnf/kf+(4/3)*R)*(diff(theta(eta), `\$`(eta, 2)))/((1-phi2)*(1-phi1+phi1*rhos1*cp1/(rhosf*cpf))+phi2*rhos2*cp2/(rhosf*cpf))+S*Pr*(f(eta)*(diff(theta(eta), eta))-eta*(diff(theta(eta), eta))-gamma*(eta^2*(diff(theta(eta), `\$`(eta, 2)))-2*eta*f(eta)*(diff(theta(eta), `\$`(eta, 2)))-eta*(diff(f(eta), eta))*(diff(theta(eta), eta))+f(eta)*(diff(f(eta), eta))*(diff(theta(eta), eta))+f(eta)^2*(diff(theta(eta), `\$`(eta, 2))))) = 0

Boundary Conditions

f(0) = 0, ((D^2)(f))(0) = 0, (D(theta))(0) = 0, f(1) = 0, (D(f))(1) = 0, theta(1) = 1

phi1 = .1, phi2 = .1, rhos1 = 2720, rhos2 = 2810, rhosf = 997.1, khnf = 1.083061737, kf = .613, cp1 = 893, cp2 = 960, cpf = 4179, Pr = 6.2, knf = .8154646474., S=0.5,R=0.5, gamma=0.5

I am tried to solve it showing following error

Error, (in dsolve/numeric/bvp) matrix is singular

## How to solve the BVP and draw a 3 D plot...

Maple

eq1 := diff(f(eta), eta, eta, eta, eta)+(2*f(eta)*(diff(f(eta), eta, eta, eta))+2*g(eta)*(diff(g(eta), eta)))*(1-phi)^2.5*(1-phi+phi*rhos/rhof)-sigmanf*M*(diff(f(eta), eta, eta))*(1-phi)^2.5/sigmaf = 0

eq2 := diff(g(eta), eta, eta)-(1-phi)^2.5*(1-phi+phi*rhos/rhof)*(2*(diff(f(eta), eta))*g(eta)-2*(diff(g(eta), eta))*f(eta))-sigmanf*M*g(eta)*(1-phi)^2.5/sigmaf = 0

eq3 := k[nf]*(diff(theta(eta), eta, eta))/(Pr*k[f])+((1-phi+phi*rhos*cps/(rhof*cpf))*2)*f(eta)*(diff(theta(eta), eta))-4*lambda*(1-phi+phi*rhos*cps/(rhof*cpf))*(f(eta)^2*(diff(theta(eta), eta, eta))+f(eta)*(diff(f(eta), eta))*(diff(theta(eta), eta)))+sigmanf*M*Ec*((diff(f(eta), eta))^2+g(eta)^2)/sigmaf = 0

eq4 := (1-phi)^2.5*(diff(chi(eta), `\$`(eta, 2)))+2*Sc*f(eta)*(diff(chi(eta), eta))-sigma*Sc*(1+delta*theta(eta))^n*exp(-E/(1+delta*theta(eta)))*chi(eta) = 0

Boundary Conditions

f(0) = 0, (D(f))(0) = A1+gamma1*((D@@2)(f))(0), f(10) = 0, (D(f))(10) = 0, g(0) = 1+gamma2*(D(g))(0), g(10) = 0, theta(0) = 1+gamma3*(D(theta))(0), theta(10) = 0, chi(0) = 1, chi(10) = 0

Parameters

lambda = 0.1e-1, sigma = .1, Ec = .2, E = .1, M = 5, delta = .1, n = .1, Sc = 3, A1 = .5, gamma1 = .5, gamma2 = .5, gamma3 = .5, Pr = 6.2, phi = 0.1e-1, rhos = 5.06*10^3, rhof = 997, cps = 397.21, cpf = 4179, k[nf] = .6358521729, k[f] = .613, sigmanf = 0.5654049962e-5, sigmaf = 5.5*10^(-6)

## Boundary Value Problem...

Maple

I am tried to solve the following problem. here is the code and boundary conditions as well as parameters used in the problem. Please help me to get the numerical solution and getting plots between Cu and eta as well as D(f)(eta) vs eta.

restart;
Digits := trunc(evalhf(Digits));
15
ODEs := [diff(f(eta), `\$`(eta, 3))+A^2+f(eta)*(diff(f(eta), `\$`(eta, 2)))-(diff(f(eta), eta))^2+beta*((diff(g(eta), eta))^2-g(eta)*(diff(g(eta), `\$`(eta, 2)))-1), lambda*(diff(g(eta), `\$`(eta, 3)))+(diff(g(eta), `\$`(eta, 2)))*f(eta)-g(eta)*(diff(f(eta), `\$`(eta, 2)))];
`<,>`(ODEs[]);
Vector[column](%id = 18446744073898822582)
LB, UB := 0, 1;
BCs := [`~`[`=`](([D(f), f, g, (D@@2)(g)])(LB), [1+B1*((D@@2)(f))(0), 0, 0, 0])[], `~`[`=`](([D(f), D(g)])(UB), [A, 0])[]];
[D(f)(0) = 1 + B1 @@(D, 2)(f)(0), f(0) = 0, g(0) = 0,

@@(D, 2)(g)(0) = 0, D(f)(1) = A, D(g)(1) = 0]
Params := Record(A = .9, B1 = .5, beta = .5, lambda = .5);
NBVs := [-((D@@2)(f))(1) = `C*__f`];
Cf := `C*__f`;
Solve := module () local nbvs_rhs, Sol, ModuleApply, AccumData, ModuleLoad; export SavedData, Pos, Init;  nbvs_rhs := `~`[rhs](:-NBVs); ModuleApply := subs(_Sys = {:-BCs[], :-NBVs[], :-ODEs[]}, proc ({ A::realcons := Params:-A, B1::realcons := Params:-B1, beta::realcons := Params:-beta, lambda::realcons := Params:-lambda }) Sol := dsolve(_Sys, _rest, numeric); AccumData(Sol, {_options}); Sol end proc); AccumData := proc (Sol::{Matrix, procedure, list({name, function} = procedure)}, params::(set(name = realcons))) local n, nbvs; if Sol::Matrix 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) 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;
colseq := [red, green, blue, brown];
Pc := B1 = .5, A = .1, beta = .5;
Ps := [B1 = .5, A = .1, beta = .5];
Pv := [lambda = [.5, 1, 1.5, 2]];
for i to nops(Ps) do plots:-display([seq(plots:-odeplot(Solve(lhs(Pv[i]) = rhs(Pv[i])[j], Ps[i][], Pc), [eta, theta(eta)], '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/process_input) invalid argument: (B1 = .5)[]