Simulating a Mathematical Model

Question

Michele Notarantonio 2020-12-3

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/676763-simulating-a-mathematical-model

评论： Star Strider 2020-12-28

My problem is this one: I have some fixed parameters taken by literature and I want simulate my model with more or less 20 equations (11 are differential other algebraic) and I have no ideas on how I can write a code for simulating this. I thought to use ode15s or od45 but I cant find a solution.

I tryed with this code but it seems that doesnt work. Both if I use ode45 or ode15s in some point of the solution the result has no physycal meaning. In particular, I'm intersted to get results about q_out and S3G, S4G and S5G because after I will take the result of the simulation and I will compare it to my experimental data. If the results is really too different from experimental data I proceed to estimate some parameters. As you can watch there are a huge quantities of parameters.

function dy = modello(t,var)
%Parameters
%      Articolo Zaher
Y1 = 1 ; %g/gbiomassa
Y2 = 0.535 ; %g/gbiomassa
Y3 = 1 ; %g/gbiomassa
Y4 = 0.0089 ; %g/gbiomassa
Y5 = 0.0152 ; %g/gbiomassa
Y6 = 0.0155 ; %g/gbiomassa
Y7 = 0.285 ; %g/gbiomassa
Y8 = 1 ; %g/gbiomassa
Y9 = 0.97 ; %g/gbiomassa
Y10 = 0.99 ; %g/gbiomassa
mu_max1 = 2.11 ; %d^-1
mu_max2 = 0.26 ; %d^-1
mu_max3 = 25.97 ; %d^-1
ks1 = 4.66 ; % g/L
ks2 = 0.26 ; % g/L
ki2 = 5.7 ; % g/L
ks3 = 0.75 ; % g/L 11.2*10^-6
kla1 = 0.294*24 ; % d^-1
kla2 = 6.6*24 ; % d^-1
kla3 = 0.261*24 ; % d^-1
%Generals
H_co2 = 7.48 ; %g/atm*L
H_h2 = 0.0016 ; %g/atm*L
H_ch4 = 0.016 ; %g/atm*L
Vg = 9.3 ; % L
Vl = 52.7 ; % L 
T = 328 ; % K
P = 1 ; % atm
pw = 0.155 ; % atm
R = 0.0821 ; % mol*L/atm*K
qh2_in = 1.2*24 ; %l/h
S1_in = 41.619 ; % g/L
S2_in = 0 ; % g/L
S3G_in = 0; % g/L
S5G_in = 0; % g/L
q_l = 0.448*24 ; % l/h
PM_co2 = 44 ; % g/mol
PM_h2 = 2 ; % g/mol
PM_ch4 = 16 ; % g/mol
S4G_in = 0.074 ; % g/L
D = 1/30 ; % d^-1
%Variables
X1=var(1) ;
X2=var(2) ;
X3=var(3) ;
S1=var(4) ;
S2=var(5) ;
S3=var(6) ;
S4=var(7) ;
S5=var(8) ;
S3G=var(9) ;
S4G=var(10) ;
S5G=var(11) ;
%Equation
mu1=mu_max1.*S1./(S1+ks1);
mu2=mu_max2.*S2./(S2+ks2+(S2^2)/ki2);
mu3=mu_max3.*S4./(S4+ks3);
S3star=(S3G./(S3G+S4G+S5G)).*P*H_co2 ;
S4star=(S4G./(S3G+S4G+S5G)).*P*H_h2 ;
S5star=(1-(S3G./(S3G+S4G+S5G))-(S4G./(S3G+S4G+S5G))).*P*H_ch4 ;
q_out=R*T/(P) *Vl*(((kla1/PM_co2)*(S3-S3star))-((kla2/PM_h2)*(S4-S4star))+(kla3/PM_ch4)*(S5-S5star)); %Portata umida
q_in=qh2_in;
%Differential equation
dx=zeros(11,1) ; %pre allocazione
dx(1,1)=(mu1-D-0.05*mu_max1).*X1 ; 
dx(2,1)=(mu2-D-0.027*mu_max2).*X2 ;
dx(3,1)=(mu3-D-0.05*mu_max3).*X3 ;
dx(4,1)=(D).*(S1_in - S1)-Y1.*mu1.*X1 ;
dx(5,1)=-(D).*S2 + Y2.*mu1.*X1-Y3.*mu2.*X2 ;
dx(6,1)=+Y4.*mu1.*X1+Y5.*mu2.*X2-Y6.*mu3.*X3-kla1.*(S3-S3star);
dx(7,1)=+Y7.*mu1.*X1-Y8.*mu3.*X3+kla2.*(S4star-S4);
dx(8,1)=+Y9.*mu2.*X2+Y10.*mu3.*X3-kla3.*(S5-S5star);
dx(9,1)=q_in.*S3G_in./Vg -q_out.*S3G./Vg +(Vl/Vg).*kla1.*(S3-S3star);
dx(10,1)=q_in.*S4G_in./Vg -q_out.*S4G./Vg +(Vl/Vg).*kla2.*(S4-S4star);
dx(11,1)=q_in.*S5G_in./Vg -q_out.*S5G./Vg +(Vl/Vg).*kla3.*(S5-S5star);
end
clc ; clear all; close all
range=[1:1:40];
IC=[2.06;0.68;0.01;41.619;0;0;0;0;0;0.0740;0];
[time, sol]=ode45(@Prova1,range,IC);

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Star Strider 2020-12-28

There may be too many parameters, depending on your data. You will have to determine that.

There are ways to fit your differential equation to your data, one example being: Parameter Estimation for a System of Differential Equations

Nonlilnear parameter estimation problems are extremely sensitive to the initial parameter estimates. One way to get around this problem is to use the ga (genetic algorithm) function to search the entire parameter space for the best set. Even then, it may take several attempts before ga converges on an acceptable set of parameters.

请先登录，再进行评论。

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