A general paradigm to model reaction-based biogeochemical processes in batch systems

[1] This paper presents the development and illustration of a numerical model of reactionbased geochemical and biochemical processes with mixed equilibrium and kinetic reactions. The objective is to provide a general paradigm for modeling reactive chemicals in batch systems, with expectations that it is applicable to reactive chemical transport problems. The unique aspects of the paradigm are to simultaneously (1) facilitate the segregation (isolation) of linearly independent kinetic reactions and thus enable the formulation and parameterization of individual rates one reaction by one reaction when linearly dependent kinetic reactions are absent, (2) enable the inclusion of virtually any type of equilibrium expressions and kinetic rates users want to specify, (3) reduce problem stiffness by eliminating all fast reactions from the set of ordinary differential equations governing the evolution of kinetic variables, (4) perform systematic operations to remove redundant fast reactions and irrelevant kinetic reactions, (5) systematically define chemical components and explicitly enforce mass conservation, (6) accomplish automation in decoupling fast reactions from slow reactions, and (7) increase the robustness of numerical integration of the governing equations with species switching schemes. None of the existing models to our knowledge has included these scopes simultaneously. This model (BIOGEOCHEM) is a general computer code to simulate biogeochemical processes in batch systems from a reaction-based mechanistic standpoint, and is designed to be easily coupled with transport models. To make the model applicable to a wide range of problems, programmed reaction types include aqueous complexation, adsorption-desorption, ionexchange, oxidation-reduction, precipitation-dissolution, acid-base reactions, and microbial mediated reactions. In addition, user-specified reaction types can be programmed into the model. Any reaction can be treated as fast/equilibrium or slow/ kinetic reaction. An equilibrium reaction is modeled with an infinite rate governed by a mass action equilibrium equation or by a user-specified algebraic equation. Programmed kinetic reaction rates include multiple Monod kinetics, nth order empirical, and elementary formulations. In addition, user-specified rate formulations can be programmed into the model. No existing models to our knowledge offer these simultaneous features. Furthermore, most available reaction-based models assume chemical components a priori so that reactions can be written in basic (canonical) forms and implicitly assume that fast equilibrium reactions occur only for homogeneous reactions. The decoupling of fast reactions from slow reactions lessens the stiffness typical of these systems. The explicit enforcement of mass conservation overcomes the mass conservation error due to numerical integration errors. The removal of redundant fast reactions alleviates the problem of singularity. The exclusion of irrelevant slow reactions eliminates the issue of exporting their problematic rate formulations/parameter estimations to different environment conditions. Taking the advantage of the nonuniqueness of components, a dynamic basis-species switching strategy is employed to make the model numerically robust. Backward basis switching allows components to freely change in the simulation of the chemistry module, while being recovered for transport simulation. Three example