Reactive transport simulation concerns the coupled solution of chemical reaction and solute transport equations. Throughout the previous lectures in this course, we have focused our attention on the relevant chemical reaction equations and looked at some of the techniques for solving them in batch systems (i.e., zero-dimensional domains without consideration of transport). In this lecture, we will see how multi-species solute transport can be coupled with solution of the chemical reaction equations to form a complete reactive transport modeling approach that will be used in the example problems of the subsequent lectures (Topic C).

In systems in which there are multiple components existing within the aqueous phase, the advective-dispersive solute transport equation may be written for each component i:

(Eq.-1)

For our purposes in this course, we will assume that the groundwater velocity field is produced entirely in response to gradients in the potentiometric surface; that is, we are not considering density-driven flow. In systems with strong thermal or chemical concentration gradients, this assumption must be studied closely.

Figure 1. The two-step approach to reactive transport modeling.

As has been our strategy in many of the prior lectures in this course, letís see how we can couple solute transport and chemical reactions using a simple example problem (with PHREEQC, of course). Consider a situation in which we have rainwater that is in equilibrium with the average CO₂ present in the earthís atmosphere (partial pressure ~ 1 x 10^-3 atm.). Now suppose this rainwater infiltrates a column of unweathered igneous rock consisting entirely of K-feldspar, KAlSi₃O₈, a major constituent of granitic rocks, at a pore velocity of 50 cm/yr.

Figure 2. Infiltration of acidic rainwater into a column of K-feldspar.

Keep in mind that CO₂, when dissolved in water, forms a weak acid:

This acid may, in turn, facilitate the dissolution of K-feldspar,

The enrichment of the aqueous phase with the constituent species comprising may lead to saturation and precipitation of the weathering product kaolinite:

Letís look at a 1-m long column, divided into 20 equal-size volume elements. Initially, at any point in the column, given a porosity of 0.20 (fractures), we can calculate that 0.004 m³ or K-feldspar exists for each kg of water, corresponding to roughly 39 moles (assuming a density of 2.7 g/cm³ and a molecular weight of 278 g/mole). Over time, as infiltrating rainwater passes through the rock, K-feldspar will convert to kaolinite, a familiar observation to geologists.

So how do we set up this problem in PHREEQC? Much of the problem definition involves key words we have already learned how to use in the past few lectures. We need now only learn how to break the system up into separate cells that define the column and specify how solutes and transferred between the cells. This is accomplished under the "TRANSPORT" key word.

Keep in mind that the number of cell shifts, multiplied by the number of cells, implies that PHREEQC will need to solve 1,000,000 separate batch chemistry models for this simulation. This may require up to an hour or more of CPU time, depending, of course, on the type of processor. This should serve as a first warning as to the perils of reactive transport modeling!

Here are some of the results of this simulation to ponder:

Figure 3. Concentration profile of the Cl- tracer in the column after 1 year.

Figure 4. Amounts of K-feldspar and kaolinite present in the column after 5000 years.

Keep in mind that this particular example problem is highly idealized and contrived. In the real world, a system such as this would exhibit far more complex behavior because (1) a number of different mineral phases besides those specified here would be present, each serving to buffer pH, Al³⁺, H₄SiO₄, etc. in the aqueous phase, and (2) rates of mineral precipitation and dissolution throughout the column would reflect the influences of reaction kinetics, which were not simulated. Nevertheless, this example clearly illustrates a central, recurring concept of reactive transport in real hydrogeochemical systems: reaction fronts will almost always lag significantly behind the water movement front. The actual movement of infiltrating rainwater through the column, as delineated by the Cl^- tracer, takes approximately 2 years (the agreement between the PHREEQC numerical model and the Ogata-Banks (1961) solution is excellent). However, only a small portion of K-feldspar has been replaced by kaolinite in the first cell of the simulated column by 5000 years. Extrapolating, the reaction front as delineated by the mineral phases would require some 240,000 years before it begins to resemble the water movement front at 1 year. The message here is clear: the relatively low solubilities of mineral phases, the large number of moles of mineral phases present per unit volume in comparison to aqueous concentrations, and low rates of groundwater movement imply that substantial mineral dissolution requires a very long time. And keep in mind that we find this to be true under equilibrium assumptions; employing reaction kinetics in our modeled systems will slow things down still further. We will find this issue raising its head again and again in forthcoming simulations in subsequent lectures.

How would the results of this simple reactive transport demonstration change if the problem were made more realistic, perhaps by introducing other mineral phases? Try, for example, adding quartz ñ SiO₂, muscovite ñ KAl₃Si₃O₁₀(OH)₂, and a finite amount of kaolinite as initial phases in the column along with K-feldspar, or try including gibbsite, Al(OH)₃, as an additional weathering product (i.e., not initially present). What are the consequences? Can you explain some of the things that are happening?

You are now ready to continue to

TOPIC C: INTRODUCTION TO REACTIVE TRANSPORT MODELING WITH PHREEQC: EXAMPLE APPLICATIONS.

You may e-mail me questions and comments.

Walt W. McNab
E-mail address: mcnab1@llnl.gov.

Copyright Notice