Back to parent page - Gas-aerosol exchange

Time integration for SOA¶

Like for sulfate, MAM4 uses the approximation that the mass transfer coefficients stay constant within one time step \(\Delta t\). However, because of the dependence of the equilibrium SOAG mixing ratio on SOA and POA (the solvent effect), the ODEs are still nonlinear even when the mass transfer coefficients are constant in time. These nonlinear ODEs are solved with a semi-implicit Euler method with adaptive sub-steps, similar to the algorithms developed by Zaveri et al. (2008) for the Model for simulating aerosol interactions and chemistry (MOSAIC).

Notation¶

A sub-step advances the solution from \(t_{old}\) to \(t_{new}\), where \(t_{old} \geqslant t_n\) and \(t_{new} \leqslant t_{n+1}\). The lengh of the sub-step is \(\Delta t_{sub} = t_{new} - t_{old}\).

The values of physical quantities at \(t_n\), \(t_{old}\), and \(t_{new}\) are indicated by superscripts \(n\), \(old\), and \(new\), respectively.

Sub-step lengths¶

The length of a sub-step is determined by requiring that, if the equations were solved with an Euler forward scheme, the condensed or evaporated SOAG does not exeed a fraction \(\alpha\) of the existing or equilibrium amount. The calculations proceed as follows.

\[ \begin{align} q^{old}_{v,SOAG,sfc-equ,i} &= \dfrac{q_{m,SOA,i}^{old}} {\max \left( q_{m,OPOA,i}^n + q_{m,SOA,i}^{old}, \,\epsilon_a \right )} \, q_{v,SOAG,equ} \,,\quad i = 1, \cdots, I\\ \Phi_{i}^{old} &= \frac{q_{v,SOAG}^{old} - q_{v,SOAG,sfc-equ,i}^{old}} {\max \left( q_{v,SOAG}^{old}, q_{v,SOAG,sfc-equ,i}^{old}, \,\epsilon_g \right)} \,,\quad i = 1, \cdots, I \\ \Delta t_{sub} &= \min \left( \frac{\alpha} {\sum\limits_{i=1}^{I} \left ( C_{SOAG,i}^n \left| \Phi_{i}^{old} \right| \right)}, t_{n+1} - t_{old} \right) \,, \end{align} \]

where \(\alpha = 0.05\); \(\epsilon_a = \epsilon_g =\) 1E-20 mol/mol are safeguard parameters used for avoiding floating point exception.

Linearizing the ODEs¶

To linearize the ODEs, we can first rewrite the evolution equations of SOAG and SOA mass mixing ratios as follows:

\[ \begin{align} \dfrac{dq_{m,SOA,i}}{dt} &= C_{SOA,i}\left( q_{v,SOAG} - \mathcal{S}_{SOAG,i}\,q_{m,SOA,i} \right) \,,\quad i = 1,\cdots,I\\ \dfrac{dq_{v,SOAG}}{dt} &= -\sum\limits_{i=1}^I \dfrac{q_{m,SOA,i}}{dt} \end{align} \]

where

\[ \mathcal{S}_{SOAG,i} = \dfrac{ q_{v,SOAG,equ} }{ q_{m,SOG,i}+q_{m,OPOA,i} } \]

Approximating \(C_{SOAG,i}\)¶

The mass transfer coefficients are fixed at their value at \(t_n\) and denoted by \(C_{SOAG,i}^n\).

Approximating \(\mathcal{S}_{SOAG,i}\)¶

An approximate constant between \(t_{old}\) and \(t_{new}\), denoted by \(\mathcal{S}_{SOAG,i}^*\), is estimate by

\[ \mathcal{S}_{SOAG,i}^* = \dfrac{ q_{v,SOAG,equ} }{ \max\left( q_{m,OPOA,i}^n + q_{m,SOA,i}^*, \,\epsilon_a \right)} \]

For a mode with \(\Phi_i \leqslant 0\), for which SOAG will likely be evaporating during \((t_{old}, t_{new})\), we let

\[q_{m,SOA,i}^* = q_{m,SOA,i}^{old} \]

For a mode with \(\Phi_i > 0\), for which SOAG will likely be condensing, we let

\[ q_{m,SOA,i}^* = q_{m,SOA,i}^{old} + \Delta t_{sub}\, C_{SOAG,i}^n \left( q_{v,SOAG}^{old} - q^{old}_{v,SOAG,sfc-equ,i} \right) \]

The linearized ODEs¶

\[ \begin{align} \dfrac{dq_{m,SOA,i}}{dt} &= C_{SOA,i}^n\left( q_{v,SOAG} - \mathcal{S}_{SOAG,i}^*\,q_{m,SOA,i} \right) \,,\quad i = 1,\cdots,I\\ \dfrac{dq_{v,SOAG}}{dt} &= -\sum\limits_{i=1}^I \dfrac{q_{m,SOA,i}}{dt} \end{align} \]

Discretizing the linearized ODEs¶

We discretize the linearized SOA mass equations using an Euler backward scheme

\[ q_{m,SOA,i}^{new} - q_{m,SOA,i}^{old} = \Delta t_{sub} C_{SOAG,i}^n \left( q_{v,SOAG}^{new} - S^*_{SOAG,i}\,q_{m,SOA,i}^{new} \right)\,,\quad i=1,\cdots,I \]

and define

\[\beta_{SOAG,i}^n = \Delta t_{sub} C_{SOAG,i}^n \]

This gives

\[ q_{m,SOA,i}^{new} = \dfrac{ q_{m,SOA,i}^{old} + \beta_{SOAG,i}^n \, q_{v,SOAG}^{new} } { 1+ \beta_{SOAG,i}^n S^*_{SOAG,i} } \,,\quad i=1,\cdots,I \]

The total mass of SOA in both gas and liquid phases is conserved, hence for the time interval \((t_{new},\,t_{new})\), the second equation in Section 3.4 can be written as

\[ q_{v,SOAG}^{new} + \sum\limits_{i=1}^I q_{m,SOA,i}^{new} = q_{v,SOAG}^{old} + \sum\limits_{i=1}^I q_{m,SOA,i}^{old} \]

Substituting the \(q_{m,SOA,i}^{new}\) in this equation by the expression obtained with the Euler backward scheme gives the following final expression for \(q_{m,SOA,i}^{new}\):

\[ \begin{align} q_{v,SOAG}^{new} &= \dfrac { q_{v,SOAG}^{old} + \sum\limits_{i=1}^I q_{m,SOA,i}^{old} - \sum\limits_{i=1}^I q_{m,SOA,i}^{old}/\left(1+ \beta_{SOAG,i}^n S^*_{SOAG,i}\right) } { 1+ \sum\limits_{i=1}^I \beta_{SOAG,i}^n /\left( 1+ \beta_{SOAG,i}^nS^*_{SOAG,i}\right) } \\ \end{align} \]

Inserting this expression to the third last equation above gives the final expression for the SOA mass mixing ratios \(q_{m,SOA,i}^{new}\,,\,i = 1,\cdots,I\).