# Resolution of Stiff Differential Systems in Biogeochemical Modeling using Some Advanced Numerical Schemes

Abbas El Hajj, Guillaume Dujardin, Christophe Preux, Benjamin Braconnier, Quang Huy Tran, Iryna Malinouskaya

► To cite this version:

Abbas El Hajj, Guillaume Dujardin, Christophe Preux, Benjamin Braconnier, Quang Huy Tran, et al. Resolution of Stiff Differential Systems in Biogeochemical Modeling using Some Advanced Numerical Schemes. 12ème Biennale de la Société des Mathématiques Appliquées et Industrielles, Jun 2025, Carcans, France. 2025. [hal-05379973](https://hal.science/hal-05379973v1)

HAL Id: hal-05379973

https://hal.science/hal-05379973v1

Submitted on 24 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Resolution of Stiff Differential Systems in Biogeochemical Modeling using Some Advanced Numerical Schemes

Abbas EL HAJJ¹,², Guillaume DUJARDIN², Christophe PREUX¹, Benjamin BRACONNIER¹, Quang-Huy TRAN¹, Iryna MALINOUSKAYA¹

¹IFP Energies Nouvelles, ²Inria Lille

## 1. Motivation and Objectives

Microbial activity plays a key role in modifying the physical and (geo)chemical conditions of natural and industrial environments. This is crucial for applications like geothermal energy, hydrogen storage, and water treatment, as part of the energy transition. For this reason, it is essential to integrate biogeochemical models of bacterial activity into porous media flow simulators. These models take into account phenomena such as the consumption of substrates contributing to the growth of organisms and their decline through biodegradation.

**Challenges:**

* Microbial activity models involve stiff differential equations.

* Numerical methods must be robust across a wide range of parameters, especially when coupled with porous media flow models.

* Ensuring positivity of concentrations.

## 2. Bacterial Reactions in the Subsurface

We focus on **sulfate-reducing bacteria (SRB)**, modeled by the simplified reaction:

$$SO_4^{2-} + 5H_2 + SRB \rightarrow H_2S + 4H_2O + \text{more SRB}$$

This describes the **consumption** of hydrogen and sulfate, the **production** of water and hydrogen sulfide (toxic and corrosive), and the **growth** of SRB.

## 3. Physical Model

We use a simplified model [1] to describe the reaction kinetics:

$$\begin{aligned} \frac{dS}{dt} &= -\gamma_S f(S, A)X, \\ \frac{dA}{dt} &= -\gamma_A f(S, A)X, \\ \frac{dX}{dt} &= Yf(S, A)X - bX \end{aligned}$$

* **S**: molar concentration of substrate [mmol/L]

* **A**: molar concentration of electron acceptor [mmol/L]

* **X**: biomass concentration [mg/L]

* $f(S, A) = \frac{\mu S}{Y\tau_S + S\tau_A + A}$: dual Monod reaction rate

* $\tau > 0$ half-saturation [mmol/L]

* $b > 0$ decay rate [1/day]

* $Y > 0$ yield coefficient [mg/mmol]

* $\gamma > 0$ stoichiometric coefficient [adim]

* $\mu > 0$ growth rate [1/day]

**Reference:** [1] de Blanc, P.C. (1998), *Development and Demonstration of a Biodegradation Model for Non-Aqueous Phase Liquids in Groundwater*, PhD Dissertation, University of Texas at Austin.

## 4. Model Configuration and Initial Conditions

* We use the same kinetic parameters studied by Robinson and Tiedje [2] and we also use their experimental hydrogen data for model validation.

* SRB (sulfate-reducing bacteria) consumes sulfate ions and $H_2$ to produce water and $H_2S$, the latter of which is both toxic and highly corrosive to metals and concrete.

<figure><img src="image_2.png" /></figure>

<table><tr><td><strong>Initial Conditions</strong></td><td>$S(0) = 0.0365$ mmol/L<br>$A(0) = 10$ mmol/L<br>$X(0) = 0.01$ mg/L</td></tr><tr><td><strong>Parameters</strong></td><td>$Y = 0.99$ mg/mmol<br>$\tau_A = 0.0001$ mmol/L<br>$\mu = 2.7$ day<sup>-1</sup><br>$b = 0.046$ day<sup>-1</sup></td></tr></table>

**Reference:** [2] Robinson J.A. and Tiedje J.M. (2005), *Kinetics of substrate depletion by mixed populations*, Appl. Environ. Microbiol.

## 5. Numerical Performances

* **System stiffness:** we set $A(0) = 60$ mmol/L and $X(0) = 0.1$ mg, and varied $S(0) \in \{1, 25, 60, 100\}$ mmol/L to study its effect on solver performance.

* **RK4** with very small time step used as a reference solution.

* **Numerical methods:** **Magnus3** ensures positivity and stability without constraints; **RK3_adaptive** needs small time steps for stability and adapts to preserve positivity; **Radau3_adaptive** and **Rosenbrock3_adaptive** are A-stable, and positivity-adaptive.

<figure><img src="image_3.png" /></figure>

## 6. Conclusion and Perspectives

* Growth of SRB biomass driven by the consumption of $H_2$, leading to the production of toxic $H_2S$.

* Increase in initial hydrogen concentration, causing greater system stiffness in our model.

* Preference for Runge-Kutta at low stiffness, Rosenbrock or Radau at high stiffness.

* Evaluating explicit Orthogonal-Runge-Kutta-Chebyshev (ROCKS) method.

* Integrating further biogeochemical processes into the current model.

* Coupling with porous media flow simulators for realistic applications.