CAMP 1.0.0
Chemistry Across Multiple Phases
|
Header file for sub model functions. More...
Go to the source code of this file.
Functions | |
void | sub_model_PDFiTE_get_used_jac_elem (int *sub_model_int_data, double *sub_model_float_data, Jacobian *jac) |
Flag Jacobian elements used by this sub model. | |
void | sub_model_PDFiTE_update_ids (int *sub_model_int_data, double *sub_model_float_data, int *deriv_ids, Jacobian jac) |
Update the time derivative and Jacbobian array indices. | |
void | sub_model_PDFiTE_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. | |
void | sub_model_PDFiTE_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. | |
void | sub_model_PDFiTE_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. | |
void | sub_model_PDFiTE_print (int *sub_model_int_data, double *sub_model_float_data) |
Print the PDFiTE Activity sub model parameters. | |
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. | |
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. | |
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. | |
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. | |
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. | |
void | sub_model_UNIFAC_print (int *sub_model_int_data, double *sub_model_float_data) |
Print the sub model data. | |
void | sub_model_ZSR_aerosol_water_get_used_jac_elem (int *sub_model_int_data, double *sub_model_float_data, Jacobian *jac) |
Flag Jacobian elements used by this sub model. | |
void | sub_model_ZSR_aerosol_water_update_ids (int *sub_model_int_data, double *sub_model_float_data, int *deriv_ids, Jacobian jac) |
Update the time derivative and Jacbobian array indices. | |
void | sub_model_ZSR_aerosol_water_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. | |
void | sub_model_ZSR_aerosol_water_calculate (int *sub_model_int_data, double *sub_model_float_data, double *sub_model_env_data, ModelData *model_data) |
Do pre-derivative calculations. | |
void | sub_model_ZSR_aerosol_water_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. | |
void | sub_model_ZSR_aerosol_water_print (int *sub_model_int_data, double *sub_model_float_data) |
Print the ZSR Aerosol Water sub model parameters. | |
Header file for sub model functions.
Definition in file sub_models.h.
void sub_model_PDFiTE_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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data including the current state and environmental conditions |
Definition at line 167 of file sub_model_PDFiTE.c.
void sub_model_PDFiTE_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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data |
J | Jacobian to be calculated |
time_step | Current time step in [s] |
Definition at line 277 of file sub_model_PDFiTE.c.
void sub_model_PDFiTE_get_used_jac_elem | ( | int * | sub_model_int_data, |
double * | sub_model_float_data, | ||
Jacobian * | jac | ||
) |
Flag Jacobian elements used by this sub model.
activity sub models are assumed to be at equilibrium
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
jac | Jacobian |
Definition at line 73 of file sub_model_PDFiTE.c.
void sub_model_PDFiTE_print | ( | int * | sub_model_int_data, |
double * | sub_model_float_data | ||
) |
Print the PDFiTE Activity sub model parameters.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
Definition at line 462 of file sub_model_PDFiTE.c.
void sub_model_PDFiTE_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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data |
Definition at line 137 of file sub_model_PDFiTE.c.
void sub_model_PDFiTE_update_ids | ( | int * | sub_model_int_data, |
double * | sub_model_float_data, | ||
int * | deriv_ids, | ||
Jacobian | jac | ||
) |
Update the time derivative and Jacbobian array indices.
activity sub models are assumed to be at equilibrium
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
deriv_ids | Indices for state array variables on the solver state array |
jac | Jacobian |
Definition at line 105 of file sub_model_PDFiTE.c.
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.
These calculations are based on equations (1)–(9) in [Marcolli2005]. Variable names used in the derivation of the Jacobian equations are adopted here for consistency. Specifically, these are:
\[ \begin{align*} \sigma & \equiv \displaystyle\sum_k r_k x_k, \\ \tau & \equiv \displaystyle\sum_k q_k x_k, \\ \mu & \equiv \displaystyle\sum_j x_j l_j, \\ c_{xX} & = \frac{1}{\displaystyle\sum_j \displaystyle\sum_i v_j^{(i)} x_i}, \textrm{ and} \\ \Pi & \equiv \displaystyle\sum_n Q_n X_n. \\ \end{align*} \]
Using these variables and subscript \(j\) for the species whose activity is being calculated for easy comparison with the Jacobian calculations, equations (1)–(9) in [Marcolli2005] become:
\[ \begin{align*} \ln{\gamma_j} & = \ln{\gamma_j^C} + \ln{\gamma_j^R} & (1) \\ \alpha_w & = \gamma_w x_w & (2) \\ \ln{\gamma_j^C} & = \ln{\frac{\Phi_j}{x_j}} + \frac{z}{2}q_j\ln{\frac{\Theta_j}{\Phi_j}} + l_j - \frac{\Phi_j}{x_j}\mu & (3) \\ \Phi_j & = \frac{r_j x_j}{\sigma}; \qquad \Theta_j = \frac{q_j x_j}{\tau} & (4) \\ l_j & = \frac{z}{2}\left(r_j - q_j\right) - r_j + 1 & (5) \\ r_j & = \displaystyle\sum_k v_k^{(j)}R_k; \qquad q_j = \displaystyle\sum_k v_k^{(j)}Q_k & (6) \\ \ln{\gamma_j^R} & = \displaystyle\sum_k v_k^{(j)} \left[ \ln{\Gamma_k} - \ln{\Gamma_k^{(j)}}\right] & (7) \\ \ln{\Gamma_k} & = Q_k \left[ 1 - \ln{\Xi_k} - \displaystyle\sum_m\frac{\Theta_m\Phi_{km}}{\Xi_m} \right] & (8) \\ \Theta_m & = \frac{Q_m X_m}{\Pi}; \qquad \Psi_{mn} = \mathrm{exp}\left[\frac{-a_{mn}}{T}\right] & (9) \\ \end{align*} \]
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data including the current state and environmental conditions |
Definition at line 226 of file sub_model_UNIFAC.c.
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.
This derivation starts from equations (1)–(9) in [Marcolli2005]. The mole fraction of species \(i\) is calculated as:
\[ x_i = \frac{m_i}{m_T}, \]
where \(m_T = \displaystyle\sum_j m_j\), \(m_i = \frac{c_i}{\mathrm{MW}_i}\) is the number in moles of species \(i\), \(c_i\) is its mass concentration ( \(\mathrm{ug} \: \mathrm{m}^{-3}\); the state variable), and \(\mathrm{MW}_i\) is its molecular weight ( \(\mathrm{ug} \: \mathrm{mol}^{-1}\)). Thus,
\[ \begin{align*} \frac{\partial x_i}{\partial c_i} & = \frac{(m_T-m_i)}{\mathrm{MW}_i m_T^2}, \textrm{ and} \\ \frac{\partial x_j}{\partial c_i} & = -\frac{m_j}{\mathrm{MW}_i m_T^2} \quad \text{for } i\neq j. \\ \end{align*} \]
The partial derivative of \(\Phi_j\) (Eq. 4) with respect to \(c_i\) is derived as follows:
\[ \begin{align*} \frac{\partial r_i x_i}{\partial c_i} & = r_i \frac{(m_T-m_i)}{\mathrm{MW}_i m_T^2}, \\ \frac{\partial r_j x_j}{\partial c_i} & = - r_j \frac{m_j}{\mathrm{MW}_i m_T^2} \quad \text{for } i\neq j, \\ \frac{\partial \displaystyle\sum_j r_j x_j}{\partial c_i} & = \frac{1}{\mathrm{MW}_i m_T}\left[r_i - \displaystyle\sum_j r_j x_j\right], \\ \sigma & \equiv \displaystyle\sum_k r_k x_k, \\ \frac{\partial \Phi_j}{\partial c_i} & = \frac{r_j}{\mathrm{MW}_i m_T \sigma^2} \left(\sigma^{\prime} - x_j r_i \right) , \qquad \sigma^{\prime} = \begin{cases} \sigma & \quad \text{if } i=j \\ 0 & \quad \text{if } i\neq j \end{cases} \\ \end{align*} \]
Similarly, the partial derivative of \(\Theta_j\) (Eq. 4) with respect to \(c_i\) is:
\[ \begin{align*} \tau & \equiv \displaystyle\sum_k q_k x_k, \\ \frac{\partial \Theta_j}{\partial c_i} & = \frac{q_j}{\mathrm{MW}_i m_T \tau^2} \left(\tau^{\prime} - x_j q_i \right) , \qquad \tau^{\prime} = \begin{cases} \tau & \quad \text{if } i=j \\ 0 & \quad \text{if } i\neq j \end{cases} \\ \end{align*} \]
From Eqs 5 and 6,
\[ \frac{\partial r_j}{\partial c_i} = \frac{\partial q_j}{\partial c_i} = \frac{\partial l_j}{\partial c_i} = 0. \]
For the last term in Eq. 3,
\[ \begin{align*} \mu & \equiv \displaystyle\sum_j x_j l_j, \\ \frac{\partial \mu}{\partial c_i} & = \frac{1}{\mathrm{MW}_i m_T} \left( l_i - \mu \right). \\ \end{align*} \]
The partial derivative of the full combinatorial term (Eq. 3) is:
\[ \begin{align*} \frac{\partial \ln{\gamma^C_j}}{\partial c_i} & = \frac{x_j}{\Phi_j}\frac{\left( \frac{\partial \Phi_j}{\partial c_i} x_j - \Phi_j\frac{\partial x_j}{\partial c_i} \right)}{x_j^2} \\ & \quad + \frac{z}{2}q_j\frac{\Phi_j}{\Theta_j}\frac{\left( \frac{\partial \Theta_j}{\partial c_i} \Phi_j - \Theta_j\frac{\partial \Phi_j}{\partial c_i}\right)}{\Phi_j^2} \\ & \quad - \frac{\left[ \left( \frac{\partial \Phi_j}{\partial c_i} \mu + \Phi_j \frac{\partial \mu}{\partial c_i} \right) x_j - \Phi_j \mu \frac{\partial x_j}{\partial c_i} \right]}{x_j^2}. \end{align*} \]
After some rearranging, this becomes:
\[ \begin{align*} \frac{\partial \ln{\gamma^C_j}}{\partial c_i} & = \frac{1}{\Phi_j}\frac{\partial \Phi_j}{\partial c_i} - \frac{1}{x_j}\frac{\partial x_j}{\partial c_i} \\ & \quad + \frac{z}{2}q_j\left( \frac{1}{\Theta_j}\frac{\partial \Theta_j}{\partial c_i} - \frac{1}{\Phi_j}\frac{\partial \Phi_j}{\partial c_i}\right) \\ & \quad - \frac{\partial \Phi_j}{\partial c_i}\frac{\mu}{x_j} - \frac{\Phi_j}{x_j}\frac{\partial \mu}{\partial c_i} + \frac{\Phi_j \mu}{x_j^2}\frac{\partial x_j}{\partial c_i} \end{align*} \]
As \(\sigma\), \(\tau\), and \(\mu\) are independent of species \(i\), these can be calculated outside the loop over the independent species. Moving to the residual term (Eq. 7), the mole fraction ( \(X_p\)) of group \(p\) in the mixture is related to the mole fraction ( \(x_j\)) of species \(j\) according to:
\[ X_p = c_{xX} \omega_p, \textrm{ where } \omega_p \equiv \displaystyle\sum_j v_p^{(j)} x_j, \\ \]
and \(c_{xX}\) is a conversion factor accounting for the difference in total species and total group number concentrations:
\[ c_{xX} = \frac{1}{\displaystyle\sum_p \omega_p}. \]
Partial derivatives of the group interaction terms in Eq 9, \(\Theta_m\) and \(\Psi_{mn}\), with respect to \(c_i\) are derived as follows:
\[ \begin{align*} \frac{\partial c_{xX}}{\partial c_i} & = - c_{xX}^2 \displaystyle\sum_p \displaystyle\sum_j v_p^{(j)} \frac{\partial x_j}{\partial c_i} = - \frac{c_{xX}^2}{\mathrm{MW}_im_T} \left( \displaystyle\sum_p v_p^{(i)} - \displaystyle\sum_p \omega_p \right), \\ & = \frac{c_{xX}}{\mathrm{MW}_im_T} \left( 1 - c_{xX} \displaystyle\sum_p v_p^{(i)} \right), \\ \frac{\partial X_n}{\partial c_i} & = c_{xX} \displaystyle\sum_j v_n^{(j)} \frac{\partial x_j}{\partial c_i} + \frac{\partial c_{xX}}{\partial c_i} \omega_n, \\ & = \frac{c_{xX}}{\mathrm{MW}_im_T} \left( v_n^{(i)} - \omega_n \right) + \frac{c_{xX}}{\mathrm{MW}_im_T} \left( \omega_n - c_{xX}\omega_n \displaystyle\sum_p v_p^{(i)} \right), \\ & = \frac{c_{xX}}{\mathrm{MW}_im_T} \left( v_n^{(i)} - X_n \displaystyle\sum_p v_p^{(i)} \right), \\ \Pi & \equiv \displaystyle\sum_n Q_n X_n, \\ \frac{\partial \Pi}{\partial c_i} & = \displaystyle\sum_n Q_n \frac{c_{xX}}{\mathrm{MW}_i m_T}\left( v_n^{(i)} - X_n \displaystyle\sum_p v_p^{(i)} \right), \\ & = \frac{c_{xX}}{\mathrm{MW}_i m_T} \left( \displaystyle\sum_n Q_n v_n^{(i)} - \Pi \displaystyle\sum_p v_p^{(i)} \right), \\ \frac{\partial \Theta_m}{\partial c_i} & = \frac{ \left(Q_m\frac{\partial X_m}{\partial c_i} \Pi - Q_m X_m \frac{\partial \Pi}{\partial c_i} \right)}{\Pi^2} , \\ & = \frac{c_{xX} Q_m}{\mathrm{MW}_im_T\Pi^2} \left[ \left( \Pi v_m^{(i)} - \Pi X_m \displaystyle\sum_p v_p^{(i)} \right) - \left( X_m \displaystyle\sum_n Q_n v_n^{(i)} - X_m \Pi \displaystyle\sum_p v_p^{(i)} \right) \right], \\ & = \frac{c_{xX} Q_m}{\mathrm{MW}_im_T\Pi^2} \left[ \Pi v_m^{(i)} - X_m\displaystyle\sum_n Q_n v_n^{(i)}\right], \textrm{ and} \\ \frac{\partial \Psi_{mn}}{\partial c_i} & = 0 \end{align*} \]
The partial derivative of the group residual activity coefficient (Eq. 8) with respect to \(c_i\) is:
\[ \begin{align*} \frac{\partial \ln{\Gamma_k}}{\partial c_i} & = Q_k \left[ - \frac{1}{\displaystyle\sum_m\Theta_m\Psi_{mk}} \displaystyle\sum_m\frac{\partial\Theta_m}{\partial c_i}\Psi_{mk} \\ \quad - \displaystyle\sum_m\left( \frac{\partial\Theta_m}{\partial c_i}\Psi_{km} \displaystyle\sum_n\Theta_n\Psi_{nm} - \Theta_m\Psi_{km}\displaystyle\sum_n \frac{\partial\Theta_n}{\partial c_i} \Psi_{nm}\right)/\left( \displaystyle\sum_n\Theta_n\Psi_{nm}\right)^2 \right] \end{align*} \]
After some rearranging, this becomes:
\[ \begin{align*} \Xi_m & \equiv \displaystyle\sum_n\Theta_n\Psi_{nm}, \\ \frac{\partial \ln{\Gamma_k}}{\partial c_i} & = - Q_k \displaystyle\sum_m \left( \frac{\Psi_{mk}}{\Xi_k}\frac{\partial\Theta_m}{\partial c_i} + \frac{\Psi_{km}}{\Xi_m} \frac{\partial\Theta_m}{\partial c_i} - \frac{\Theta_m\Psi_{km}} {\Xi_m^2} \displaystyle\sum_n\frac{\partial\Theta_n}{\partial c_i}\Psi_{nm} \right) \end{align*} \]
The left side of the three bracketed terms in the above equation are independent of species \(i\) and can be calculated outside of the loop over the independent species. The partial derivative of the full residual term with respect to \(c_i\) is:
\[ \begin{align*} \frac{\partial \ln{\Gamma_k^{(j)}}}{\partial c_i} & = 0 \\ \frac{\partial \ln{\gamma_j^R}}{\partial c_i} & = \displaystyle\sum_k v_k^{(j)} \frac{\partial \ln{\Gamma_k}}{\partial c_i} \end{align*} \]
The overall equation for the partial derivative of \(\gamma_j\) with respect to species \(i\) is:
\[ \frac{\partial\gamma_j}{\partial c_i} = \gamma_j \left( \frac{\partial\ln{\gamma_j^C}}{\partial c_i} + \frac{\partial\ln{\gamma_j^R}}{\partial c_i} \right) \]
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data |
J | Jacobian to be calculated |
time_step | Current time step [s] |
Definition at line 575 of file sub_model_UNIFAC.c.
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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
jac | Jacobian |
Definition at line 82 of file sub_model_UNIFAC.c.
void sub_model_UNIFAC_print | ( | int * | sub_model_int_data, |
double * | sub_model_float_data | ||
) |
Print the sub model data.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
Definition at line 771 of file sub_model_UNIFAC.c.
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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data |
Definition at line 129 of file sub_model_UNIFAC.c.
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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
deriv_ids | Indices for state array variables on the solver state array |
jac | Jacobian |
Definition at line 106 of file sub_model_UNIFAC.c.
void sub_model_ZSR_aerosol_water_calculate | ( | int * | sub_model_int_data, |
double * | sub_model_float_data, | ||
double * | sub_model_env_data, | ||
ModelData * | model_data | ||
) |
Do pre-derivative calculations.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data, including the state array |
Definition at line 216 of file sub_model_ZSR_aerosol_water.c.
void sub_model_ZSR_aerosol_water_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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data |
J | Jacobian to be calculated |
time_step | Current time step [s] |
Definition at line 314 of file sub_model_ZSR_aerosol_water.c.
void sub_model_ZSR_aerosol_water_get_used_jac_elem | ( | int * | sub_model_int_data, |
double * | sub_model_float_data, | ||
Jacobian * | jac | ||
) |
Flag Jacobian elements used by this sub model.
ZSR aerosol water sub models are assumed to be at equilibrium
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
jac | Jacobian |
Definition at line 81 of file sub_model_ZSR_aerosol_water.c.
void sub_model_ZSR_aerosol_water_print | ( | int * | sub_model_int_data, |
double * | sub_model_float_data | ||
) |
Print the ZSR Aerosol Water sub model parameters.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
Definition at line 445 of file sub_model_ZSR_aerosol_water.c.
void sub_model_ZSR_aerosol_water_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.
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
sub_model_env_data | Pointer to the sub model environment-dependent data |
model_data | Pointer to the model data |
Definition at line 189 of file sub_model_ZSR_aerosol_water.c.
void sub_model_ZSR_aerosol_water_update_ids | ( | int * | sub_model_int_data, |
double * | sub_model_float_data, | ||
int * | deriv_ids, | ||
Jacobian | jac | ||
) |
Update the time derivative and Jacbobian array indices.
ZSR aerosol water sub models are assumed to be at equilibrium
sub_model_int_data | Pointer to the sub model integer data |
sub_model_float_data | Pointer to the sub model floating-point data |
deriv_ids | Indices for state array variables on the solver state array |
jac | Jacobian |
Definition at line 133 of file sub_model_ZSR_aerosol_water.c.