can it be possible to use odeVectorField

1 次查看(过去 30 天)
syms f(x) g(x) h(x) Pr Nb Nt Le
[V] = odeToVectorField( diff(f,3) - diff(f)^2 + f * diff(f,2) == 0,diff(g,2) + Pr * ( f * diff(g) + Nb * diff(g)* diff(h)+ Nt * diff(g)^2) == 0,...
diff(h,2) + Le * f * diff(h) + ( Nt/Nb ) * Pr * ( f * diff(g) + Nb * diff(g)* diff(h)+ Nt * diff(g)^2)== 0)
%% I want to incorporate the expression of [V] into below code, but I am getting wrong expression (not as dydx below)
%%% can it be used in bvp 4c code and in which style
function main
Pr=10; Le=10; Nt=.1;Nb=.1;
xa=0;xb=6;solinit=bvpinit(linspace(xa,xb,1000),[0 1 1 1 0 0 0]);options=[];
sol=bvp4c(@ode,@bc,solinit,options);xint=linspace(xa,xb,100);S=deval(sol,xint);
function res=bc(ya,yb)
res=[ya(1); ya(2)-1; ya(4)-1; ya(6)-1; yb(2); yb(4); yb(6)];
end
function dydx=ode(x,y)
dydx=[y(2); y(3); y(2)^2 - y(3)*y(1); y(5); - Pr*(y(1)*y(5) + Nb*y(5)*y(7) + Nt*y(5)^2);
y(7); - Le*y(7)*y(1) + (Nt/Nb)*Pr*(y(1)*y(5) + Nb*y(5)*y(7) + Nt*y(5)^2)];
end
end
  3 个评论
MINATI PATRA
MINATI PATRA 2020-5-23
Dear Darova that M gives wrong code as you can compare and another question is how to include in bvp 4c code
MINATI
MINATI 2020-5-24
Please anybody invest some time here to help me.

请先登录,再进行评论。

回答(1 个)

darova
darova 2020-5-24
Try this
clc,clear
syms f(x) g(x) h(x) x Y
Pr=10;
Le=10;
Nt=.1;
Nb=.1;
[V,y1] = odeToVectorField( diff(f,3) - diff(f)^2 + f * diff(f,2) == 0,...
diff(g,2) + Pr * ( f * diff(g) + Nb * diff(g)* diff(h)+ Nt * diff(g)^2) == 0,...
diff(h,2) + Le * f * diff(h) + ( Nt/Nb ) * Pr * ( f * diff(g) + Nb * diff(g)* diff(h)+ Nt * diff(g)^2)== 0)
F = matlabFunction(V,'vars', [x Y]);
bc = @(ya,yb) [ya(1) % g
ya(2)-1 % Dg
ya(4)-1 % f
ya(6)-1 % Df
yb(2) % D2f
yb(4) % h
yb(6)]; % dh
xa=0;
xb=6;
solinit = bvpinit(linspace(xa,xb,1000),[0 1 1 1 0 0 0]);
options = [];
sol = bvp4c(F,bc,solinit,options);
plot(sol.x,sol.y)
Change boundary conditions if needed

类别

Help CenterFile Exchange 中查找有关 Calculus 的更多信息

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by