Question

Brewing   beer   involves   putting   sugar   (S)   and   yeast   (Y)   together   in   a   vessel   so   that   alcohol   (A)   is   produced   as   a   by-product   of   the   metabolism   of   sugar   by   yeast.   Two   important   aspects   of   fermentation   are   the   limited   amount   of   sugar   and   the   formation   of   excessive   alcohol   that   will   kill   yeast   cells.   The   model   has   some   analogy   to   enzyme   kinetics   (sugar   is   a   substrate   for   yeast)   and   predator   –   prey   interaction   (alcohol   preys   yeast).       The   set   of   differential   equations,   initial   conditions   and   constants   are:  
  
S=1000 (rel.   conc.)

Y   =100   (rel.   conc.)  

A   =   0   (rel.   conc.)  
  
dS/dt  = -(a*b*S*Y) - (a*f*S*Y)  

dY/dt  = a*c*S*Y  - (d*Y*A)   d

A/dt = a*b*S*Y  
  
Constants:   a=0.1   b=0.01   c=0.05   d=0.01   f=0.01  
  
1) Use   Euler’s   method   with   a   step   size   of   h   =   0.2   to   find   approximate   values   of   the   state   variables   S,   Y   and   A   from   t=0   up   to   t=5   hours.   Produce   a   graph   showing   the   change   of   state   variables.   Do   not   use   Matlab   ODE   solvers,   but   develop   your   own   Euler   solver.       2) Develop   a   Runga-Kutta   solver   (RK4)   with   a   step   size   of   h   =   0.2   .   Graph   the   functions   for   all   three   state   variables.       3) Determine   the   absolute   and   relative   errors   for   all   three   predicted   state   variables   at   t=5   between   Euler   and   RK4.       4) Summarize   your   results   with   2-3   sentences.   What   do   your   predictions   tell   us   about   the   process   of   fermentation?   To   what   degree   is   RK4   better   than   Euler   and   does   it   depend   on   the   behavior   of   the   state   variable?  
  
>>>   Include   your   Matlab   code   at   the   end   of   your   word-file   or   attach   matlab   m.file   >>>   Consult   publicly   available   Matlab   codes   for   Euler   and   RK4  

Answer 1

Code

solution

euller method b 0.01; c 0.05; d0.01 f=0.01; hvalue-0.2; % step size t-0 : 0 . 1 : 5; % time num=length (t) ; s1 1000 zeros (1, num-1); Y1 =[ 100 zeros (1, num-1)]; A1 0 zeros (1, num-1)1 % euler method : yl (num+ 1,-yl (num)+hvalue*dy for 1Ξ1 : num-1 ds 1--(a*b*S1 (i) *Y1 (i) ) - 81 (i+1)-S 1 (i) +hvalue * ds 1; (a*f*S1 (i) *Y1 (i) ) ; end Splot hold on plot (t,si, "b' plot (t,Y1, 'r') plot (t,A1, 'g title(' euler method)

xlabel ('time') ylabel ('Rel.Conc') Legend('Bugar', 'Yeast', 'aleohol')

//output

euler method 200 Sugar Yesst Alcohol 600_ 400

//copyable code

euler method

a=0.1;

b=0.01;

c=0.05;

d=0.01;

f=0.01;

hvalue=0.2; % step size

t=0:0.1:5; % time

num=length(t);

S1 =[1000 zeros(1,num-1)];

Y1 =[100 zeros(1,num-1)];

A1 =[0 zeros(1,num-1)];

% euler method: y1(num+1)=y1(num)+hvalue*dy

for i=1:num-1

dS1 = -(a*b*S1(i)*Y1(i)) - (a*f*S1(i)*Y1(i));

S1(i+1)=S1(i)+hvalue*dS1;

dY1 = a*c*S1(i)*Y1(i) - (d*Y1(i)*A1(i));

Y1(i+1)=Y1(i)+hvalue*dY1;

dA1 = a*b*S1(i)*Y1(i);

A1(i+1)=A1(i)+hvalue*dA1;

end

  

%plot

hold on

plot(t,S1,'b')

plot(t,Y1,'r')

plot(t,A1,'g')

title(' euler method')

xlabel('time')

ylabel('Rel.Conc')

legend('Sugar','Yeast','Alcohol')

Runga kutta Method clc;clear Constants : b 0.01; с 0.05; d=0.01; f-0.01 h-0 . 2; % size step t-0 : 0 . 1 : 5 ; % time vector num-length (t); % vector initialization s [ 1000 zeros (1, num-1) ] ; Y 1100 zeros (1, num-1) A [0 zeros (1, num-1); Runge Kutta 4 for i-1: num-1 k 1=h-ds 1 (t (i) ,S (i) ,y(i) ) ; k-2-h * ds 1 (t (i) +0 . 5th, s(i)+0.5+k-1 , Y (İ) +0.5* k-1 ) ; k 3-h*dsl (t (i) +0.5*h, s(i)+0.5 k_2, Y(i)+0.5*k_2); k 4-h*dS1 (t (i)+h, S(i)+k 3, Y (i) +k 3);

1-h * dY1 (t (i) ,S (i) ,Y(1) ,A (i) ) ; k k 1h*dA1 (t (i),S(i) , Y (i)); k 2-h-dA1 (t (i)+0.5+h, s(i)+0.5*k 1,y(i)+0.5*k 1); end %plotting graphs hold on plot (t,s,'r' plot (t, Y, 'b') plot (t,A, 'g') title('Runge-Kutta method') xlabel 't') ylabel ( ' Rel . Conc ' ) legend ('Sugar, 'Yeast, "Alcohol')

//output

Fermentation modelling with Runge-Kutta method 200 Sugar Yesst 1000 Alcohol 800 600_

//copyable code

euler method

a=0.1;

b=0.01;

c=0.05;

d=0.01;

f=0.01;

hvalue=0.2; % step size

t=0:0.1:5; % time

num=length(t);

S1 =[1000 zeros(1,num-1)];

Y1 =[100 zeros(1,num-1)];

A1 =[0 zeros(1,num-1)];

% euler method: y1(num+1)=y1(num)+hvalue*dy

for i=1:num-1

dS1 = -(a*b*S1(i)*Y1(i)) - (a*f*S1(i)*Y1(i));

S1(i+1)=S1(i)+hvalue*dS1;

dY1 = a*c*S1(i)*Y1(i) - (d*Y1(i)*A1(i));

Y1(i+1)=Y1(i)+hvalue*dY1;

dA1 = a*b*S1(i)*Y1(i);

A1(i+1)=A1(i)+hvalue*dA1;

end

  

%plot

hold on

plot(t,S1,'b')

plot(t,Y1,'r')

plot(t,A1,'g')

title(' euler method')

xlabel('time')

ylabel('Rel.Conc')

legend('Sugar','Yeast','Alcohol')