[scilab-Users] No reaction
jaundreventer at gmail.com
jaundreventer at gmail.com
Tue Nov 1 16:18:22 CET 2011
Hi adrien
Thank you very much for your help. I will try your tips.
Thank you once again.
Regards.
Sent via my BlackBerry from Vodacom - let your email find you!
-----Original Message-----
From: Adrien Vogt-Schilb <vogt at centre-cired.fr>
Date: Mon, 31 Oct 2011 23:40:04
To: <users at lists.scilab.org>
Reply-To: users at lists.scilab.org
Subject: Re: [scilab-Users] No reaction
Hi Jaundre Venter
your code runs just fine
you just ploted y(9) and y(8) instead of y(9,:) and y(8,:) !
it happens that when M is a matrix, say a 3x4 matrix, things like M($),
M(1), M(7) or M([1 3 5 9]) are legal scilab statements, so you don't
have any error message.
see the documentation: help extraction
*scilab team*, maybe you could put an option for displaying warnings
when tables are accessed with less arguments than their dimensions, what
do you think?
On 31/10/2011 21:10, Jaundre Venter wrote:
> Hi Adrien
>
> Sorry for late response but had some issues on other projects that i
> had to attend.
>
> i am still getting no reactions what so ever for the last three
> curves. see attached the new code. i don't know it is showing 0 or at
> the initial values! i have checked the mass balance if some input
> values are wrong that is causing the no reaction but the constant etc
> is correct.
>
> On Thu, Oct 6, 2011 at 2:07 PM, Adrien Vogt-Schilb
> <vogt at centre-cired.fr <mailto:vogt at centre-cired.fr>> wrote:
>
> Hi
>
> Maybe you should look at x (the output of ode) without ploting it.
> That way you'll make sure if the problem comes from the plot or not
>
> Can you attach a file with the new code?
>
>
> On 06/10/2011 10:00, Jaundre Venter wrote:
>> Hi
>>
>> okay i have made the changes and the reaction curves are looking
>> better but for some reason i still get no reaction on the last
>> three curves. It is just showing a straight line. can it be that
>> i have choosen to plot the wrong type of graph? Maybe because of
>> a time value or to much ODE"s?
>>
>> On Wed, Oct 5, 2011 at 10:24 PM, Jaundre Venter
>> <jaundreventer at gmail.com <mailto:jaundreventer at gmail.com>> wrote:
>>
>> Thank you very much Adrien.
>>
>> always nice if someone can explain to you where your problems
>> are and why. Thanks
>>
>> Was there any other problems you saw that i have to be aware of?
>>
>>
>> On Wed, Oct 5, 2011 at 5:38 PM, Adrien Vogt-Schilb
>> <vogt at centre-cired.fr <mailto:vogt at centre-cired.fr>> wrote:
>>
>> Hi
>>
>> When you use ode, it's ok, if say, dx(1) depends on dx(4).
>> but you still have say that to scilab properly, something
>> like:
>>
>> function dx = f(t,x)
>> dx(6)=((F+Fab+Floss)*(x(2))), // culture Volume V
>> dx(1)=(((mu)*(x(1)))-(((x(1))/(x(6)))*((dx(6))))*(CO)*(x(2)))),
>> //biomass concentration X
>> and so one
>>
>> note that because i had to know dx(6) to compute dx(1) i
>> just computed dx(6) before dx(1): no problem. and note
>> that i used x(2). The idea of the ode is to compute dx
>> from x!
>>
>> make sure you understand that using dx_6 instead of
>> dx(6), your ODE solver is not updating dx_6 at each time
>> step, it is using the initial and only dx_6 forever.
>> That's why your last varaibles do not move, somehow their
>> speeds are never updated.
>> for instance, dx(6)=((F+Fab+Floss)*(HION)), // culture
>> Volume V is constance in time (i guess)
>>
>>
>>
>> On 05/10/2011 15:00, Jaundre Venter wrote:
>>> Hi Adrien
>>>
>>> i am new to SCILAB! I just want to say that.
>>>
>>> Yes dx_1 is equal to dx1 but the only reason why i have
>>> programmed it like that is becasue the ODE's looks as
>>> follows - (see word file attached that can explain the
>>> ODE"s better. with regards to the HION and CO it
>>> actually refers to [H+] and CO2 as you said. the only
>>> reason why n multiplied HION and CO wit hsome ODE's is
>>> because some does have a influence on some ODE's.
>>>
>>> This is the first time i am working with SCILAB thus i
>>> am struggling to understand how SCILAB wants the code so
>>> that all 9 ODE's are shown and so that the ODE's that is
>>> having a effect on other does happen. I though if you
>>> refer to dx(1) for example in a other ODE it means that
>>> SCILAb will know the dx(1) has a influence on the other ODE.
>>>
>>> The main goal of my assesment is to deliver similar
>>> results obtained from MATLAB on SCILAB. all i got was
>>> how the grpahs should look like and the ODE's.
>>>
>>> On Wed, Oct 5, 2011 at 2:43 PM, Adrien Vogt-Schilb
>>> <vogt at centre-cired.fr <mailto:vogt at centre-cired.fr>> wrote:
>>>
>>> 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*
>>>
>>>
>>
>>
>> --
>> Adrien Vogt-Schilb (Cired)
>> Tel: (+33) 1 43 94 *73 77*
>>
>>
>>
>
>
> --
> Adrien Vogt-Schilb (Cired)
> Tel: (+33) 1 43 94 *73 77*
>
>
--
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/20111101/fb165fcf/attachment.htm>
More information about the users
mailing list