source: branches/newlanguage/sample/miscellaneous/tenprobs/prob07.mso @ 228

Last change on this file since 228 was 228, checked in by Rodolfo Rodrigues, 15 years ago

Add more detailed header

File size: 4.7 KB
Line 
1#*---------------------------------------------------------------------
2* EMSO Model Library (EML) Copyright (C) 2004 - 2007 ALSOC.
3*
4* This LIBRARY is free software; you can distribute it and/or modify
5* it under the therms of the ALSOC FREE LICENSE as available at
6* http://www.enq.ufrgs.br/alsoc.
7*
8* EMSO Copyright (C) 2004 - 2007 ALSOC, original code
9* from http://www.rps.eng.br Copyright (C) 2002-2004.
10* All rights reserved.
11*
12* EMSO is distributed under the therms of the ALSOC LICENSE as
13* available at http://www.enq.ufrgs.br/alsoc.
14*
15*----------------------------------------------------------------------
16* 7. Diffusion with chemical reaction in a one dimensional slab
17*----------------------------------------------------------------------
18*
19*   Description:
20*      This problem is part of a collection of 10 representative
21*       problems in Chemical Engineering for solution by numerical methods
22*       developed for Cutlip (1998).
23*
24*   Subject:
25*       * Transport Phenomena
26*               * Reaction Engineering
27*
28*       Concepts utilized:
29*               Methods for solving second order ODEs with 2 point boundary
30*       values typically used in transport phenomena and reaction kinetics.
31*
32*       Numerical method:
33*               * Simultaneous ODEs with split boundary conditions
34*               * Resolved by finite difference method
35*
36*   Reference:
37*               * CUTLIP et al. A collection of 10 numerical problems in
38*       chemical engineering solved by various mathematical software
39*       packages. Comp. Appl. in Eng. Education. v. 6, 169-180, 1998.
40*       * More informations and a detailed description of all problems
41*       is available online in http://www.polymath-software.com/ASEE
42*
43*----------------------------------------------------------------------
44* Author: Rodolfo Rodrigues
45* GIMSCOP/UFRGS - Group of Integration, Modeling, Simulation,
46*                                       Control, and Optimization of Processes
47* $Id$
48*--------------------------------------------------------------------*#
49using "types";
50
51
52
53Model problem
54        PARAMETERS
55outer N as Integer      (Brief="Number of discrete points", Lower=3);
56       
57        Co              as conc_mol             (Brief="Constant concentration at the surface");
58        D               as diffusivity  (Brief="Binary diffusion coefficient");
59        k               as Real                 (Brief="Homogeneous reaction rate constant", Unit='1/s');
60        L               as length               (Brief="Bottom surface");
61       
62       
63        VARIABLES
64        C(N+2)  as conc_mol             (Brief="Concentration of reactant");
65        z(N+2)  as length               (Brief="Distance", Default=1e-3);
66        dz              as length_delta (Brief="Distance increment");   
67
68
69        EQUATIONS
70        "Discrete interval"
71        dz = (z(N+2) - z(1))/(N+1);
72       
73        for i in [2:(N+1)]
74        "Concentration of reactant"     
75                (C(i+1) - 2*C(i)+ C(i-1))/(z(i) - z(i-1))^2
76                        = (k/D)*C(i);
77               
78        "Discrete length"
79                z(i) = z(i-1) + dz;
80        end
81
82        # Boundary conditions
83        "Initial and boundary condition" # z = 0
84        C(1) = Co;
85       
86        "Upper boundary" # z = L
87        (C(N+2) - C(N+1))/(z(N+2) - z(N+1)) = 0*'kmol/m^4';
88       
89       
90        SET
91        Co= 0.2*'kmol/m^3';
92        D = 1.2e-9*'m^2/s';
93        k = 1e-3/'s';
94        L = 1e-3*'m';
95end
96
97
98
99FlowSheet numerical_solution
100        PARAMETERS
101        N               as Integer;
102       
103       
104        DEVICES
105        reac    as problem;
106       
107       
108        SET
109        N = 10; # Number of discrete points
110       
111       
112        SPECIFY
113        "Lower boundary"
114        reac.z(1) = 0*'m';
115
116        "Upper boundary"
117        reac.z(N+2) = reac.L;
118
119
120        OPTIONS
121        Dynamic = false;
122end
123
124
125
126FlowSheet comparative
127        PARAMETERS
128        N               as Integer;
129       
130       
131        VARIABLES
132        C_(N+2) as conc_mol     (Brief="Concentration of reactant by analytical solution");
133       
134        r_              as Real         (Brief="Pearson product-moment correlation coefficient");
135       
136        Cm              as conc_mol     (Brief="Arithmetic mean of calculated C");
137        C_m             as conc_mol     (Brief="Arithmetic mean of analytical C");
138       
139       
140        DEVICES
141        reac    as problem;
142       
143       
144        SET
145        N = 10; # Number of discrete points
146       
147       
148        EQUATIONS
149        "Analytical solution"
150        C_ = reac.Co*cosh(reac.L*sqrt(reac.k/reac.D)*(1 - reac.z/reac.L))
151        /cosh(reac.L*sqrt(reac.k/reac.D));
152
153        "Pearson correlation coefficient" # used by softwares like MS Excel
154        r_ = (sum((reac.C - Cm)*(C_ - C_m)))/sqrt(sum((reac.C - Cm)^2)*sum((C_ - C_m)^2));
155       
156        "Arithmetic mean of C"
157        Cm = sum(reac.C)/(N+1);
158       
159        "Arithmetic mean of C_"
160        C_m = sum(C_)/(N+1);
161       
162       
163        SPECIFY
164        "Lower boundary"
165        reac.z(1) = 0*'m';
166
167        "Upper boundary"
168        reac.z(N+2) = reac.L;
169
170
171        OPTIONS
172        Dynamic = false;
173end
174
175
176FlowSheet analytical_solution
177        PARAMETERS
178        Co              as conc_mol;
179        L               as length;
180        k               as Real(Unit='1/s');
181        D               as diffusivity;
182       
183       
184        VARIABLES
185        C               as conc_mol             (Default=0.2);
186        z               as length               (Default=1e-3);
187       
188       
189        EQUATIONS
190        "Change time in z"
191        z = time*'m/s';
192       
193        "Analytical solution"
194        C = Co*cosh(L*sqrt(k/D)*(1 - z/L))/cosh(L*sqrt(k/D));
195       
196       
197        SET
198        Co= 0.2*'kmol/m^3';
199        D = 1.2e-9*'m^2/s';
200        k = 1e-3/'s';
201        L = 1e-3*'m';
202       
203       
204        OPTIONS
205        TimeStart = 0;
206        TimeStep = 1e-6;
207        TimeEnd = 1e-3;
208end
Note: See TracBrowser for help on using the repository browser.