Page 89 - Fundamentals of Computational Geoscience Numerical Methods and Algorithms

76 4 Algorithm for Simulating Fluid-Rock Interaction Problems

Lasaga 1994, Raffensperger and Garven 1995). Considering this fact, we have suc-
cessfully used the finite element method to model coupled reactive multi-species
transport problems with homogeneous reactions (Zhao et al. 1999a). Here we will
extend the numerical method developed to model coupled reactive multi-species
transport problems with heterogeneous reactions.

4.1 Key Issues Associated with the Numerical Modelling
of Fluid-Rock Interaction Problems

Of central importance to the numerical modelling of fluid-rock interaction problems
in pore-fluid saturated hydrothermal/sedimentary basins is the appropriate consid-
eration of the heterogeneous chemical reaction which takes place slowly between
aqueous reactive species in the pore-fluid and solid minerals in fluid-saturated
porous rock masses. From the rock alteration point of view, the pore-fluid flow,
which carries reactive aqueous species, is the main driving force causing the hetero-
geneous chemical reactions at the interface between the pore-fluid and solid miner-
als so that the rock can be changed from one type into another. The reason for this
is that these reactions between the aqueous species and solid minerals may result
in dissolution of one kind of mineral and precipitation of another and therefore, the
reactive minerals may be changed from one type into another type. Due to the dis-
solution and precipitation of minerals, the porosity of the porous rock mass evolves
during the rock alteration. Since a change in porosity can result in a change in the
pore-fluid flow path, a feed-back loop is formed between the pore-fluid flow and
the transport of reactive chemical species involved in heterogeneous chemical reac-
tions in fluid-rock interaction systems and this porosity change in the pore-fluid
flow needs to be considered in the numerical modelling of fluid-rock interaction
problems. This implies that the average linear velocity of pore-fluid flow varies with
time due to this porosity evolution. Since both the mesh Peclet number and Courant
number are dependent on the average linear velocity, we have to overcome a dif-
ficulty in dealing with the problem of variable Peclet and Courant numbers in the
transient analysis of fluid-rock interaction problems. Generally, there are two ways
to overcome this difficulty, from the computational point of view. The first one is to
use a very fine mesh and very small time step of integration so that the requirements
for both the mesh Peclet number and Courant number can be satisfied at every time
step of the computation. The second one is to regenerate the mesh and re-determine
the time step of integration at every time step of computation so that the finite ele-
ment method can be used effectively. Since either the use of a very fine mesh in the
whole process of computation or the regeneration of the mesh at every time step of
computation is computationally inefficient, there is a definite need for developing
new numerical algorithms to deal with this kind of problem.
Another important issue related to the numerical modelling of fluid-rock interac-
tion problems is that the dissolution rates of minerals are dependent on the existence
of the dissolving minerals. Once a dissolving mineral is exhausted in the rock mass,