Guest User

pha_cont

a guest
Sep 9th, 2019
34
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. function f=pha_cont(t,y)
  2.  
  3. %Pure culture PHA accumulation model written for CHEE4001
  4. %Damien Batstone 27/8/2019
  5. %Mass basis describes growth of relevant biomass and accumulation
  6. %In fixed volume two-stage chemostat system
  7.  
  8. %local variable assignment
  9. %reactor 1
  10. Ss_1 = y(1) ; %substrate
  11. Snh_1 = y(2); %ammonia
  12. Xb_1 = y(3); %bacteria
  13. Xphb_1 = y(4); %PHB
  14.  
  15. %reactor 2
  16. Ss_2 = y(5) ; %substrate
  17. Snh_2 = y(6); %ammonia
  18. Xb_2 = y(7); %bacteria
  19. Xphb_2 = y(8); %PHB
  20.  
  21. %Inputs
  22. qin_1 = 78161; %L/d fresh feed to R1
  23. qin_2 = 16332; %L/d fresh feed to R2
  24. Ss_in = 206; %g/L glucose feed concentration
  25. Snh_in_1 = 3; %g/L ammonia feed concentration to R1
  26.  
  27. %Parameters
  28. %physical
  29. V=[78161,94493] ; %L
  30. %V=[10000,15000] ; %L
  31.  
  32. %Biochemical original units
  33. mumax_gs = 0.4321*24; %d-1 Maximum biomass growth rate
  34. mumax_ps = 0.2628*24; %d-1 Maximum PHA production rate
  35. mumax_gp = 0.1566*24; %d-1 Maximum biomass growth rate
  36. KS_s_g = 1.2; %kg/m3 Sat coeff for growth
  37. KI_s_g = 16.728; %kg/m3 Sat inhib coeff for growth
  38. KS_s_p = 4.13; %kg/m3 Sat coeff for PHA production
  39. KS_nh_g = 0.254; %kg/m3 Ammonia sat coeff for growth
  40. KI_nh_g = 1.691; %kg/m3 Ammonia inhib coeff for growth
  41. %KI_nh_p = 80; %kg/m3 Ammonia inhib coeff for PHA uptake
  42. KI_nh_p = 0.2576; %kg/m3 Ammonia inhib coeff for PHA uptake
  43. Xmax = 68; %g/L Maximum concentration of Biomass
  44. PHAmax = 0.75; %max fraction of PHA content
  45. theta = 5.8; %Inhibitory effect index
  46.  
  47. %local constitutive eqs
  48. frac_PHA_1 = Xphb_1/(Xb_1+Xphb_1);
  49. frac_PHA_2 = Xphb_2/(Xb_2+Xphb_2);
  50.  
  51. %stoichiometry from Excel
  52. %Masters 1 values (fed batch)
  53. %v = [-2.083333333 -0.806293019 -0.150442478 1 0 1.108652901
  54. %-1.886792453 -0.338160012 0 0 1 0.72078397];
  55.  
  56. %2S-CSTR values
  57. v = [-1.779878826 -0.482608211 -0.150442478 1 0 0.66358629
  58. -2.777777778 -1.288544358 0 0 1 2.027562446];
  59.  
  60.  
  61. %rates
  62. % reactor 1
  63. r1_1 = mumax_gs*Xb_1*(Ss_1/(KS_s_g+Ss_1))*(Snh_1/(KS_nh_g+Snh_1))
  64. %r1_1 = mumax_gs * (1-((Xb_1/Xmax)^theta)) * (Ss_1/(KS_s_g+Ss_1+((Ss_1^2)/KI_s_g))) * (Snh_1/(KS_nh_g+Snh_1+((Snh_1^2)/KI_nh_g)))
  65. r2_1 = mumax_ps * Xb_1*Ss_1/(KS_s_p+Ss_1)*KI_nh_p/(KI_nh_p+Snh_1)*PHAmax/(PHAmax+frac_PHA_1)
  66.  
  67. % reactor 2
  68. r1_2 = mumax_gs*Xb_2*Ss_2/(KS_s_g+Ss_2)*Snh_2/(KS_nh_g+Snh_2)
  69. %r1_2 = mumax_gs * (1-((Xb_2/Xmax)^theta)) * (Ss_2/(KS_s_g+Ss_1+((Ss_1^2)/KI_s_g))) * (Snh_2/(KS_nh_g+Snh_2+((Snh_2^2)/KI_nh_g)))
  70. r2_2 = mumax_ps*Xb_2*Ss_2/(KS_s_p+Ss_2)*KI_nh_p/(KI_nh_p+Snh_2)*PHAmax/(PHAmax+frac_PHA_2)
  71.  
  72. % SS SNH XB XPHB Add SO2 SCO2 if you need them.
  73.  
  74. %derivs
  75. f(1) = qin_1/V(1)*(Ss_in-Ss_1) + v(1,1)*r1_1 + v(2,1)*r2_1;
  76. %skip O2
  77. %f(9) = KLa_O*(O_in-O)*(mumax_gs
  78. f(2) = qin_1/V(1)*(Snh_in_1-Snh_1) + v(1,3)*r1_1 + v(2,3)*r2_1;
  79. f(3) = qin_1/V(1)*(0-Xb_1) + v(1,4)*r1_1 + v(2,4)*r2_1;
  80. f(4) = qin_1/V(1)*(0-Xphb_1) + v(1,5)*r1_1 + v(2,5)*r2_1;
  81. %skip co2
  82. f(5) = qin_1/V(2)*(Ss_1-Ss_2)+qin_2/V(2)*(Ss_in-Ss_2)+v(1,1)*r1_2+v(2,1)*r2_2;
  83. %skip O2
  84. f(6) = qin_1/V(2)*(Snh_1-Snh_2)+qin_2/V(2)*(0-Snh_2)+v(1,3)*r1_2+v(2,3)*r2_2;
  85. f(7) = qin_1/V(2)*(Xb_1-Xb_2)+qin_2/V(2)*(0-Xb_2)+v(1,4)*r1_2+v(2,4)*r2_2;
  86. f(8) = qin_1/V(2)*(Xphb_1-Xphb_2)+qin_2/V(2)*(0-Xphb_2)+v(1,5)*r1_2+v(2,5)*r2_2;
  87. %skip co2
  88.  
  89. f=f';
  90. end
RAW Paste Data