I am getting an error of 'Subscript indices must either be real positive integers or logicals' while solving a runge kutta problem

3 次查看(过去 30 天)
clf;clc;
%parameters
s=0.01;
t_final=10;
global Re Pr R phi b eps Mn Bi alp rhof cpf kf rhos cps ks
Re= 2.5;
Pr= 6.8;
R= 1;
phi= 0.1;
b=1.0;
eps=0.1;
Mn= 1.0;
Bi= 1.0;
alp= 1;
rhof= 997.1;
cpf= 4179.0;
kf= 0.613;
rhos= 8933.0;
cps= 385.0;
ks= 400;
phi1= (1-phi)^2.5*((1-phi)+phi*(rhos/rhof));
phi2= (1-phi)+phi((rhos*cps)/(rhof*cpf));
k = (ks+2*kf-2*phi*(kf-ks))/(ks+2*kf+phi*(kf-ks));
%initial conditions
t(1)=0;
h(1)=1;
%function handle
f=@(t,h) -(1+b)*h+eps*((1+b^2)*phi1*Re+(1+b)*((1-phi)^2.5)*Mn^2)*h^3/3+...
2*eps*alp*(1-phi)^2.5*(k/(k+Bi*h))*h^2-(2/3)*eps^2*(1-phi)^2.5*alp*Pr*(Re*phi2/k+R)*...
((2*k*Bi^2*h^3*h_t)/(k+Bi*h)^3+(2*(1+b)*Bi*k*h^3)/(k+Bi*h)^3+((1+b)*(3*k+Bi*h)*k*h^2)/(k+Bi*h)^2)*h^2/2
% RK4 loop
for i=1:ceil(t_final/s)
t(i+1)=t(i)+s;
k1= f(t(i) , h(i));
k2= f(t(i)+0.5*s , h(i)+0.5*k1*s);
k3= f(t(i)+0.5*s , h(i)+0.5*k2*s);
k4= f(t(i)+s , h(i)+k3*s);
h(i+1)= h(i)+s/6*(k1+2*k2+2*k3+k4);
end

采纳的回答

Birdman
Birdman 2018-3-20
Change this line
phi2= (1-phi)+phi((rhos*cps)/(rhof*cpf));
to
phi2= (1-phi)+phi*((rhos*cps)/(rhof*cpf));

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by