 Timestamp:
 Nov 26, 2007, 9:50:55 AM (15 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/eml/reactors/equil.mso
r403 r415 17 17 * 18 18 * 19 * 19 * Description: 20 * Thermodynamic equilibrium modeling of a reactor based on 21 * equilibrium constants approach. 22 * 23 * Assumptions: 24 * * thermodynamic equilibrium 25 * * steadystate 26 * 27 * Specify: 28 * * inlet stream 29 * * stoichiometric matrix 30 * * equilibrium temperature 20 31 * 21 32 * … … 30 41 * only vapour phase 31 42 **# 32 33 43 Model equil_vap as tank_vap 34 44 ATTRIBUTES 35 45 Pallete = true; 36 46 Icon = "icon/cstr"; 37 Brief = "Model of a generic vapourphase Equilibrium CSTR";47 Brief = "Model of a generic vapourphase equilibrium CSTR"; 38 48 Info = " 39 49 Requires the information of: 40 50 * number of reactions 41 * matrix of stoichiometric coefficients (compo unds, reaction)51 * matrix of stoichiometric coefficients (components by reactions) 42 52 "; 43 53 … … 47 57 Rg as Real (Brief="Universal gas constant", Unit='J/mol/K', Default=8.314); 48 58 fs(NComp) as pressure (Brief="Fugacity in standard state", Default=1, DisplayUnit='atm'); 49 T 0as temperature (Brief="Reference temperature", Default=298.15);59 To as temperature (Brief="Reference temperature", Default=298.15); 50 60 51 61 VARIABLES 52 62 out Outlet as vapour_stream; # Outlet stream 53 63 54 G s(NComp) as energy_mol (Brief="Gibbs energy in standard state");64 G(NComp) as energy_mol (Brief="Gibbs freeenergy of formation"); 55 65 K(NReac) as Real (Brief="Equillibrium constant",Default=1.5); 56 66 activ(NComp)as Real (Brief="Activity",Default=0.2); 57 phi(NComp) as fugacity (Brief="Fugacity coefficient", Default=1);58 67 59 68 rate(NComp) as reaction_mol (Brief="Overall component rate of reaction"); … … 72 81 73 82 "Steadystate" 74 Outlet.F = Inlet.F+sum(sumt(stoic*extent)); 75 76 "Gibbs Energy of Formation" 77 Gs = PP.IdealGasGibbsOfFormation(Outlet.T); 78 79 # "Gibbs Energy of Formation without Cp correction" 80 # Gs = PP.IdealGasGibbsOfFormationAt25C()*Outlet.T/T0+PP.IdealGasEnthalpyOfFormationAt25C()*(1Outlet.T/T0); 81 83 Outlet.F = Inlet.F + sum(sumt(stoic*extent)); 84 85 "Gibbs freeenergy of formation" 86 G = PP.IdealGasGibbsOfFormation(Outlet.T); 87 88 # "Gibbs freeenergy of formation without Cp correction" 89 # G = PP.IdealGasGibbsOfFormationAt25C()*Outlet.T/To+PP.IdealGasEnthalpyOfFormationAt25C()*(1Outlet.T/To); 90 91 "Gibbs free energy of reaction" 92 # sumt(G*stoic) = Rg*Outlet.T*ln(K); 93 K = exp(sumt(G*stoic)/(Rg*Outlet.T)); 94 95 for j in [1:NReac] 82 96 "Equilibrium constant" 83 K = exp(sumt(Gs*stoic)/(Rg*Outlet.T));84 85 for j in [1:NReac]86 "Equilibrium"87 97 K(j) = prod(activ^stoic(:,j)); 88 98 end … … 107 117 end 108 118 109 "Fugacity coefficient"110 phi = PP.VapourFugacityCoefficient(Outlet.T,Outlet.P,Outlet.z);111 112 119 "Activity" 113 activ = phi*Outlet.P*Outlet.z/fs;120 activ = PP.VapourFugacityCoefficient(Outlet.T,Outlet.P,Outlet.z)*Outlet.P*Outlet.z/fs; 114 121 end 122 123 124 #* 125 * only liquidphase 126 **# 127 Model equil_liq as tank_liq 128 ATTRIBUTES 129 Pallete = true; 130 Icon = "icon/cstr"; 131 Brief = "Model of a generic liquidphase equilibrium CSTR"; 132 Info = " 133 Requires the information of: 134 * number of reactions 135 * matrix of stoichiometric coefficients (components by reactions) 136 "; 137 138 PARAMETERS 139 NReac as Integer (Brief="Number of reactions", Default=1); 140 stoic(NComp,NReac) as Real (Brief="Stoichiometric matrix"); 141 Rg as Real (Brief="Universal gas constant", Unit='J/mol/K', Default=8.314); 142 Ps as pressure (Brief="Standard pressure", Default=1, DisplayUnit='bar'); 143 To as temperature (Brief="Reference temperature", Default=298.15); 144 145 VARIABLES 146 out Outlet as liquid_stream; # Outlet stream 147 148 G(NReac) as enth_mol (Brief="Gibbs freeenergy of formation"); 149 K(NReac) as fraction (Brief="Equillibrium constant"); 150 activ(NComp)as Real (Brief="Activity"); 151 152 rate(NComp) as reaction_mol (Brief="Overall component rate of reaction"); 153 extent(NReac)as flow_mol (Brief="Extent of reaction"); 154 conv(NComp) as Real (Brief="Fractional conversion of component", Default=0); 155 156 EQUATIONS 157 "Outlet stream" 158 Outlet.F*Outlet.z = Outletm.F*Outletm.z + rate*V; 159 160 "Mechanical equilibrium" 161 Outlet.P = Outletm.P; 162 163 "Energy balance" 164 Outlet.F*Outlet.h = Outletm.F*Outletm.h; 165 166 "Steadystate" 167 Outlet.F = Inlet.F + sum(sumt(stoic*extent)); 168 169 "Gibbs freeenergy of formation" 170 G = PP.IdealGasGibbsOfFormation(Outlet.T); 171 172 # "Gibbs freeenergy of formation without Cp correction" 173 # G = PP.IdealGasGibbsOfFormationAt25C()*Outlet.T/To+PP.IdealGasEnthalpyOfFormationAt25C()*(1Outlet.T/To); 174 175 "Gibbs free energy of reaction" 176 # sumt(G*stoic) = Rg*Outlet.T*ln(K); 177 K = exp(sumt(G*stoic)/(Rg*Outlet.T)); 178 179 for j in [1:NReac] 180 "Equilibrium constant" 181 K(j) = prod(activ^stoic(:,j)); 182 end 183 184 for i in [1:NComp] 185 "Outlet molar fraction" 186 Outlet.F*Outlet.z(i) = (Inlet.F*Inlet.z(i) + sumt(stoic(i,:)*extent)); 187 end 188 189 for i in [1:NComp] 190 if (Outletm.z(i) > 0) then 191 "Molar conversion" 192 Outlet.F*Outlet.z(i) = Outletm.F*Outletm.z(i)*(1  conv(i)); 193 else if (Outlet.z(i) > 0) then 194 "Molar conversion" 195 conv(i) = 1; # ? 196 else 197 "Molar conversion" 198 conv(i) = 0; # ? 199 end 200 end 201 end 202 203 "Activity" 204 activ = PP.LiquidFugacityCoefficient(Outlet.T,Outlet.P,Outlet.z)*Outlet.z 205 *exp(PP.LiquidVolume(Outlet.T,Outlet.P,Outlet.z)*(Outlet.P  Ps)/Rg/Outlet.T); 206 end
Note: See TracChangeset
for help on using the changeset viewer.