|
NEML2 2.0.0
|
This example demonstrates the use of the chemical reactions physics module to compose models for reactive infiltration kinetics, as a solid \((S)\) material is infiltrated by a liquid \((L)\) to create a product \((P)\). The framework and usage of the model is explained below, with both mathematical formulations and example input files.
The infiltration process for a material considers a cylindrical Representative Volume Element (RVE) with radius \(R\) and span \(H\), depicting a solid capillary of radius \(r_o\), corresponding to the initial porosity of \(\varphi_0\). The mathematical description of the model uses cylindrical coordinates and assumes axisymmetry around \(r = 0\).
As the liquid enters the cylinder, it reacts with a solid to form a product with thickness \(\delta_P\), as shown schematically in the figure below.
Throughout this process, the solid radius, \(r_o\), decreases while the product thickness, \(\delta_P\), increases. Let \(r_i\) denote the inner radius of the product. In addition, define \(\alpha_i\) (in units of mole per volume) as the amount of substance in the RVE, and \(\Omega_i = \dfrac{M_i}{\rho_i}\) as the molar volume of a material with molar mass \(M_i\) (in units of amu, or gram per mole) and density \(\rho_i\) (mass per volume), with subscripts taking \(L\), \(S\), and \(P\), respectively.
The volume fraction, \(\varphi_i\) of each material is then
\[ \varphi_i = \alpha_i \Omega_i \]
and the RVE porosity (void) is
\[ \varphi_v = 1 - \varphi_L - \varphi_P - \varphi_S \]
Key assumptions made throughout the derivation are:
The following nondimensionalization is applied to the constitutive model for the reactive infiltration process.
\begin{align*} &\ \bar{\delta}_P = \frac{\delta_P}{R}, \quad \bar{r}_o = \dfrac{r_o}{R} = \sqrt{1-\varphi_S}, \quad \bar{r}_i = \dfrac{r_i}{R} = \sqrt{1-\varphi_S-\varphi_P}. \end{align*}
The complete state of the RVE is denoted by the tuple \(\left( \varphi_L, \varphi_S, \varphi_P\right)\), with \(\alpha_L\) as the prescribed force.
Mathematically, it is possible that \( \varphi_L + \varphi_P + \varphi_S \ge 1 \). Physically, this implies "overflow", aka the prescribed \(\alpha_L\) is larger than the available voids. Care must be taken at the macroscopic model to avoid or resolve this issue.
The initial-value problem (IVP) corresponding to the constitutive model is
\[ \mathbf{r} = \begin{Bmatrix} r_{\varphi_P} \\ r_{\varphi_S} \end{Bmatrix} = \mathbf{0}, \]
where \(r_{\varphi_i}\) is the residual corresponding to the volume fraction of the liquid, solid and product. The residuals are defined below. Note an implicit backward-Euler time intergration scheme is used.
\[ r_{\varphi_P} = \dfrac{\alpha_{P, n+1} - \alpha_{P, n}}{t_{n+1}-t_n} - \dfrac{2}{\Omega_L}\dfrac{D}{l_c^2}\dfrac{\bar{r}_o+\bar{r}_i}{\bar{r}_o-\bar{r}_i} R_L R_S \]
Here \(R_L(\varphi_L), R_S(\varphi_S) \in [0, 1]\) are the reactivity of liquid and solid, represented by a step function of the void fraction, \(\varphi_i\).\[ r_{\varphi_S} = \dfrac{\alpha_{S, n+1} - \alpha_{S, n}}{t_{n+1}-t_n} + \dfrac{k_P}{k_S} \dfrac{\alpha_{P, n+1} - \alpha_{P, n}}{t_{n+1}-t_n} \]
The above governing equations are fairly general for a wide range of reactive infiltration systems within this framework given the chemical reaction kinetics with appropriate reaction coefficient ( \(k_i\)) and the corresponding molar volume \(\Omega_i\) of the materials.
Example for determining \(k_i\) and \(\Omega_i\): Consider a reactive infiltration process of liquid (chemical formula \(L_x\)) and solid ( \(S_y\)), with the molar mass \(M_L\) and \(M_S\). The chemical reaction creates a product with chemical formula \(L_m S_n\), with reaction coefficients \(k_S\) and \(k_P\).
\[ L_x + k_S S_y \rightarrow k_P L_mS_n \]
For stoichiometric balance, the reaction coefficients are:
\begin{align*} k_S = \dfrac{xn}{ym}, \quad k_P = \dfrac{x}{m} \end{align*}
The densities of the liquid, solid, and product are respectively \(\rho_L, \rho_S, \rho_P\).
Then, the molar volume is:
\[ \Omega_L = \dfrac{x M_L}{\rho_L}. \quad \Omega_S = \dfrac{y M_S}{\rho_S}, \quad \Omega_P = \dfrac{mM_L+nM_S}{\rho_P}. \]
The instantaneous balance also implies that
\[ \dot{\alpha}_S = \dfrac{k_P}{k_S} \dot{\alpha}_P, \]
Since this constitutive model considers \(\alpha_L\) as the forcing function for the IVP problem, \(\dot{\alpha}_L \ne k_P \dot{\alpha}_P\).
The following table summarizes the relationship between the mathematical expressions and the NEML2 models.
| Expression | Syntax |
|---|---|
| \( \bar{r}_i = \sqrt{1-\varphi_S-\varphi_P}; \quad \bar{r}_o = \sqrt{1-\varphi_S}\) | CylindricalChannelGeometry |
| \( \dfrac{2}{\Omega_L}\dfrac{D}{l_c^2}\dfrac{\bar{r}_o+\bar{r}_i}{\bar{r}_o-\bar{r}_i} R_L R_S\) | DiffusionLimitedReaction |
| \( R_L; \quad R_S \) | HermiteSmoothStep |
| \( \varphi_i = \alpha_i \Omega_i \) | Linear Combination |
| \( r_{\varphi_P} = \dfrac{\alpha_{P, n+1} - \alpha_{P, n}}{t_{n+1}-t_n} - \dfrac{2}{\Omega_L}\dfrac{D}{l_c^2}\dfrac{\bar{r}_o+\bar{r}_i}{\bar{r}_o-\bar{r}_i} R_L R_S\) | Linear Combination, Variable Rate |
| \(r_{\varphi_S} = \dfrac{\alpha_{S, n+1} - \alpha_{S, n}}{t_{n+1}-t_n} + \dfrac{k_P}{k_S} \dfrac{\alpha_{P, n+1} - \alpha_{P, n}}{t_{n+1}-t_n}\) | Linear Combination, Variable Rate |
A complete example input file for the reactive infiltration is shown below with the appropriate model composition.