Page 92 - Fundamentals of Computational Geoscience Numerical Methods and Algorithms

4.2 Development of the Term Splitting Algorithm 79

the vertical cross-sections, which are perpendicular to the direction of the pore-fluid
flow, may have different values of the areal porosity (Zhao 1998e, Phillips 1991,
Nield and Bejan 1992, Zhao et al. 1994c). However, the linear averaged velocity,
which is the velocity averaged over the pore area of a representative elementary
area, varies if the vertical cross-sections have different values of the areal porosity.
It is also noted that Eq. (4.9) is a typical transport equation for species i.This
equation describes the mass conservation of species i in the pore-fluid saturated
porous medium for both constant and variable porosity cases (Zhao et al. 1998e,
1994c).
In order to consider the feedback effect of porosity due to heterogeneous chemi-
cal reactions, we need a porosity evolution equation in the numerical analysis. Using
the generalized concentration of solid minerals, the evolution equation for the poros-
ity of the rock mass can be expressed as follows:
m
W i

φ = 1 − C Gi , (4.11)
ρ s
i=1
where W i is the molecular weight of mineral i; ρ s is the unit weight of solid minerals.
For the transient analysis, C Gi are the functions of space and time variables, so that
the porosity of the rock mass varies with both space and time. This is the reason
why Eq. (4.11) is called the evolution equation of porosity with time.
From the chemical kinetics point of view, the reaction rates involved in Eqs. (4.8)
and (4.9) can be expressed as follows (Lasaga 1984):

N

R i =− τ ij r j , (4.12)
j=1
where τ ij is the contribution factor of mineral/species j to the reaction rate of min-
eral/species i; N is the number of the heterogeneous chemical reactions to be con-
sidered in the rock mass; r j is the overall dissolution rate of mineral j.

A j Q j
r j = k j 1 − , (4.13)
V f K j
where A j is the surface area of mineral j; k j and K j are the reaction constant and the
equilibrium constant of the jth heterogeneous chemical reaction; V f is the volume of
the solution; Q j is the chemical affinity of the jth heterogeneous chemical reaction.
It needs to be pointed out that the surface area of a dissolving mineral varies in
the process of mineral dissolution. Once the dissolving mineral is depleted in the
rock mass, the surface area of the mineral must be equal to zero and therefore, the
chemical reaction dissolving this particular mineral should stop. If this area is sim-
ply assumed to be a constant in the numerical modelling, then the resulting numer-
ical solutions are incorrect because the real mechanism of chemical kinetics cannot
be simulated during the consumption of the dissolving mineral in the rock mass.