solving pohlhausen equation,how I can write true function?

2019 6 6
1 个回答
更新时间：2019 6 7
6 次查看（30 天）

0 个投票

I write this function for solve pohlhausen equation.(T''+Pr/2 F T'=0 ) That F is output of blasius equation( 2F'''+FF''=0).but my code has error.how I can write true function?
function [Tetadot,fdot]=pohlhausen(eta,Teta,f);
fdot=blasius(eta,f);
Tetadot(1)=Teta(2);
Tetadot(2)=-(0.7/2)*f(1)*Teta(2);
Tetadot=Tetadot';
function fdot=blasius(eta,f);
fdot(1)=f(2);
fdot(2)=f(3);
fdot(3)=-(1/2)*f(1)*f(3);
fdot=fdot';
end
end

2 个评论
显示 无 隐藏 无

mahdi etminan
mahdi etminan 2019-6-6
I solve the blasius equation numerically.
clear all;
clc;
p1=0.01;
p2=1;
Err=0.00001;
i=1;
err=1;
while err>Err
if i==1
[eta, f]=ode45('blasius1',[0:0.2: 10],[0 0 p1]);
m1=f(end,2);
[eta,f]=ode45('blasius1',[0:0.2: 10],[0 0 p2]);
m2=f(end,2);
else
p2=p1+(p2-p1)*(1-m1)/(m2-m1);
[eta,f]=ode45('blasius1',[0:0.2: 10],[0 0 p2]);
m2=f(end,2);
err=abs(1-m2);
end
i=i+1;
end
that the function of blasius is:
function fdot=blasius1(eta,f);
fdot(1)=f(2);
fdot(2)=f(3);
fdot(3)=-(1/2)*f(1)*f(3);
fdot=fdot';
I do not know how replace F data in pohlhausen's function to earn Teta numerically.
n1=0.01;
n2=1;
Er=0.00001;
j=1;
er=1;
while er>Er
if j==1
[eta,Teta]=ode45('pohlhausen',[0:0.2: 10],[0 n1]);
h1=Teta(end,1);
[eta,Teta]=ode45('pohlhausen',[0:0.2: 10],[0 n2]);
h2=Teta(end,1);
else
n2=n1+(n2-n1)*(1-h1)/(h2-h1);
[eta,Teta]=ode45('pohlhausen',[0:0.2: 10],[0 n2]);
h2=Teta(end,1);
er=abs(1-h2);
end
j=j+1;
pohlhausen function that I write is:
function [Tetadot,fdot]=pohlhausen(eta,Teta,f);
fdot=bla(eta,f);
Tetadot(1)=Teta(2);
Tetadot(2)=-(0.7/2)*f(1)*Teta(2);
Tetadot=Tetadot';
function fdot=bla(eta,f);
fdot(1)=f(2);
fdot(2)=f(3);
fdot(3)=-(1/2)*f(1)*f(3);
fdot=fdot';
end
end
end
Torsten
Torsten 2019-6-6
You'll have to solve both equations (Blasius and Pohlhausen) together, not one after the other.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Star Strider
Star Strider 2019-6-6

0 个投票

If you have the Symbolic Math Toolbox, let it do the programming for you:
% % 2F'''+FF''=0
% % T''+Pr/2 F T'=0
syms F(t) Pr T(t) t Y
Blasius = 2*diff(F,3)+F*diff(F,2) == 0;
Polhausen = diff(T,2) + Pr/2 * F * diff(T) == 0;
[VF,Sbs] = odeToVectorField(Blasius, Polhausen) % Vector Field Representation & Substitutions
polhfcn = matlabFunction(VF, 'Vars',{t,Y,Pr}) % Anonymous Function
varcell = sym2cell(Sbs); % Cell Array Of Substitutions
varstr = sprintfc('%s',[varcell{:}]); % Cell Array Of Substitution Strings (For ‘legend’ Or Other Uses), Can Also Use The ‘compose’ Function
producing:
polhfcn = @(t,Y,Pr) [Y(2);Pr.*Y(2).*Y(3).*(-1.0./2.0);Y(4);Y(5);Y(3).*Y(5).*(-1.0./2.0)];
and other information.
You would define a value for ‘Pr’ in your workspace, then use the function (in ode45, for example, although I have no idea what solver is most appropriate here):
[t,y] = ode45(@(t,Y)polhfcn(t,Y,Pr), tspan, ics);
Check to be sure I entered the equations correctly. I believe I did, although it is always appropriate to check.

类别

在 帮助中心 和 File Exchange 中查找有关 Data Type Identification 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by