/[ascend]/trunk/models/kinetics.a4l
ViewVC logotype

Annotation of /trunk/models/kinetics.a4l

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2651 - (hide annotations) (download) (as text)
Thu Dec 13 07:29:48 2012 UTC (11 years, 6 months ago) by jpye
File MIME type: text/x-ascend
File size: 9982 byte(s)
Fixing GPL header, removing postal address (rpmlint incorrect-fsf-address)
1 aw0a 1 REQUIRE "atoms.a4l";
2     (* => atoms.a4l, measures.a4l, system.a4l, basemodel.a4l *)
3     PROVIDE "kinetics.a4l";
4    
5     (*
6     * kinetics.a4l
7     * by Duncan Coffey
8     * Part of the ASCEND Library
9     * $Date: 1998/06/17 19:10:29 $
10     * $Revision: 1.3 $
11     * $Author: mthomas $
12     * $Source: /afs/cs.cmu.edu/project/ascend/Repository/models/kinetics.a4l,v $
13     *
14     * This file is part of the ASCEND Modeling Library.
15     *
16     * Copyright (C) 1998 Duncan Coffey
17     *
18     * The ASCEND Modeling Library is free software; you can redistribute
19     * it and/or modify it under the terms of the GNU General Public
20     * License as published by the Free Software Foundation; either
21     * version 2 of the License, or (at your option) any later version.
22     *
23     * The ASCEND Modeling Library is distributed in hope that it
24     * will be useful, but WITHOUT ANY WARRANTY; without even the implied
25     * warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
26     * See the GNU General Public License for more details.
27     *
28     * You should have received a copy of the GNU General Public License
29 jpye 2651 * along with this program. If not, see <http://www.gnu.org/licenses/>.
30 aw0a 1 *)
31    
32    
33     (* ************************************************************************ *)
34     (* ******************** Kinetics ******************** *)
35     (* ************************************************************************ *)
36     MODEL base_kinetics(
37     components IS_A set OF symbol_constant;
38     nr IS_A set OF symbol_constant;
39     Kr[nr] IS_A constant;
40     active_energy[nr] IS_A constant;
41     reac_T WILL_BE temperature;
42     species[nr] IS_A set OF symbol_constant;
43     nu[components][nr] WILL_BE integer_constant;
44     conc[components] WILL_BE molar_density;
45     )WHERE(
46    
47     )REFINES cmumodel;
48     NOTES
49     'purpose' SELF {
50     This MODEL computes the production rates per volume of the reactant
51     speices given. There are three possible models to use.
52    
53     base_kinetics - No rate law is constructed and is expected to be given
54     by the user. Production rates are then calculated from the given rate laws.
55    
56     nr - names of reactions
57     Kr[nr] - rate constants
58     species[nr] - set of components involved in each reaction
59     nu[components][nr] - stoichiometry for each components and reaction
60     conc[components] - Concentrations used in rate laws.
61    
62     element_kinetics - Here the rate laws are constructed from the
63     stoichiometry (nu) given.
64    
65     specify_kinetics - Here the rate laws are constructed from the orders (order) given.
66    
67     }
68     'developer-Duncan' SELF {From the element_kinetics and
69     specify_kinetics most rate laws can be constructed. However in many cases
70     there may need to be intermediate steps and it might just be easier to
71     use base_kinetics.
72     }
73     END NOTES;
74    
75     rate[nr] IS_A conc_rate;
76     reactions[components] IS_A set OF symbol_constant;
77     production[components] IS_A conc_rate;
78     Ftot_scale IS_A scaling_constant;
79 jose 2071 R IS_A molar_gas_constant;
80 aw0a 1
81     (* Define the production rate for a component *)
82     FOR i IN components CREATE
83     production[i] / Ftot_scale = SUM[nu[i][j]*rate[j] | j IN reactions[i]] / Ftot_scale;
84     END FOR;
85    
86     (* define rxns where each component is present *)
87     FOR i IN components CREATE
88     reactions[i] :== [j IN nr | i IN species[j]];
89     END FOR;
90    
91     METHODS
92    
93     METHOD default_self;
94     rate[nr].lower_bound :=-1e100 {mole/m^3/s};
95     production[components].lower_bound :=-1e100 {mole/m^3/s};
96     Ftot_scale :=1000;
97     END default_self;
98    
99     METHOD default_all;
100     RUN reac.default_self;
101     END default_all;
102    
103     METHOD check_self;
104     END check_self;
105    
106     METHOD check_all;
107     RUN check_self;
108     END check_all;
109    
110     METHOD bound_self;
111     END bound_self;
112    
113     METHOD bound_all;
114     RUN bound_self;
115     END bound_all;
116    
117     METHOD scale_self;
118     END scale_self;
119    
120     METHOD scale_all;
121     RUN scale_self;
122     END scale_all;
123    
124     METHOD seqmod;
125 johnpye 576 FIX conc[components];
126     FIX reac_T;
127 aw0a 1 END seqmod;
128    
129     METHOD specify;
130     RUN seqmod;
131     END specify;
132     END base_kinetics;
133    
134    
135     MODEL element_kinetics(
136     components IS_A set OF symbol_constant;
137     nr IS_A set OF symbol_constant;
138     Kr[nr] IS_A constant;
139     active_energy[nr] IS_A constant;
140     reac_T WILL_BE temperature;
141     species[nr] IS_A set OF symbol_constant;
142     nu[components][nr] WILL_BE integer_constant;
143     conc[components] WILL_BE molar_density;
144     )WHERE(
145    
146     )REFINES base_kinetics;
147    
148     (* define the rate equations for each component *)
149     FOR j IN nr CREATE
150     rate[j] / Ftot_scale = (Kr[j]*exp(-active_energy[j]/R/reac_T)*PROD[
151     PROD[ conc[i] | m IN [1..-(nu[i][j])]] | i IN species[j]]) / Ftot_scale;
152     END FOR;
153    
154    
155    
156     METHODS
157    
158     METHOD default_self;
159     rate[nr].lower_bound :=-1e100 {mole/m^3/s};
160     production[components].lower_bound :=-1e100 {mole/m^3/s};
161     Ftot_scale :=1000;
162     END default_self;
163    
164     METHOD default_all;
165     RUN reac.default_self;
166     END default_all;
167    
168     METHOD check_self;
169     END check_self;
170    
171     METHOD check_all;
172     RUN check_self;
173     END check_all;
174    
175     METHOD bound_self;
176     END bound_self;
177    
178     METHOD bound_all;
179     RUN bound_self;
180     END bound_all;
181    
182     METHOD scale_self;
183     END scale_self;
184    
185     METHOD scale_all;
186     RUN scale_self;
187     END scale_all;
188    
189     METHOD seqmod;
190 johnpye 576 FIX conc[components];
191     FIX reac_T;
192 aw0a 1 END seqmod;
193    
194     METHOD specify;
195     RUN seqmod;
196     END specify;
197    
198     END element_kinetics;
199    
200    
201     MODEL specify_kinetics(
202     components IS_A set OF symbol_constant;
203     nr IS_A set OF symbol_constant;
204     Kr[nr] IS_A constant;
205     active_energy[nr] IS_A constant;
206     reac_T WILL_BE temperature;
207     species[nr] IS_A set OF symbol_constant;
208     nu[components][nr] WILL_BE integer_constant;
209     conc[components] WILL_BE molar_density;
210     order[components][nr] WILL_BE constant;
211     )WHERE(
212    
213     )REFINES base_kinetics;
214    
215    
216     (* define the rate equations for each component *)
217     FOR j IN nr CREATE
218     rate[j] / Ftot_scale = (Kr[j]*exp(-active_energy[j]/R/reac_T)*
219     PROD[conc[i]^order[i][j] | i IN species[j]]) / Ftot_scale;
220     END FOR;
221    
222     METHODS
223    
224     METHOD default_self;
225     rate[nr].lower_bound :=-1e100 {mole/m^3/s};
226     production[components].lower_bound :=-1e100 {mole/m^3/s};
227     Ftot_scale :=1000;
228     END default_self;
229    
230     METHOD seqmod;
231 johnpye 576 FIX conc[components];
232     FIX reac_T;
233 aw0a 1 END seqmod;
234    
235     METHOD specify;
236     RUN seqmod;
237     END specify;
238    
239     END specify_kinetics;
240    
241     MODEL test_elm_kinetics;
242     components IS_A set OF symbol_constant;
243     nr IS_A set OF symbol_constant;
244     Kr[nr] IS_A constant;
245     active_energy[nr] IS_A constant;
246     reac_T IS_A temperature;
247     species[nr] IS_A set OF symbol_constant;
248     nu[components][nr] IS_A integer_constant;
249     conc[components] IS_A molar_density;
250    
251    
252     components :==['methane','ethane','n_butane'];
253     nr :==['reaction_1','reaction_2'];
254     species['reaction_1'] :==['methane','ethane'];
255     species['reaction_2'] :==['methane','ethane'];
256     Kr['reaction_1'] :==2 {m^3/mol/s};
257     Kr['reaction_2'] :==3 {1/s};
258     active_energy['reaction_1'] :==70000 {J};
259     active_energy['reaction_2'] :==70000 {J};
260    
261     nu['methane']['reaction_1'] :==-2;
262     nu['ethane']['reaction_1'] :==1;
263     nu['methane']['reaction_2'] :==2;
264     nu['ethane']['reaction_2'] :==-1;
265     (* define nu=0 FOR components not IN reactions*)
266     FOR i IN nr CREATE
267     FOR j IN components - [species[i]] CREATE
268     nu[j][i] :==0;
269     END FOR;
270     END FOR;
271    
272    
273    
274     reac IS_A element_kinetics(components, nr, Kr, active_energy,
275     reac_T, species, nu, conc);
276    
277     METHODS
278    
279     METHOD default_self;
280     RUN reac.default_self;
281     END default_self;
282    
283     METHOD default_all;
284     RUN reac.default_self;
285     END default_all;
286    
287     METHOD check_self;
288     END check_self;
289    
290     METHOD check_all;
291     RUN check_self;
292     END check_all;
293    
294     METHOD bound_self;
295     END bound_self;
296    
297     METHOD bound_all;
298     RUN bound_self;
299     END bound_all;
300    
301     METHOD scale_self;
302     END scale_self;
303    
304     METHOD scale_all;
305     RUN scale_self;
306     END scale_all;
307     END test_elm_kinetics;
308    
309     MODEL test_ord_kinetics;
310     components IS_A set OF symbol_constant;
311     nr IS_A set OF symbol_constant;
312     Kr[nr] IS_A constant;
313     active_energy[nr] IS_A constant;
314     reac_T IS_A temperature;
315     species[nr] IS_A set OF symbol_constant;
316     nu[components][nr] IS_A integer_constant;
317     order[components][nr] IS_A constant;
318     conc[components] IS_A molar_density;
319    
320    
321     components :==['methane','ethane','n_butane'];
322     nr :==['reaction_1','reaction_2'];
323     species['reaction_1'] :==['methane','ethane'];
324     species['reaction_2'] :==['methane','ethane'];
325     Kr['reaction_1'] :==2 {m^3/mol/s};
326     Kr['reaction_2'] :==3 {1/s};
327     active_energy['reaction_1'] :==70000 {J};
328     active_energy['reaction_2'] :==70000 {J};
329    
330     nu['methane']['reaction_1'] :==-2;
331     nu['ethane']['reaction_1'] :==1;
332     nu['methane']['reaction_2'] :==2;
333     nu['ethane']['reaction_2'] :==-1;
334     (* define nu=0 FOR components not IN reactions*)
335     FOR i IN nr CREATE
336     FOR j IN components - [species[i]] CREATE
337     nu[j][i] :==0;
338     END FOR;
339     END FOR;
340    
341     order['methane']['reaction_1'] :==2;
342     order['ethane']['reaction_1'] :==0;
343     order['methane']['reaction_2'] :==3;
344     order['ethane']['reaction_2'] :==1;
345     FOR i IN nr CREATE
346     FOR j IN components - [species[i]] CREATE
347     order[j][i] :==0;
348     END FOR;
349     END FOR;
350    
351     reac IS_A specify_kinetics(components, nr, Kr, active_energy,
352     reac_T, species, nu, conc, order);
353    
354     METHODS
355     METHOD default_self;
356     RUN reac.default_self;
357     END default_self;
358    
359     METHOD default_all;
360     RUN reac.default_self;
361     END default_all;
362    
363     METHOD check_self;
364     END check_self;
365    
366     METHOD check_all;
367     RUN check_self;
368     END check_all;
369    
370     METHOD bound_self;
371     END bound_self;
372    
373     METHOD bound_all;
374     RUN bound_self;
375     END bound_all;
376    
377     METHOD scale_self;
378     END scale_self;
379    
380     METHOD scale_all;
381     RUN scale_self;
382     END scale_all;
383     END test_ord_kinetics;

john.pye@anu.edu.au
ViewVC Help
Powered by ViewVC 1.1.22