can it be possible to use odeVectorField

1 回表示 (過去 30 日間)
MINATI
MINATI 2020 年 5 月 22 日
コメント済み: MINATI 2020 年 5 月 25 日
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
  1 件のコメント
MINATI
MINATI 2020 年 5 月 25 日
No dear Darova, its not matching. You can compare also

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeAnimation についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by