[scilab-Users] No reaction
Adrien Vogt-Schilb
vogt at centre-cired.fr
Wed Oct 5 14:43:50 CEST 2011
Hi
Try to isalote your problem
if i understood well, the following code
// initial values
x0=[0.1, 1e-5, 0, 15, 1.16, 100,0,297,0.5]';
t=0:0.005:400;
y=ode(x0, 0, t, f);
returns y such that sum(y(6:9,:)>x0) == 0 ?
if this is true, we do not need the plots to solve the problem
can you check that ?
I believe the f function is erroneous.
It seems that dx_1 should be equal to dx(1) at each time step, and that
HION should be equal to x(2) at each time step, etc.
in other terms, some of your phisical variables seem to be represented
by to variables (i am guessing HION=[H+] and x(2)=[H+] also) but scilab
does not have any chance to know that.
if my guess is right, you have to rewrite the f function in a way that
eliminates all references to HION, dx_1, dx_6 and so on
On 05/10/2011 14:27, Jaundre Venter wrote:
> Hi all
>
> Can someone please explain to me the following:
>
> I am busy with a project of simulation the production of penicillin in
> a bio reactor. Now i have 9 ODE's which i want to simulate.
>
> now for some reason the last three graphs i am getting doesn't show
> any response what so ever. i am using the following code.
>
> dx(1)=(((mu)*(x(1)))-(((x(1))/(x(6)))*((dx_6)))*(CO)*(HION)),
> //biomass concentration X
> dx(2)=((z*(((mu)*(x(1)))-(((F)*(x(1)))/(x(6)))))+(QQ)), //hydrogen ion
> concentration H+
> dx(3)=((((mupp)*(x(1)))-((K)*(x(3)))-((x(3))/(x(6)))*(dx_6))*(HION)),
> //Penicilin concentration P
> dx(4)=((-((mu)/(Yxs))*(x(1)))-(((mupp)/(Yps))*(x(1)))-((mx)*(x(1)))+((Fsf)/(x(6)))-((x(4)/(x(6)))*(dx_6))),
> //Substrate concentration S
> dx(5)=(-(((mu)/(Yxo))*(x(1)))-(((mupp)/(Ypo))*(x(1)))-(((mo))*(x(1)))+((kla)*(cll-(x(5))))-(((x(5))/(x(6)))*(dx_6))),
> //dissolved oxygen
> dx(6)=((F+Fab+Floss)*(HION)), // culture Volume V
> dx(7)=(((rq1)*(dx_1)*(x(6)))+(rq2)*(x(1))*(x(6))), //Heat generation Qrxn
> dx(8)=((((F)/(sf))*(Tf-(x(8))))+(1/((x(6))*(pcp)))*(QT)), //
> Temperature T
> dx(9)=(((a1)*(dx_1))+((a2)*(x(1)))+(a3)), // CO2 evolution, CO2
> endfunction
>
> now when i ask for plotting the graphs i am using the following.:
>
> // initial values
> x0=[0.1, 1e-5, 0, 15, 1.16, 100,0,297,0.5]';
> t=0:0.005:400;
> y=ode(x0, 0, t, f);
>
> // the plots of each variable
> da.title.text="BIOMASS CONCENTRATION"
> da.x_label.text="Time, hours";
> da.y_label.text="X,g/l ";
> scf(1);clf; //Opens and clears figure 1
> plot(t,y(1,:))
>
> da.title.text="HYDROGEN ION H+ CONCENTRATION"
> da.y_label.text="H+,mol/l ";
> scf(2);clf; //Opens and clears figure 2
> plot(t,y(2,:))
>
> da.title.text="PENICILLIN CONCENTRATION"
> da.y_label.text="P,g/l ";
> scf(3);clf; //Opens and clears figure 3
> plot(t,y(3,:))
>
> da.title.text="SUBSTRATE CONCENTRATION"
> da.y_label.text="S,g/l ";
> scf(4);clf; //Opens and clears figure 4
> plot(t,y(4,:))
>
> da.title.text="DISSOLVED OXYGEN CONCENTRATION"
> da.y_label.text="C_l,g/l ";
> scf(5);clf; //Opens and clears figure 5
> plot(t,y(5,:))
>
> da.title.text="CULTURE VOLUME"
> da.y_label.text="V,l";
> scf(6);clf; //Opens and clears figure 6
> plot(t,y(6,:))
>
> da.title.text="HEAT OF REACTION"
> da.y_label.text="Qrxn,cal";
> scf(7);
> clf; //Opens and clears figure 7
> plot(t,y(7),:)
>
> da.title.text="TEMPERATURE"
> da.y_label.text="T,Kelvin";
> scf(8);
> clf; //Opens and clears figure 8
> plot(t,y(8),:)
>
> da.title.text="CO2 EVOLUTION"
> da.y_label.text="CO2,mmol/l/";
> scf(9);
> clf; //Opens and clears figure 9
> plot(t,y(9),:)
>
> Am i doing something wrong? before the ODE's i have just programmed
> the initial values and constants :
>
> funcprot(0);
> function dx = f(t,x)
> K1=1.0e-10 //mol/l
> K2=7.0e-05 //mol/l
> Kx=0.15 // Contois saturation constant, g/l
> Kox=2e-02 // oxygen limitation constant
> mux=0.092 // maitenance coefficient on subsrate
> p=3 //constant
> Kp=0.0002 // inhibition constant
> Kop=2e-02 // oxygen limitation constant
> K=0.04 // Penicillin hydrolysis constant, per h
> Yxs=0.45 // yield constant,g biomass/g glucose = dimensionless
> Yps=0.90 // yield constant, g pinicillin/ g glucose = dimensionless
> mx= 0.014 // Maintenance coefficient on substrate, per h
> Yxo=0.04 // yield constant, g biomass/g oxygen = dimensionless
> Ypo=0.20 // yield constant, g penicillin/g oxygen= dimensionless
> mo= 0.467 // maintenance coefficient of oxygen, per h
> mup=0.0005 // specific rate of penicilline production (per h)
> sf=600 // Feed substrate concentration, g/l
> kla=23 // function of agitation power input and oxugen flow rate,
> dimensional
> cll=1.16 // dissolved oxygen concentration, g/l
> Cab=3 // concentrations in both solutions
> Fa=5 // acid flow rate, l/h !!
> Fb=5 // base flow rate, l/h !!
> delta_t=0.01 // time step in digital PID controller - arbitrary
> value!!!
> z=10e-5 // constant
> F=0.042 // feed substrate flow rate l/h
> T0=273 // temperature at freezing, K
> Tv=373 // temperature at boiling
> T=298 // feed temp of substrate
> h=(2.5e-4) // constant
> Floss=(x(6)*(h)*(exp(5)*((T-T0)/(Tv-T0))))
> Fab=Fa+Fb // volume increase due to influx of acid Fa and base Fb
> Fsf=((sf)*(F))
> kg= 7e-3 // Arrhenius constant for growth
> kd=10e33 // Arrhenius constant for cell death
> Eg= 5100 // Activation energy for growth, cal/mol
> Ed= 50000 // Activation energy for cell death, cal/mol
> R= 1.987 // gas constant, cal/mol k
> T= 297 // Temperature
> RT= R*T
> alpa= 70 // constant in Kla
> betha= 0.4 // constant in Kla
> Pw= 30 // Agitation power input, W
> fg= 8.6 // Flow rate of oxygen
> V=100 // Volume
> QE= ((kg*exp(-(Eg/RT)))-(kd*exp(-(Ed/RT))))
> kla= alpa*((sqrt(fg)*(Pw/x(6)))^betha)
> mu
> =(((mux)/(1+((K1)/(x(2)))+((x(2))/(K2))))*((x(3))/(((Kx)*(x(1)))+(x(3))))*((x(5))/(((Kox)*(x(1)))+(x(5))))*(QE))
> // Specific growth rate
> mupp =
> ((mup)*((x(4))/((Kp)+(x(4))+(x(4)^2)/(K1)))*((x(5)^p)/((Kop)*(x(1)))+(x(5)^p)))
> // Specific penicillin production rate
> B =(((1e-14/x(2)-x(2))*x(6)-Cab*(Fa+Fb)*delta_t)/(x(6)+(Fa+Fb)*delta_t))
> QQ =((-B+sqrt(B^2+4e-14))/2-(x(2)))*(1/delta_t)
> dx_6 = (F+Fab+Floss) //Culture Volume V
> dx_1 = (((mu)*(x(1)))-((x(1))/(x(6)))*(dx_6)) //biomass concentration X
> rq1 = 60 // yield of heat generation, cal/g biomass
> rq2 = 1.6783e-4 // Constant, cal/g biomass h
> Tf = 296 // substrate feed temperature, Kelvin
> a = 1000 // heat transfer coefficient of cooling/heating liquid,
> cal/h degree C
> b = 0.60 // constant
> Fc=0.1 // Cooling water flow rate, not sure about value, l/h
> pcCpc = 1/2000 // Density times heat capacity of cooling liquid,
> per l degree C
> pcp = 1/1500 // density times heat capacity of medium
> QT = ((x(7)-(((a)*(Fc^b+1))/((Fc)+((a)*(Fc^b))/2*pcCpc))))
> a1=0.143 // constant relating CO2 to growth, mmol CO2/g biomass
> a2=4e-7 // Constant relating CO2 to mainteneance energy, mmol CO2/g
> biomass h
> a3=1e-4 // Constant relating CO2 to penicillin production, mmol CO2/l h
> CO= (((a1)*(dx_1))+((a2)*(x(1)))+(a3)), // CO2 evolution, CO2
> HION=((z*(((mu)*(x(1)))-(((F)*(x(1)))/(x(6)))))+(QQ))
>
> Thanks.
--
Adrien Vogt-Schilb (Cired)
Tel: (+33) 1 43 94 *73 77*
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.scilab.org/pipermail/users/attachments/20111005/356be67a/attachment.htm>
More information about the users
mailing list