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

1 visualizzazione (ultimi 30 giorni)
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

Risposta accettata

Birdman
Birdman il 20 Mar 2018
Change this line
phi2= (1-phi)+phi((rhos*cps)/(rhof*cpf));
to
phi2= (1-phi)+phi*((rhos*cps)/(rhof*cpf));
  2 Commenti
naygarp
naygarp il 20 Mar 2018
How would I add another for loop if I want to vary the value of 'phi' in the above code and obtain 3 graphs for each value of 'phi'. It would be a nested for loop. say like
qq=[0 0.05 0.1];
for phi=qq

Accedi per commentare.

Più risposte (0)

Categorie

Scopri di più su Mathematics in Help Center e File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by