CAMP 1.0.0
Chemistry Across Multiple Phases
sub_model_UNIFAC.c
Go to the documentation of this file.
1/* Copyright (C) 2021 Barcelona Supercomputing Center and University of
2 * Illinois at Urbana-Champaign
3 * SPDX-License-Identifier: MIT
4 *
5 * UNIFAC activity coefficient calculation
6 *
7 */
8/** \file
9 * \brief UNIFAC activity coefficient calculation
10 *
11 * For more info see the camp_sub_model_UNIFAC module
12 *
13 * Equation references are to Marcolli and Peter, ACP 5(2), 1501-1527, 2005.
14 *
15 */
16#include <math.h>
17#include <stdio.h>
18#include <stdlib.h>
19#include "../Jacobian.h"
20#include "../sub_models.h"
21
22// TODO Lookup environmental indicies during initialization
23#define TEMPERATURE_K_ env_data[0]
24#define PRESSURE_PA_ env_data[1]
25
26#define NUM_UNIQUE_PHASE_ (int_data[0])
27#define NUM_GROUP_ (int_data[1])
28#define TOTAL_INT_PROP_ (int_data[2])
29#define TOTAL_FLOAT_PROP_ (int_data[3])
30#define NUM_INT_PROP_ 4
31#define NUM_FLOAT_PROP_ 0
32#define PHASE_INT_LOC_(p) (int_data[NUM_INT_PROP_ + p] - 1)
33#define PHASE_FLOAT_LOC_(p) \
34 (int_data[NUM_INT_PROP_ + NUM_UNIQUE_PHASE_ + p] - 1)
35#define PHASE_ENV_LOC_(p) \
36 (int_data[NUM_INT_PROP_ + 2 * NUM_UNIQUE_PHASE_ + p] - 1)
37#define NUM_PHASE_INSTANCE_(p) (int_data[PHASE_INT_LOC_(p)])
38#define NUM_SPEC_(p) (int_data[PHASE_INT_LOC_(p) + 1])
39#define PHASE_INST_ID_(p, c) (int_data[PHASE_INT_LOC_(p) + 2 + c] - 1)
40#define SPEC_ID_(p, i) \
41 (int_data[PHASE_INT_LOC_(p) + 2 + NUM_PHASE_INSTANCE_(p) + i])
42#define GAMMA_ID_(p, i) \
43 (int_data[PHASE_INT_LOC_(p) + 2 + NUM_PHASE_INSTANCE_(p) + NUM_SPEC_(p) + i])
44#define JAC_ID_(p, c, j, i) \
45 int_data[PHASE_INT_LOC_(p) + 2 + NUM_PHASE_INSTANCE_(p) + \
46 (2 + c * NUM_SPEC_(p) + j) * NUM_SPEC_(p) + i]
47#define V_IK_(p, i, k) \
48 (int_data[PHASE_INT_LOC_(p) + 2 + NUM_PHASE_INSTANCE_(p) + \
49 (k + 2 + NUM_PHASE_INSTANCE_(p) * NUM_SPEC_(p)) * NUM_SPEC_(p) + \
50 i])
51
52#define Q_K_(k) (float_data[k])
53#define R_K_(k) (float_data[NUM_GROUP_ + k])
54#define X_K_(m) (float_data[2 * NUM_GROUP_ + m])
55#define DTHETA_M_DC_I_(m) (float_data[3 * NUM_GROUP_ + m])
56#define XI_M_(m) (float_data[4 * NUM_GROUP_ + m])
57#define LN_GAMMA_K_(m) (float_data[5 * NUM_GROUP_ + m])
58#define A_MN_(m, n) (float_data[(m + 6) * NUM_GROUP_ + n])
59#define R_I_(p, i) (float_data[PHASE_FLOAT_LOC_(p) + i])
60#define Q_I_(p, i) (float_data[PHASE_FLOAT_LOC_(p) + NUM_SPEC_(p) + i])
61#define L_I_(p, i) (float_data[PHASE_FLOAT_LOC_(p) + 2 * NUM_SPEC_(p) + i])
62#define MW_I_(p, i) (float_data[PHASE_FLOAT_LOC_(p) + 3 * NUM_SPEC_(p) + i])
63#define X_I_(p, i) (float_data[PHASE_FLOAT_LOC_(p) + 4 * NUM_SPEC_(p) + i])
64
65#define THETA_M_(m) (sub_model_env_data[m])
66#define PSI_MN_(m, n) (sub_model_env_data[(m + 1) * NUM_GROUP_ + n])
67#define LN_GAMMA_IK_(p, i, k) \
68 (sub_model_env_data[PHASE_ENV_LOC_(p) + i * NUM_GROUP_ + k])
69
70#define INT_DATA_SIZE_ (TOTAL_INT_PROP_)
71#define FLOAT_DATA_SIZE_ (TOTAL_FLOAT_PROP_)
72
73// Update types (These must match values in sub_model_UNIFAC.F90)
74// (none right now)
75
76/** \brief Get the Jacobian elements used for a particular row of the matrix
77 *
78 * \param sub_model_int_data Pointer to the sub model integer data
79 * \param sub_model_float_data Pointer to the sub model floating-point data
80 * \param jac Jacobian
81 */
82void sub_model_UNIFAC_get_used_jac_elem(int *sub_model_int_data,
83 double *sub_model_float_data,
84 Jacobian *jac) {
85 int *int_data = sub_model_int_data;
86 double *float_data = sub_model_float_data;
87
88 for (int i_phase = 0; i_phase < NUM_UNIQUE_PHASE_; ++i_phase)
89 for (int i_inst = 0; i_inst < NUM_PHASE_INSTANCE_(i_phase); ++i_inst)
90 for (int i_gamma = 0; i_gamma < NUM_SPEC_(i_phase); ++i_gamma)
91 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
92 jacobian_register_element(
93 jac,
94 PHASE_INST_ID_(i_phase, i_inst) + GAMMA_ID_(i_phase, i_gamma),
95 PHASE_INST_ID_(i_phase, i_inst) + SPEC_ID_(i_phase, i_spec));
96}
97
98/** \brief Update stored ids for elements used within a row of the Jacobian
99 * matrix
100 *
101 * \param sub_model_int_data Pointer to the sub model integer data
102 * \param sub_model_float_data Pointer to the sub model floating-point data
103 * \param deriv_ids Indices for state array variables on the solver state array
104 * \param jac Jacobian
105 */
106void sub_model_UNIFAC_update_ids(int *sub_model_int_data,
107 double *sub_model_float_data, int *deriv_ids,
108 Jacobian jac) {
109 int *int_data = sub_model_int_data;
110 double *float_data = sub_model_float_data;
111
112 for (int i_phase = 0; i_phase < NUM_UNIQUE_PHASE_; ++i_phase)
113 for (int i_inst = 0; i_inst < NUM_PHASE_INSTANCE_(i_phase); ++i_inst)
114 for (int i_gamma = 0; i_gamma < NUM_SPEC_(i_phase); ++i_gamma)
115 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
116 JAC_ID_(i_phase, i_inst, i_gamma, i_spec) = jacobian_get_element_id(
117 jac,
118 PHASE_INST_ID_(i_phase, i_inst) + GAMMA_ID_(i_phase, i_gamma),
119 PHASE_INST_ID_(i_phase, i_inst) + SPEC_ID_(i_phase, i_spec));
120}
121
122/** \brief Update sub-model data for new environmental conditions
123 *
124 * \param sub_model_int_data Pointer to the sub model integer data
125 * \param sub_model_float_data Pointer to the sub model floating-point data
126 * \param sub_model_env_data Pointer to the sub model environment-dependent data
127 * \param model_data Pointer to the model data
128 */
129void sub_model_UNIFAC_update_env_state(int *sub_model_int_data,
130 double *sub_model_float_data,
131 double *sub_model_env_data,
132 ModelData *model_data) {
133 int *int_data = sub_model_int_data;
134 double *float_data = sub_model_float_data;
135 double *env_data = model_data->grid_cell_env;
136
137 // Update the interaction parameters
138 for (int m = 0; m < NUM_GROUP_; m++)
139 for (int n = 0; n < NUM_GROUP_; n++)
140 PSI_MN_(m, n) = exp(-A_MN_(m, n) / TEMPERATURE_K_);
141
142 // Calculate the pure liquid residual acitivity ceofficient ln(GAMMA_k^(i))
143 // terms. Eq. 7 & 8
144 for (int i_phase = 0; i_phase < NUM_UNIQUE_PHASE_; i_phase++) {
145 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); i_spec++) {
146 // Calculate the sum (Q_n * X_n) in the denominator of Eq. 9 for the
147 // pure liquid
148 double sum_Qn_Xn = 0.0;
149 double total_group_moles = 0.0;
150 for (int i_group = 0; i_group < NUM_GROUP_; i_group++) {
151 sum_Qn_Xn += Q_K_(i_group) * V_IK_(i_phase, i_spec, i_group);
152 total_group_moles += V_IK_(i_phase, i_spec, i_group);
153 }
154 sum_Qn_Xn *= 1.0 / total_group_moles;
155
156 // Calculate THETA_m Eq. 9
157 for (int m = 0; m < NUM_GROUP_; m++)
158 THETA_M_(m) =
159 Q_K_(m) * V_IK_(i_phase, i_spec, m) / sum_Qn_Xn / total_group_moles;
160
161 // Calculate ln(GAMMA_k^(i))
162 for (int k = 0; k < NUM_GROUP_; k++) {
163 double sum_m_A = 0.0; // ln( sum_m_A ) term in Eq. 8
164 double sum_m_B = 0.0; // last term in Eq. 8
165 for (int m = 0; m < NUM_GROUP_; m++) {
166 sum_m_A += THETA_M_(m) * PSI_MN_(m, k);
167 double sum_n = 0.0;
168 for (int n = 0; n < NUM_GROUP_; n++)
169 sum_n += THETA_M_(n) * PSI_MN_(n, m);
170 sum_m_B += THETA_M_(m) * PSI_MN_(k, m) / sum_n;
171 }
172 LN_GAMMA_IK_(i_phase, i_spec, k) =
173 Q_K_(k) * (1.0 - log(sum_m_A) - sum_m_B);
174 }
175 }
176 }
177}
178
179/** \brief Perform the sub-model calculations for the current model state
180 *
181 * These calculations are based on equations (1)--(9) in \cite Marcolli2005.
182 * Variable names used in the derivation of the Jacobian equations are
183 * adopted here for consistency. Specifically, these are:
184 * \f[
185 * \begin{align*}
186 * \sigma & \equiv \displaystyle\sum_k r_k x_k, \\
187 * \tau & \equiv \displaystyle\sum_k q_k x_k, \\
188 * \mu & \equiv \displaystyle\sum_j x_j l_j, \\
189 * c_{xX} & = \frac{1}{\displaystyle\sum_j
190 * \displaystyle\sum_i v_j^{(i)} x_i}, \textrm{ and} \\
191 * \Pi & \equiv \displaystyle\sum_n Q_n X_n. \\
192 * \end{align*}
193 * \f]
194 * Using these variables and subscript \f$j\f$ for the species whose activity
195 * is being calculated for easy comparison with the Jacobian calculations,
196 * equations (1)--(9) in \cite Marcolli2005 become:
197 * \f[
198 * \begin{align*}
199 * \ln{\gamma_j} & = \ln{\gamma_j^C} + \ln{\gamma_j^R} & (1) \\
200 * \alpha_w & = \gamma_w x_w & (2) \\
201 * \ln{\gamma_j^C} & = \ln{\frac{\Phi_j}{x_j}}
202 * + \frac{z}{2}q_j\ln{\frac{\Theta_j}{\Phi_j}}
203 * + l_j
204 * - \frac{\Phi_j}{x_j}\mu & (3) \\
205 * \Phi_j & = \frac{r_j x_j}{\sigma}; \qquad
206 * \Theta_j = \frac{q_j x_j}{\tau} & (4) \\
207 * l_j & = \frac{z}{2}\left(r_j - q_j\right) - r_j + 1 & (5) \\
208 * r_j & = \displaystyle\sum_k v_k^{(j)}R_k; \qquad
209 * q_j = \displaystyle\sum_k v_k^{(j)}Q_k & (6) \\
210 * \ln{\gamma_j^R} & = \displaystyle\sum_k v_k^{(j)} \left[
211 * \ln{\Gamma_k} - \ln{\Gamma_k^{(j)}}\right] & (7) \\
212 * \ln{\Gamma_k} & = Q_k \left[ 1 - \ln{\Xi_k} -
213 * \displaystyle\sum_m\frac{\Theta_m\Phi_{km}}{\Xi_m}
214 * \right] & (8) \\
215 * \Theta_m & = \frac{Q_m X_m}{\Pi}; \qquad
216 * \Psi_{mn} = \mathrm{exp}\left[\frac{-a_{mn}}{T}\right] & (9) \\
217 * \end{align*}
218 * \f]
219 *
220 * \param sub_model_int_data Pointer to the sub model integer data
221 * \param sub_model_float_data Pointer to the sub model floating-point data
222 * \param sub_model_env_data Pointer to the sub model environment-dependent data
223 * \param model_data Pointer to the model data including the current state and
224 * environmental conditions
225 */
226void sub_model_UNIFAC_calculate(int *sub_model_int_data,
227 double *sub_model_float_data,
228 double *sub_model_env_data,
229 ModelData *model_data) {
230 int *int_data = sub_model_int_data;
231 double *float_data = sub_model_float_data;
232
233 double *state = model_data->grid_cell_state;
234
235 // Loop through each instance of each phase to calculate activity
236 for (int i_phase = 0; i_phase < NUM_UNIQUE_PHASE_; ++i_phase) {
237 for (int i_instance = 0; i_instance < NUM_PHASE_INSTANCE_(i_phase);
238 ++i_instance) {
239 // Get the total number of moles of species in this phase instance
240 double m_T = 0.0;
241 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
242 m_T += state[PHASE_INST_ID_(i_phase, i_instance) +
243 SPEC_ID_(i_phase, i_spec)] /
244 MW_I_(i_phase, i_spec);
245 }
246
247 // If there are no species present, skip the calculations
248 if (m_T <= 0.0) {
249 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
250 state[PHASE_INST_ID_(i_phase, i_instance) +
251 GAMMA_ID_(i_phase, i_spec)] = 0.0;
252 }
253 continue;
254 }
255
256 // Pre-calculate mole fractions and variables that do not depend
257 // on species j
258
259 double sigma = 0.0;
260 double tau = 0.0;
261 double mu = 0.0;
262 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
263 // Mole fractions x_i
264 double x_i = state[PHASE_INST_ID_(i_phase, i_instance) +
265 SPEC_ID_(i_phase, i_spec)] /
266 MW_I_(i_phase, i_spec) / m_T;
267 X_I_(i_phase, i_spec) = x_i;
268
269 // Combinatorial variables
270 sigma += R_I_(i_phase, i_spec) * x_i;
271 tau += Q_I_(i_phase, i_spec) * x_i;
272 mu += L_I_(i_phase, i_spec) * x_i;
273 }
274
275 // Residual variables
276 double c_xX = 0.0;
277 double Pi = 0.0;
278 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
279 X_K_(i_group) = 0.0;
280 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
281 double v_ik_x_i =
282 V_IK_(i_phase, i_spec, i_group) * X_I_(i_phase, i_spec);
283 c_xX += v_ik_x_i;
284 X_K_(i_group) += v_ik_x_i;
285 }
286 }
287 c_xX = 1.0 / c_xX;
288 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
289 X_K_(i_group) *= c_xX;
290 Pi += Q_K_(i_group) * X_K_(i_group);
291 }
292 // Group surface area fractions Theta_m (Eq. 9)
293 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group)
294 THETA_M_(i_group) = Q_K_(i_group) * X_K_(i_group) / Pi;
295
296 // Xi_m
297 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
298 XI_M_(i_group) = 0.0;
299 for (int j_group = 0; j_group < NUM_GROUP_; ++j_group) {
300 XI_M_(i_group) += THETA_M_(j_group) * PSI_MN_(j_group, i_group);
301 }
302 }
303
304 // Group residual activity coefficients ln(Gamma_k) (Eq. 8)
305 for (int k = 0; k < NUM_GROUP_; ++k) {
306 LN_GAMMA_K_(k) = 1.0 - log(XI_M_(k));
307 for (int m = 0; m < NUM_GROUP_; ++m)
308 LN_GAMMA_K_(k) -= THETA_M_(m) * PSI_MN_(k, m) / XI_M_(m);
309 LN_GAMMA_K_(k) *= Q_K_(k);
310 }
311
312 // Calculate activity coefficients for each species in the phase instance
313 for (int j = 0; j < NUM_SPEC_(i_phase); ++j) {
314 // Mole fraction for species j
315 double x_j = X_I_(i_phase, j);
316
317 // Phi_j (Eq. 4)
318 double Phi_j;
319 Phi_j = R_I_(i_phase, j) * x_j / sigma;
320
321 // Theta_j (Eq. 4)
322 double Theta_j;
323 Theta_j = Q_I_(i_phase, j) * x_j / tau;
324
325 // Calculate the combinatorial term ln(gamma_i^C) Eq. 3
326 double lngammaC_j = 0.0;
327 if (X_I_(i_phase, j) > 0.0) {
328 lngammaC_j = log(Phi_j / x_j) +
329 5.0 * Q_I_(i_phase, j) * log(Theta_j / Phi_j) +
330 L_I_(i_phase, j) - Phi_j / x_j * mu;
331 }
332
333 // Full residual term
334 double lngammaR_j = 0.0;
335 for (int k = 0; k < NUM_GROUP_; ++k) {
336 lngammaR_j += V_IK_(i_phase, j, k) *
337 (LN_GAMMA_K_(k) - LN_GAMMA_IK_(i_phase, j, k));
338 }
339
340 // Calculate gamma_i Eq. 1
341 state[PHASE_INST_ID_(i_phase, i_instance) + GAMMA_ID_(i_phase, j)] =
342 exp(lngammaC_j + lngammaR_j) * X_I_(i_phase, j);
343 }
344 }
345 }
346}
347
348/** \brief Add contributions to the Jacobian from derivates calculated using the
349 * output of this sub model
350 *
351 * This derivation starts from equations (1)--(9) in \cite Marcolli2005.
352 * The mole fraction of species \f$i\f$ is calculated as:
353 * \f[
354 * x_i = \frac{m_i}{m_T},
355 * \f]
356 * where \f$m_T = \displaystyle\sum_j m_j\f$,
357 * \f$m_i = \frac{c_i}{\mathrm{MW}_i}\f$
358 * is the number in moles of species \f$i\f$,
359 * \f$c_i\f$ is its mass concentration
360 * (\f$\mathrm{ug} \: \mathrm{m}^{-3}\f$; the state variable), and
361 * \f$\mathrm{MW}_i\f$ is its molecular weight
362 * (\f$\mathrm{ug} \: \mathrm{mol}^{-1}\f$). Thus,
363 * \f[
364 * \begin{align*}
365 * \frac{\partial x_i}{\partial c_i} & =
366 * \frac{(m_T-m_i)}{\mathrm{MW}_i m_T^2}, \textrm{ and} \\
367 * \frac{\partial x_j}{\partial c_i} & =
368 * -\frac{m_j}{\mathrm{MW}_i m_T^2} \quad \text{for } i\neq j. \\
369 * \end{align*}
370 * \f]
371 * The partial derivative of \f$\Phi_j\f$ (Eq. 4) with respect to \f$c_i\f$
372 * is derived as follows:
373 * \f[
374 * \begin{align*}
375 * \frac{\partial r_i x_i}{\partial c_i} & =
376 * r_i \frac{(m_T-m_i)}{\mathrm{MW}_i m_T^2}, \\
377 * \frac{\partial r_j x_j}{\partial c_i} & =
378 * - r_j \frac{m_j}{\mathrm{MW}_i m_T^2}
379 * \quad \text{for } i\neq j, \\
380 * \frac{\partial \displaystyle\sum_j r_j x_j}{\partial c_i} & =
381 * \frac{1}{\mathrm{MW}_i m_T}\left[r_i -
382 * \displaystyle\sum_j r_j x_j\right], \\
383 * \sigma & \equiv \displaystyle\sum_k r_k x_k, \\
384 * \frac{\partial \Phi_j}{\partial c_i} & =
385 * \frac{r_j}{\mathrm{MW}_i m_T \sigma^2}
386 * \left(\sigma^{\prime} - x_j r_i \right)
387 * , \qquad
388 * \sigma^{\prime} =
389 * \begin{cases}
390 * \sigma & \quad \text{if } i=j \\
391 * 0 & \quad \text{if } i\neq j
392 * \end{cases} \\
393 * \end{align*}
394 * \f]
395 * Similarly, the partial derivative of \f$\Theta_j\f$ (Eq. 4) with respect
396 * to \f$c_i\f$ is:
397 * \f[
398 * \begin{align*}
399 * \tau & \equiv \displaystyle\sum_k q_k x_k, \\
400 * \frac{\partial \Theta_j}{\partial c_i} & =
401 * \frac{q_j}{\mathrm{MW}_i m_T \tau^2}
402 * \left(\tau^{\prime} - x_j q_i \right)
403 * , \qquad
404 * \tau^{\prime} =
405 * \begin{cases}
406 * \tau & \quad \text{if } i=j \\
407 * 0 & \quad \text{if } i\neq j
408 * \end{cases} \\
409 * \end{align*}
410 * \f]
411 * From Eqs 5 and 6,
412 * \f[
413 * \frac{\partial r_j}{\partial c_i} =
414 * \frac{\partial q_j}{\partial c_i} =
415 * \frac{\partial l_j}{\partial c_i} = 0.
416 * \f]
417 * For the last term in Eq. 3,
418 * \f[
419 * \begin{align*}
420 * \mu & \equiv \displaystyle\sum_j x_j l_j, \\
421 * \frac{\partial \mu}{\partial c_i} & = \frac{1}{\mathrm{MW}_i m_T}
422 * \left( l_i - \mu \right). \\
423 * \end{align*}
424 * \f]
425 * The partial derivative of the full combinatorial term (Eq. 3) is:
426 * \f[
427 * \begin{align*}
428 * \frac{\partial \ln{\gamma^C_j}}{\partial c_i} & =
429 * \frac{x_j}{\Phi_j}\frac{\left( \frac{\partial \Phi_j}{\partial c_i}
430 * x_j - \Phi_j\frac{\partial x_j}{\partial c_i} \right)}{x_j^2} \\
431 * & \quad + \frac{z}{2}q_j\frac{\Phi_j}{\Theta_j}\frac{\left(
432 * \frac{\partial \Theta_j}{\partial c_i} \Phi_j -
433 * \Theta_j\frac{\partial \Phi_j}{\partial c_i}\right)}{\Phi_j^2} \\
434 * & \quad - \frac{\left[ \left( \frac{\partial \Phi_j}{\partial c_i} \mu +
435 * \Phi_j \frac{\partial \mu}{\partial c_i} \right) x_j
436 * - \Phi_j \mu \frac{\partial x_j}{\partial c_i} \right]}{x_j^2}.
437 * \end{align*}
438 * \f]
439 * After some rearranging, this becomes:
440 * \f[
441 * \begin{align*}
442 * \frac{\partial \ln{\gamma^C_j}}{\partial c_i} & =
443 * \frac{1}{\Phi_j}\frac{\partial \Phi_j}{\partial c_i}
444 * - \frac{1}{x_j}\frac{\partial x_j}{\partial c_i} \\
445 * & \quad + \frac{z}{2}q_j\left(
446 * \frac{1}{\Theta_j}\frac{\partial \Theta_j}{\partial c_i}
447 * - \frac{1}{\Phi_j}\frac{\partial \Phi_j}{\partial c_i}\right) \\
448 * & \quad - \frac{\partial \Phi_j}{\partial c_i}\frac{\mu}{x_j}
449 * - \frac{\Phi_j}{x_j}\frac{\partial \mu}{\partial c_i}
450 * + \frac{\Phi_j \mu}{x_j^2}\frac{\partial x_j}{\partial c_i}
451 * \end{align*}
452 * \f]
453 * As \f$\sigma\f$, \f$\tau\f$,
454 * and \f$\mu\f$ are independent of species \f$i\f$,
455 * these can be calculated outside the loop over the independent species.
456 * Moving to the residual term (Eq. 7), the mole fraction (\f$X_p\f$)
457 * of group \f$p\f$ in the mixture is related to the mole fraction
458 * (\f$x_j\f$) of species \f$j\f$ according to:
459 * \f[
460 * X_p = c_{xX} \omega_p, \textrm{ where }
461 * \omega_p \equiv \displaystyle\sum_j v_p^{(j)} x_j, \\
462 * \f]
463 * and \f$c_{xX}\f$ is a conversion factor accounting for the difference
464 * in total species and total group number concentrations:
465 * \f[
466 * c_{xX} = \frac{1}{\displaystyle\sum_p \omega_p}.
467 * \f]
468 * Partial derivatives of the group interaction terms in Eq 9,
469 * \f$\Theta_m\f$ and \f$\Psi_{mn}\f$, with respect to \f$c_i\f$ are
470 * derived as follows:
471 * \f[
472 * \begin{align*}
473 * \frac{\partial c_{xX}}{\partial c_i} & = - c_{xX}^2
474 * \displaystyle\sum_p \displaystyle\sum_j v_p^{(j)}
475 * \frac{\partial x_j}{\partial c_i}
476 * = - \frac{c_{xX}^2}{\mathrm{MW}_im_T} \left(
477 * \displaystyle\sum_p v_p^{(i)} -
478 * \displaystyle\sum_p \omega_p \right), \\
479 * & = \frac{c_{xX}}{\mathrm{MW}_im_T} \left(
480 * 1 - c_{xX} \displaystyle\sum_p v_p^{(i)} \right), \\
481 * \frac{\partial X_n}{\partial c_i} & =
482 * c_{xX} \displaystyle\sum_j v_n^{(j)}
483 * \frac{\partial x_j}{\partial c_i}
484 * + \frac{\partial c_{xX}}{\partial c_i} \omega_n, \\
485 * & = \frac{c_{xX}}{\mathrm{MW}_im_T}
486 * \left( v_n^{(i)} - \omega_n \right)
487 * + \frac{c_{xX}}{\mathrm{MW}_im_T}
488 * \left( \omega_n - c_{xX}\omega_n
489 * \displaystyle\sum_p v_p^{(i)} \right), \\
490 * & = \frac{c_{xX}}{\mathrm{MW}_im_T} \left(
491 * v_n^{(i)} - X_n
492 * \displaystyle\sum_p v_p^{(i)} \right), \\
493 * \Pi & \equiv \displaystyle\sum_n Q_n X_n, \\
494 * \frac{\partial \Pi}{\partial c_i} & =
495 * \displaystyle\sum_n Q_n
496 * \frac{c_{xX}}{\mathrm{MW}_i m_T}\left(
497 * v_n^{(i)} - X_n \displaystyle\sum_p v_p^{(i)} \right), \\
498 * & = \frac{c_{xX}}{\mathrm{MW}_i m_T} \left(
499 * \displaystyle\sum_n Q_n v_n^{(i)} -
500 * \Pi \displaystyle\sum_p v_p^{(i)} \right), \\
501 * \frac{\partial \Theta_m}{\partial c_i} & = \frac{
502 * \left(Q_m\frac{\partial X_m}{\partial c_i} \Pi -
503 * Q_m X_m \frac{\partial \Pi}{\partial c_i} \right)}{\Pi^2} , \\
504 * & = \frac{c_{xX} Q_m}{\mathrm{MW}_im_T\Pi^2} \left[
505 * \left( \Pi v_m^{(i)} - \Pi X_m
506 * \displaystyle\sum_p v_p^{(i)} \right) -
507 * \left( X_m \displaystyle\sum_n Q_n v_n^{(i)} -
508 * X_m \Pi \displaystyle\sum_p v_p^{(i)} \right) \right], \\
509 * & = \frac{c_{xX} Q_m}{\mathrm{MW}_im_T\Pi^2} \left[
510 * \Pi v_m^{(i)} - X_m\displaystyle\sum_n Q_n v_n^{(i)}\right],
511 * \textrm{ and} \\
512 * \frac{\partial \Psi_{mn}}{\partial c_i} & = 0
513 * \end{align*}
514 * \f]
515 * The partial derivative of the group residual activity coefficient (Eq. 8)
516 * with respect to \f$c_i\f$ is:
517 * \f[
518 * \begin{align*}
519 * \frac{\partial \ln{\Gamma_k}}{\partial c_i} & = Q_k \left[
520 * - \frac{1}{\displaystyle\sum_m\Theta_m\Psi_{mk}}
521 * \displaystyle\sum_m\frac{\partial\Theta_m}{\partial c_i}\Psi_{mk} \\
522 * \quad - \displaystyle\sum_m\left(
523 * \frac{\partial\Theta_m}{\partial c_i}\Psi_{km}
524 * \displaystyle\sum_n\Theta_n\Psi_{nm}
525 * - \Theta_m\Psi_{km}\displaystyle\sum_n
526 * \frac{\partial\Theta_n}{\partial c_i}
527 * \Psi_{nm}\right)/\left(
528 * \displaystyle\sum_n\Theta_n\Psi_{nm}\right)^2
529 * \right]
530 * \end{align*}
531 * \f]
532 * After some rearranging, this becomes:
533 * \f[
534 * \begin{align*}
535 * \Xi_m & \equiv \displaystyle\sum_n\Theta_n\Psi_{nm}, \\
536 * \frac{\partial \ln{\Gamma_k}}{\partial c_i} & =
537 * - Q_k \displaystyle\sum_m \left(
538 * \frac{\Psi_{mk}}{\Xi_k}\frac{\partial\Theta_m}{\partial c_i}
539 * + \frac{\Psi_{km}}{\Xi_m}
540 * \frac{\partial\Theta_m}{\partial c_i}
541 * - \frac{\Theta_m\Psi_{km}}
542 * {\Xi_m^2}
543 * \displaystyle\sum_n\frac{\partial\Theta_n}{\partial c_i}\Psi_{nm}
544 * \right)
545 * \end{align*}
546 * \f]
547 * The left side of the three bracketed terms in the above equation are
548 * independent of species \f$i\f$ and can be calculated outside of the loop
549 * over the independent species. The partial derivative of the full residual
550 * term with respect to \f$c_i\f$ is:
551 * \f[
552 * \begin{align*}
553 * \frac{\partial \ln{\Gamma_k^{(j)}}}{\partial c_i} & = 0 \\
554 * \frac{\partial \ln{\gamma_j^R}}{\partial c_i} & =
555 * \displaystyle\sum_k v_k^{(j)}
556 * \frac{\partial \ln{\Gamma_k}}{\partial c_i}
557 * \end{align*}
558 * \f]
559 * The overall equation for the partial derivative of \f$\gamma_j\f$ with
560 * respect to species \f$i\f$ is:
561 * \f[
562 * \frac{\partial\gamma_j}{\partial c_i} = \gamma_j \left(
563 * \frac{\partial\ln{\gamma_j^C}}{\partial c_i}
564 * + \frac{\partial\ln{\gamma_j^R}}{\partial c_i} \right)
565 * \f]
566 *
567 * \param sub_model_int_data Pointer to the sub model integer data
568 * \param sub_model_float_data Pointer to the sub model floating-point data
569 * \param sub_model_env_data Pointer to the sub model environment-dependent data
570 * \param model_data Pointer to the model data
571 * \param J Jacobian to be calculated
572 * \param time_step Current time step [s]
573 */
574#ifdef CAMP_USE_SUNDIALS
575void sub_model_UNIFAC_get_jac_contrib(int *sub_model_int_data,
576 double *sub_model_float_data,
577 double *sub_model_env_data,
578 ModelData *model_data, realtype *J,
579 double time_step) {
580 int *int_data = sub_model_int_data;
581 double *float_data = sub_model_float_data;
582 double *state = model_data->grid_cell_state;
583 double *env_data = model_data->grid_cell_env;
584
585 // Loop through each instance of each phase to calculate activity
586 for (int i_phase = 0; i_phase < NUM_UNIQUE_PHASE_; ++i_phase) {
587 for (int i_instance = 0; i_instance < NUM_PHASE_INSTANCE_(i_phase);
588 ++i_instance) {
589 // Get the total number of moles of species in this phase instance
590 // (in umoles)
591 double m_T = 0.0;
592 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
593 m_T += state[PHASE_INST_ID_(i_phase, i_instance) +
594 SPEC_ID_(i_phase, i_spec)] /
595 MW_I_(i_phase, i_spec);
596 }
597
598 // If there are no species present, skip the calculations
599 if (m_T <= 0.0) {
600 continue;
601 }
602
603 // Pre-calculate mole fractions and variables that do not depend on
604 // dependent species j or independent species i
605
606 double sigma = 0.0;
607 double tau = 0.0;
608 double mu = 0.0;
609 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
610 // Mole fractions x_i for this phase instance
611 double x_i = state[PHASE_INST_ID_(i_phase, i_instance) +
612 SPEC_ID_(i_phase, i_spec)] /
613 MW_I_(i_phase, i_spec) / m_T;
614 X_I_(i_phase, i_spec) = x_i;
615
616 // Combinatorial variables
617 sigma += R_I_(i_phase, i_spec) * x_i;
618 tau += Q_I_(i_phase, i_spec) * x_i;
619 mu += L_I_(i_phase, i_spec) * x_i;
620 }
621
622 // Residual variables
623 double c_xX = 0.0;
624 double Pi = 0.0;
625 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
626 X_K_(i_group) = 0.0;
627 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec) {
628 double v_ik_x_i =
629 V_IK_(i_phase, i_spec, i_group) * X_I_(i_phase, i_spec);
630 c_xX += v_ik_x_i;
631 X_K_(i_group) += v_ik_x_i;
632 }
633 }
634 c_xX = 1.0 / c_xX;
635 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
636 X_K_(i_group) *= c_xX;
637 Pi += Q_K_(i_group) * X_K_(i_group);
638 }
639
640 // Group surface area fractions Theta_m (Eq. 9)
641 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group)
642 THETA_M_(i_group) = Q_K_(i_group) * X_K_(i_group) / Pi;
643
644 // Xi_m
645 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
646 XI_M_(i_group) = 0;
647 for (int j_group = 0; j_group < NUM_GROUP_; ++j_group)
648 XI_M_(i_group) += THETA_M_(j_group) * PSI_MN_(j_group, i_group);
649 }
650
651 // Group residual activity coefficients ln(Gamma_k) (Eq. 8)
652 for (int k = 0; k < NUM_GROUP_; ++k) {
653 LN_GAMMA_K_(k) = 1.0 - log(XI_M_(k));
654 for (int m = 0; m < NUM_GROUP_; ++m)
655 LN_GAMMA_K_(k) -= THETA_M_(m) * PSI_MN_(k, m) / XI_M_(m);
656 LN_GAMMA_K_(k) *= Q_K_(k);
657 }
658
659 // Loop over dependent species j
660 for (int j = 0; j < NUM_SPEC_(i_phase); ++j) {
661 // Mole fraction for species j
662 double x_j = X_I_(i_phase, j);
663
664 // Phi_j (Eq. 4)
665 double Phi_j;
666 Phi_j = R_I_(i_phase, j) * x_j / sigma;
667
668 // Theta_j (Eq. 4)
669 double Theta_j;
670 Theta_j = Q_I_(i_phase, j) * x_j / tau;
671
672 // Full combinatorial term
673 double lngammaC_j = 0.0;
674 if (X_I_(i_phase, j) > 0.0) {
675 lngammaC_j = log(Phi_j / x_j) +
676 5.0 * Q_I_(i_phase, j) * log(Theta_j / Phi_j) +
677 L_I_(i_phase, j) - Phi_j / x_j * mu;
678 }
679
680 // Full residual term
681 double lngammaR_j = 0.0;
682 for (int k = 0; k < NUM_GROUP_; ++k) {
683 lngammaR_j += V_IK_(i_phase, j, k) *
684 (LN_GAMMA_K_(k) - LN_GAMMA_IK_(i_phase, j, k));
685 }
686
687 // Activity coefficient gamma_j (in m3/ug)
688 double gamma_j = exp(lngammaC_j + lngammaR_j);
689
690 // Loop over independent species i
691 double dx_j_dc_i, dPhi_j_dc_i, dTheta_j_dc_i; // combinatorial
692 double dmu_dc_i, dlngammaC_j_dc_i; // terms
693 for (int i = 0; i < NUM_SPEC_(i_phase); ++i) {
694 // combinatorial partial derivatives
695 if (i == j) {
696 dx_j_dc_i = 1.0 - x_j;
697 dPhi_j_dc_i = sigma - x_j * R_I_(i_phase, i);
698 dTheta_j_dc_i = tau - x_j * Q_I_(i_phase, i);
699 } else {
700 dx_j_dc_i = -x_j;
701 dPhi_j_dc_i = -x_j * R_I_(i_phase, i);
702 dTheta_j_dc_i = -x_j * Q_I_(i_phase, i);
703 }
704 double MWimT_inv = 1.0 / (MW_I_(i_phase, i) * m_T);
705 dx_j_dc_i *= MWimT_inv;
706 dPhi_j_dc_i *= R_I_(i_phase, j) * MWimT_inv / sigma / sigma;
707 dTheta_j_dc_i *= Q_I_(i_phase, j) * MWimT_inv / tau / tau;
708
709 dmu_dc_i = MWimT_inv * (L_I_(i_phase, i) - mu);
710
711 // partial derivative of full combinatorial term ln(gammaC_j)
712 dlngammaC_j_dc_i =
713 dPhi_j_dc_i / Phi_j - dx_j_dc_i / x_j +
714 5.0 * Q_I_(i_phase, j) *
715 (dTheta_j_dc_i / Theta_j - dPhi_j_dc_i / Phi_j) -
716 dPhi_j_dc_i * mu / x_j - Phi_j / x_j * dmu_dc_i +
717 Phi_j * mu / x_j / x_j * dx_j_dc_i;
718
719 // Residual partial derivatives
720
721 // partial derivative of group surface area fractions Theta_m
722 double sumQnvn_i = 0.0;
723 for (int n = 0; n < NUM_GROUP_; ++n)
724 sumQnvn_i += Q_K_(n) * V_IK_(i_phase, i, n);
725 for (int m = 0; m < NUM_GROUP_; ++m) {
726 DTHETA_M_DC_I_(m) =
727 Q_K_(m) * c_xX * MWimT_inv / Pi / Pi *
728 (Pi * V_IK_(i_phase, i, m) - X_K_(m) * sumQnvn_i);
729 }
730
731 // partial derivative of full residual terms ln(gammaR_j)
732 double dlngammaR_j_dc_i = 0.0;
733 for (int k = 0; k < NUM_GROUP_; ++k) {
734 // partial derivative of group residual activity coefficients
735 // ln(Gamma_k)
736 double dlnGamma_k_dc_i = 0.0;
737 for (int m = 0; m < NUM_GROUP_; ++m) {
738 double sum_n = 0.0;
739 for (int n = 0; n < NUM_GROUP_; ++n)
740 sum_n += DTHETA_M_DC_I_(n) * PSI_MN_(n, m);
741 dlnGamma_k_dc_i +=
742 DTHETA_M_DC_I_(m) *
743 (PSI_MN_(m, k) / XI_M_(k) + PSI_MN_(k, m) / XI_M_(m)) -
744 THETA_M_(m) * PSI_MN_(k, m) / pow(XI_M_(m), 2) * sum_n;
745 }
746 dlnGamma_k_dc_i *= -Q_K_(k);
747
748 dlngammaR_j_dc_i += V_IK_(i_phase, j, k) * dlnGamma_k_dc_i;
749 }
750
751 // partial derivative of activity coefficient gamma_j
752 double dgamma_j_dc_i =
753 gamma_j * x_j * (dlngammaC_j_dc_i + dlngammaR_j_dc_i) +
754 gamma_j * dx_j_dc_i;
755
756 // Add the partial derivative contribution to the
757 // Jacobian
758 J[JAC_ID_(i_phase, i_instance, j, i)] += dgamma_j_dc_i;
759 }
760 }
761 }
762 }
763}
764#endif
765
766/** \brief Print the sub model data
767 *
768 * \param sub_model_int_data Pointer to the sub model integer data
769 * \param sub_model_float_data Pointer to the sub model floating-point data
770 */
771void sub_model_UNIFAC_print(int *sub_model_int_data,
772 double *sub_model_float_data) {
773 int *int_data = sub_model_int_data;
774 double *float_data = sub_model_float_data;
775
776 printf("\n\nUNIFAC sub model\n\n");
777 printf("\nint_data");
778 for (int i = 0; i < INT_DATA_SIZE_; i++) printf(" %d", int_data[i]);
779 printf("\nfloat_data");
780 for (int i = 0; i < FLOAT_DATA_SIZE_; i++) printf(" %le", float_data[i]);
781 printf("\nNumber of unique phases %d", NUM_UNIQUE_PHASE_);
782 printf("\nNumber of UNIFAC groups %d", NUM_GROUP_);
783 printf("\n*** General group data ***");
784 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
785 printf("\nGroup %d:", i_group);
786 printf("\n Q_K: %le R_K: %le X_k: %le", Q_K_(i_group), R_K_(i_group),
787 X_K_(i_group));
788 printf("\n Xi_M: %le ln(Gamma_k): %le", XI_M_(i_group),
789 LN_GAMMA_K_(i_group));
790 printf("\n dTheta_n / dc_i): %le", DTHETA_M_DC_I_(i_group));
791 printf("\n A_MN (by group):");
792 for (int j_group = 0; j_group < NUM_GROUP_; ++j_group)
793 printf(" %le", A_MN_(i_group, j_group));
794 }
795 printf("\n\n*** Phase data ***");
796 for (int i_phase = 0; i_phase < NUM_UNIQUE_PHASE_; ++i_phase) {
797 printf("\nPhase %d:", i_phase);
798 printf("\n Number of instances: %d", NUM_PHASE_INSTANCE_(i_phase));
799 printf("\n Number of species: %d", NUM_SPEC_(i_phase));
800 printf("\n Species ids:");
801 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
802 printf(" %d", SPEC_ID_(i_phase, i_spec));
803 printf("\n ** Instance data **");
804 for (int i_inst = 0; i_inst < NUM_PHASE_INSTANCE_(i_phase); ++i_inst) {
805 printf("\n Instance %d:", i_inst);
806 printf("\n State id: %d", PHASE_INST_ID_(i_phase, i_inst));
807 printf("\n Gamma ids:");
808 for (int i_gamma = 0; i_gamma < NUM_SPEC_(i_phase); ++i_gamma)
809 printf(" %d",
810 PHASE_INST_ID_(i_phase, i_inst) + GAMMA_ID_(i_phase, i_gamma));
811 printf("\n Species ids:");
812 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
813 printf(" %d",
814 PHASE_INST_ID_(i_phase, i_inst) + SPEC_ID_(i_phase, i_spec));
815 printf("\n Jacobian ids:");
816 for (int i_gamma = 0; i_gamma < NUM_SPEC_(i_phase); ++i_gamma)
817 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
818 printf(" J[%d][%d] = %d", i_gamma, i_spec,
819 JAC_ID_(i_phase, i_inst, i_gamma, i_spec));
820 }
821 printf("\n R_I (by species):");
822 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
823 printf(" %le", R_I_(i_phase, i_spec));
824 printf("\n Q_I (by species):");
825 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
826 printf(" %le", Q_I_(i_phase, i_spec));
827 printf("\n L_I (by species):");
828 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
829 printf(" %le", L_I_(i_phase, i_spec));
830 printf("\n MW_I (by species):");
831 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
832 printf(" %le", MW_I_(i_phase, i_spec));
833 printf("\n X_I (by species):");
834 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
835 printf(" %le", X_I_(i_phase, i_spec));
836 printf("\n ** Phase-specific group data **");
837 for (int i_group = 0; i_group < NUM_GROUP_; ++i_group) {
838 printf("\n Group %d:", i_group);
839 printf("\n v_ik (by species):");
840 for (int i_spec = 0; i_spec < NUM_SPEC_(i_phase); ++i_spec)
841 printf(" %d", V_IK_(i_phase, i_spec, i_group));
842 }
843 }
844}
jacobian_get_element_id
unsigned int jacobian_get_element_id(Jacobian jac, unsigned int dep_id, unsigned int ind_id)
Get an element id in the Jacobian data arrays.
Definition Jacobian.c:200
jacobian_register_element
void jacobian_register_element(Jacobian *jac, unsigned int dep_id, unsigned int ind_id)
Adds an element to the sparse matrix.
Definition Jacobian.c:105
ModelData::grid_cell_env
double * grid_cell_env
Definition camp_common.h:104
ModelData::grid_cell_state
double * grid_cell_state
Definition camp_common.h:101
sub_model_UNIFAC_get_jac_contrib
void sub_model_UNIFAC_get_jac_contrib(int *sub_model_int_data, double *sub_model_float_data, double *sub_model_env_data, ModelData *model_data, realtype *J, double time_step)
Add contributions to the Jacobian from derivates calculated using the output of this sub model.
Definition sub_model_UNIFAC.c:575
X_K_
#define X_K_(m)
Definition sub_model_UNIFAC.c:54
MW_I_
#define MW_I_(p, i)
Definition sub_model_UNIFAC.c:62
sub_model_UNIFAC_update_ids
void sub_model_UNIFAC_update_ids(int *sub_model_int_data, double *sub_model_float_data, int *deriv_ids, Jacobian jac)
Update stored ids for elements used within a row of the Jacobian matrix.
Definition sub_model_UNIFAC.c:106
SPEC_ID_
#define SPEC_ID_(p, i)
Definition sub_model_UNIFAC.c:40
R_I_
#define R_I_(p, i)
Definition sub_model_UNIFAC.c:59
LN_GAMMA_IK_
#define LN_GAMMA_IK_(p, i, k)
Definition sub_model_UNIFAC.c:67
sub_model_UNIFAC_update_env_state
void sub_model_UNIFAC_update_env_state(int *sub_model_int_data, double *sub_model_float_data, double *sub_model_env_data, ModelData *model_data)
Update sub-model data for new environmental conditions.
Definition sub_model_UNIFAC.c:129
FLOAT_DATA_SIZE_
#define FLOAT_DATA_SIZE_
Definition sub_model_UNIFAC.c:71
DTHETA_M_DC_I_
#define DTHETA_M_DC_I_(m)
Definition sub_model_UNIFAC.c:55
XI_M_
#define XI_M_(m)
Definition sub_model_UNIFAC.c:56
X_I_
#define X_I_(p, i)
Definition sub_model_UNIFAC.c:63
PHASE_INST_ID_
#define PHASE_INST_ID_(p, c)
Definition sub_model_UNIFAC.c:39
Q_I_
#define Q_I_(p, i)
Definition sub_model_UNIFAC.c:60
sub_model_UNIFAC_get_used_jac_elem
void sub_model_UNIFAC_get_used_jac_elem(int *sub_model_int_data, double *sub_model_float_data, Jacobian *jac)
Get the Jacobian elements used for a particular row of the matrix.
Definition sub_model_UNIFAC.c:82
THETA_M_
#define THETA_M_(m)
Definition sub_model_UNIFAC.c:65
sub_model_UNIFAC_calculate
void sub_model_UNIFAC_calculate(int *sub_model_int_data, double *sub_model_float_data, double *sub_model_env_data, ModelData *model_data)
Perform the sub-model calculations for the current model state.
Definition sub_model_UNIFAC.c:226
NUM_PHASE_INSTANCE_
#define NUM_PHASE_INSTANCE_(p)
Definition sub_model_UNIFAC.c:37
NUM_SPEC_
#define NUM_SPEC_(p)
Definition sub_model_UNIFAC.c:38
INT_DATA_SIZE_
#define INT_DATA_SIZE_
Definition sub_model_UNIFAC.c:70
GAMMA_ID_
#define GAMMA_ID_(p, i)
Definition sub_model_UNIFAC.c:42
TEMPERATURE_K_
#define TEMPERATURE_K_
Definition sub_model_UNIFAC.c:23
NUM_UNIQUE_PHASE_
#define NUM_UNIQUE_PHASE_
Definition sub_model_UNIFAC.c:26
Q_K_
#define Q_K_(k)
Definition sub_model_UNIFAC.c:52
R_K_
#define R_K_(k)
Definition sub_model_UNIFAC.c:53
sub_model_UNIFAC_print
void sub_model_UNIFAC_print(int *sub_model_int_data, double *sub_model_float_data)
Print the sub model data.
Definition sub_model_UNIFAC.c:771
PSI_MN_
#define PSI_MN_(m, n)
Definition sub_model_UNIFAC.c:66
LN_GAMMA_K_
#define LN_GAMMA_K_(m)
Definition sub_model_UNIFAC.c:57
JAC_ID_
#define JAC_ID_(p, c, j, i)
Definition sub_model_UNIFAC.c:44
NUM_GROUP_
#define NUM_GROUP_
Definition sub_model_UNIFAC.c:27
A_MN_
#define A_MN_(m, n)
Definition sub_model_UNIFAC.c:58
V_IK_
#define V_IK_(p, i, k)
Definition sub_model_UNIFAC.c:47
L_I_
#define L_I_(p, i)
Definition sub_model_UNIFAC.c:61
CAMP 1.0.0 documentation generated by doxygen 1.9.8