Geochemical modeling of diagenetic reactions in Snorre Field reservoir sandstones : a comparative study of computer codes

Manuscript ID: 30145. Received: 07/31/2014. Approved: 06/22/2015. ABSTRACT: Diagenetic reactions, characterized by the dissolution and precipitation of minerals at low temperatures, control the quality of sedimentary rocks as hydrocarbon reservoirs. Geochemical modeling, a tool used to understand diagenetic processes, is performed through computer codes based on thermodynamic and kinetic parameters. In a comparative study, we reproduced the diagenetic reactions observed in Snorre Field reservoir sandstones, Norwegian North Sea. These reactions had been previously modeled in the literature using DISSOL-THERMAL code. In this study, we modeled the diagenetic reactions in the reservoirs using Geochemist’s Workbench (GWB) and TOUGHREACT software, based on a convective-diffusive-reactive model and on the thermodynamic and kinetic parameters compiled for each reaction. TOUGHREACT and DISSOLTHERMAL modeling showed dissolution of quartz, K-feldspar and plagioclase in a similar temperature range from 25 to 80°C. In contrast, GWB modeling showed dissolution of albite, plagioclase and illite, as well as precipitation of quartz, K-feldspar and kaolinite in the same temperature range. The modeling generated by the different software for temperatures of 100, 120 and 140°C showed similarly the dissolution of quartz, K-feldspar, plagioclase and kaolinite, but differed in the precipitation of albite and illite. At temperatures of 150 and 160°C, GWB and TOUGHREACT produced different results from the DISSOL-THERMAL, except for the dissolution of quartz, plagioclase and kaolinite. The comparative study allows choosing the numerical modeling software whose results are closer to the diagenetic reactions observed in the petrographic analysis of the modeled reservoirs.


INTRODUCTION
Petroleum systems are dynamic physical-chemical environments of generation, migration and storage of oil.Initially, oil is generated in a rock rich in organic matter that suffers a maturation process in specific depth and temperature (source rock).The oil then migrates to a porous lithology that will store it (the reservoir rock), where it is trapped by another rock that has low permeability (the seal) (Magoon & Dow 1994;Gluyas & Swarbrick 2004).
Unlike the seal rocks, reservoir rocks are characterized by high porosity and permeability.The reservoirs are normally carbonate rocks (Roehl & Choquette 1985) or siliciclastic rocks (Barwis et al. 1990), formed as the result of several transformations that occur in the deposited sediments, via modification of the original texture, composition and, consequently, in the porosity.Such processes are called diagenesis.
The most common diagenetic processes that affect the porosity and permeability of reservoirs are compaction and precipitation of various types of carbonates, silicates and other minerals (Tucker 2001).If precipitated as cements in significant amounts in the space between the primary depositional constituents, these diagenetic minerals compromise the porosity of the rocks, making them unfeasible as a reservoir.Nevertheless, there is also the possibility of subsequent dissolution of diagenetic cements, generating secondary porosity, as occurs in many reservoirs (Bjorlykke 1984;Schmidt & McDonald 1979;Surdam et al. 1984).
The diagenetic processes allow a variety of minerals to occur as cement in sedimentary rocks.The most common cements are: 1. quartz in the form of secondary overgrowths around detrital grains and intergranular microcrystals; 2. K-feldspar as secondary overgrowths or discrete crystals filling pores; 3. carbonates, like calcite and dolomite, siderite and more rarely ankerite, rhodochrosite and magnesite; 4. clay minerals, like kaolinite, illite, chlorite, smectite and interstratified mixed-layers (especially illite-smectite).
Clay minerals are volumetrically less significant, but extremely important, since they commonly occur as rims or coatings around grains, exerting enormous influence on the permeability of rocks, damaging the quality of reservoirs (Houseknecht & Pittman 1992).
Chemical phenomena are the most common during diagenesis, comprising reactions that result in precipitation of authigenic minerals or secondary dissolution.From the geochemical point of view, these reactions occur because the constituents of the sediments are always trying to reach equilibrium, and, therefore, tend to interact with the interstitial fluids through diagenetic processes (Burley et al. 1985).The diagenetic processes are controlled not only by temperature (T), pressure (P) and the original mineral assemblage, but also by the activities of the ions dissolved in the formation waters, including organic solutes, Eh and pH conditions.Galloway (1984) establishes three main hydrological regimes for fluids within sedimentary basins, related to their movement and occurrence: 1. meteoric regime, affecting the shallower portions of basins, characterized by infiltration of surface to groundwaters; 2. compactional regime, linked to expulsion of pore waters due to burial; 3. thermobaric or abyssal regime, associated to the deeper parts of the basins, where the interstitial fluids are affected by clay minerals dehydration and oil generation.
Geochemical modeling is a tool used to understand the diagenetic processes.It is performed through computer codes based on thermodynamic and kinetic parameters.However, since there are several mutually-influencing parameters, each code must be used carefully in order to better simulate the natural conditions.In this study, we reproduced the diagenetic reactions observed in the Snorre Field sandstones reservoir, Norwegian North Sea using Geochemist's Workbench (GWB) and TOUGHREACT software, based on a convective-diffusive-reactive model, and the thermodynamic and kinetic parameters compiled for each reaction.The Snorre Field was chosen considering the availability of petrographic and stratigraphic data to build a conceptual model of the diagenetic reactions (Baccar & Fritz 1993;Morad et al. 1990).Nevertheless, such study has potential application to the modeling of Brazilian reservoirs with similar composition, depositional and burial evolution.
Baccar and Fritz (1993) modeled the diagenetic reactions in Snorre Field reservoir sandstones using DISSOL-THERMAL code.We reproduced the same study using GWB and THOUGHREACT codes, so that we could evaluate similarities and divergences in the results of each system according to a wide temperature range.

GEOCHEMICAL MODELING APPROACH
Although the fundamental concepts of chemical modeling have evolved much before, geochemical modeling begins in the 1960s, with the pioneering efforts of researchers such as Garrels and Mackenzie (1967), Helgeson (1967aHelgeson ( , 1967b)), and Helgeson and James (1968).The first simulations were originally applied to the understanding of the basic issues involving chemical reactions in aquatic environments related to pollution of surface waters and assessment of diagenetic processes (involving natural formation and alteration of rocks).
A geochemical model is really only useful as a tool for prediction if there is the possibility to validate the results.In reality, this is the goal that most often is unattainable, due to the complexity of natural systems, insufficient field data and uncertainties related to changes in the system over time (Nordstrom 1994).
In the following sections, the main characteristics of TOUGHREACT and GWB, software used in this study, are described.We also give some characteristics of the DISSOL-THERMAL code used by Baccar and Fritz (1993).

TOUGHREACT
TOUGHREACT was developed to simulate fluid flow non-isothermal multi-component and geochemical transport, in order to investigate problems involving simple and complex geological environments (Xu & Pruess 1998).TOUGHREACT can also use "batch models", simulating the reactions in a closed system.A number of thermo-physical and chemical processes, which occur in the subsurface, are considered under the aspect of hydrological and geochemical conditions of pressure, temperature, water saturation and ionic strength.The code uses non-isothermal equations (Pruess 1987;Pruess et al. 1999).The numerical solution to the water-rock interaction is based on finite difference (Narasimhan & Witherspoon 1976).These equations are solved by the Newton-Raphson interaction (Pruess 1991).The activity coefficient is calculated using the Debye-Huckel equation (Helgeson et al. 1981).
Kinetic rates (mineral dissolution/precipitation) can be expressed as function of non-basis species as well.Usually the species appearing in rate laws happen to be basis species.In this model, we use a rate expression given by Lasaga et al. (1994): where positive values of r n indicate dissolution, and negative values indicate precipitation; k n is the rate constant (moles per mineral surface area unit and time unit), which is temperature dependent; A n is the specific reactive surface area per kg; and Q/K is the kinetic mineral saturation ratio.The temperature dependence of the reaction rate constant can be expressed reasonably well via an Arrhenius equation (Lasaga 1984;Steefel and Lasaga 1994).Because many rate constants are reported at 25°C, it is convenient to approximate rate constant dependency as a function of temperature, thus: where e is Euler's number; Ea is the activation energy; k 25 is the reaction rate at 25°C; R is gas constant; and T is absolute temperature.
The kinetic rate constant k n in equations 1 and 2 only considers the most well-studied mechanisms in pure H 2 O (at neutral pH).Dissolution and precipitation of minerals are often catalyzed by H + (acid mechanism) and OH -(base mechanism).For many minerals, the kinetic rate constant k includes each of these three mechanisms (Lasaga et al. 1994;Palandri & Kharaka 2004).

Geochemist's Workbench
The GWB system was developed by the Department of Geology at the University of Illinois at Urbana-Champaign in 1978.It is currently marketed by "Aqueous Solutions LLC" (http://www.gwb.com/).The software allows fast assessment of the problems most commonly encountered in geochemical modeling: balance reactions speciation in aqueous solution, water-rock interaction models, and fluid mixture in conditions of "batch model" (Bethke 2008).
GWB can integrate thermodynamic and kinetic parameters with species processing of geochemical reactions in redox environments (Drever 1997).The batch modeling is a process by which a series of irreversible reactions (no equilibrium) progress to a state of thermodynamic equilibrium (Helgeson et al. 1969).These reactions lead to chemical evolution of a system to reach the equilibrium state.To solve the equilibrium of an initial system, GWB determines which reactions are more likely to happen.
GWB solves simultaneously a set of non-linear algebraic equations, composed by the equation of mass action and the equilibrium constant, corresponding to aqueous species, gas or minerals that are contained in a database.While GWB solves the equation of mass action, it is solving all the mass balance of the system.The GWB solves the equations using Newton-Raphson method.
The GWB allows the user to build the rate law for the reactions kinetics.Settings are available for kinetic dissolution/precipitation reactions, redox reactions and microbial metabolism.The pattern that GWB uses for the equation of the rate law for the dissolution and precipitation of minerals takes the following form: where r k is reaction rate; A S is mineral surface area; k + is rate constant; Q is activity product; and K is equilibrium constant.The ratio Q/K is known as mineral saturation index.In many cases, changes in the temperature will affect the rates of reactions in progress.In GWB, the activation energy and pre-exponential factor can be used to specify a temperature-dependent rate, based on the formulation of Arrhenius: where A is pre-exponential factor; Ea is activation energy; R is gas constant; and T is absolute temperature.
Diagenesis is defined as the sum of physical and chemical changes during progressive burial and heating, reflected in compaction owing to known physical processes and to reactions where minerals are dissolved or precipitated.A wide variety of diagenetic reactions occur in sandstones, and some have significant effects on porosity and permeability.

Mathematical correlations between TOUGHREACT and Geochemist's Workbench, and kinetic parameters
The mathematical correlations of TOUGHREACT resemble in part those of GWB, but there is an important difference between them.Solving equation 2, we have: Thus, the correlation of GWB with TOUGHREACT is: The major difference observed in TOUGHREACT is the numerical expression: The kinetic parameters used in the simulations are shown in Tab. 1.Each software has different numerical methods to solve the diagenetic reactions observed in the geological environment.Thus, the kinetic parameters adopted by each simulator are submitted to some adjustment factors.
After the kinetic parameters are inserted in each software interface, simulators calculate the "global rate constant", whose values are in Tab. 2.
If we look at two different reference books on the equilibrium constant value of a chemical reaction, there is a big chance that the values found are different (sometimes by a factor of 10 or more).This discrepancy is because the constant value may have been determined in different conditions and perhaps using different techniques.A common source of variation in published values of K is the ionic composition of the solution (Harris 2012).
The diagenetic reactions used by GWB and TOUGHREACT are shown in Tab. 3.For GWB, the equations are written with a view to all possible chemical equilibria, as a function of hydrogen ion.The GWB uses hydrogen ions (H + ) to assist in the dissociation of minerals, thus enabling the adjustment of the stoichiometric coefficients.

DISSOL-THERMAL
The DISSOL-THERMAL simulator (Fritz 1975(Fritz , 1981) is a geochemical code able to predict whether the mineral constituents are dissolved or precipitated, as well as the composition of the solution after the water-rock interaction.Through the attributed concentration of dissolved species, temperature, pH, and redox potential, DISSOL-THERMAL code provides the diagenetic reactions and calculates the saturation index of minerals present in the geological environment.
DISSOL-THERMAL is divided into two geochemical codes: 1. DISSOL -calculates a model for simulating water-rock interactions at a given temperature between 0° and 300°C; and 2. THERMAL -calculates a model for simulating the effect of temperature variation on water-rock equilibrium.
These two codes both operate by combining an initial solution saturation test, followed by path calculation (dissolution or temperature variation).
In this study, we used the DISSOL-THERMAL modeling results from Baccar and Fritz (1993) and compared to our modeling using GWB and TOUGHREACT simulations.

Geological setting and mineralogy of Snorre Field reservoirs
The Snorre Field is located in the northern part of the North Sea (Fig. 1).The oil field is situated on the Tampen Spur, which is a platform high on the western side of the Viking Graben (Jurassic-Cretaceous).The main reservoir horizons are fluvial sandstones from the upper member of the Upper Triassic Lunde Formation and the Upper Triassic to Lower Jurassic Statfjord Formation (Hollander 1987).
The Lunde Formation, which is the subject of this paper, is the uppermost formation of the Triassic Hegre Group (Vollset & Dore 1984).It was deposited during the thermal subsidence phase following a Late Permian to Earliest Triassic rifting episode (Badley et al. 1988;Nystuen et al. 1989).The Triassic to Middle Jurassic sequence on the Snorre Field was uplifted, tilted and partly eroded during the Kimmeridgian rifting episode in Late Jurassic to Early Cretaceous times.The uplift was followed by subsidence and burial during crustal cooling from Early Cretaceous to Recent (Morad et al. 1990).
Subordinate constituents include extrabasinal quartz-feldspar-mica rock fragments that probably represent granitic rocks and/or schists or gneisses.They also include mud and carbonate intraclasts.Micas and detrital clay minerals seldom make up more than 2% of the total mineral content in the sandstones.Clay minerals of detrital origin include smectite, mixed-layer clay minerals, chlorite and subordinate amounts of kaolinite and illite.

Minerals GWB TOUGHREACT
chlorite.Diagenetic kaolinite also replaces detrital micas and feldspars.Other cements include carbonates, which play a significant role in porosity reduction in some of the sandstones.The most dominant are calcite and ferroan-calcite.Dolomite and ankerite, however, are also present in many of the samples, especially above the oil-water contact, and within 100 m of the Kimmeridgian unconformity (Fig. 2).Calcite occurs mainly as nodular to patchy pore-filling or replacive poikilotopic cement.Replacement of calcite by ferroan-dolomite was observed by Morad et al. (1990).Ferroandolomite and ankerite are commonly associated with detrital biotite.Siderite is a subordinate cement in these sandstones, but does occur in the flood plain sediments, predominantly as nodules.Authigenic quartz overgrowths (< 1%) are rarely observed.
The diagenetic evolution of the upper member of the Lunde Formation has been controlled by water-rock reactions acting during deposition, burial through the cooling phase succeeding the Permo-Triassic rifting episode, Kimmeridgian rifting and uplift, and post-Kimmeridgian burial.Over this time interval, four major diagenetic regimes (cf. Schmidt & McDonald 1979)   Table 3. Dissociation reaction of the minerals according to the software.

Minerals Dissociation reaction
4. Mesodiagenesis in the mid-Cretaceous to recent.

Geochemical modeling of diagenetic events
Baccar and Fritz (1993) studies were focused on the mesodiagenesis, in which there were several interactions between the rocks and modified waters in different temperatures.Our studies are referred to the same conditions.The scenario that will be modeled in this study is the interaction of Lunde Formation with a modified seawater in different conditions of temperature, which correspond to mesodiagenesis conditions (Morad et al. 1990).The primary rock is a sandstone with an average simplified composition of 70% quartz, 20% feldspars and 10% clays.The complete mineral composition of this rock is detailed in Tab. 4 and the quantity of each mineral is represented by volume fraction dissolved per kg of H 2 O. Seawater composition was compiled from Nordstrom et al. (1979), and the major element composition is shown in Tab. 5. Diagenetic events are modeled in the following temperatures: ■ Event 1: 25, 40, 60 and 80°C; ■ Event 2: 100, 120 and 140°C; ■ Event 3: 150 and 160°C.

RESULTS
Diagenetic processes, as predicted by modeling the waterrock interaction between 25 and 160°C using TOUGHREAT and GWB, in comparison to DISSOL-THERMAL, are summarized in Fig. 3.The modeling with TOUGHREACT and DISSOL-THERMAL showed dissolution of quartz, K-feldspar, plagioclase, kaolinite and illite between 25 and 160°C.When GWB was used, quartz, K-feldspar, albite, illite and kaolinite were precipitated.
At temperatures of 100, 120 and 140°C, the 3 different modeling systems indicate the dissolution of quartz, K-feldspar, plagioclase and kaolinite, converging regarding their results.
At temperatures of 150 and 160°C, GWB and TOUGHREACT produce different results from those of DISSOL-THERMAL with dissolution of K-feldspar.
Above 100°C, we observed similar behavior in GWB and TOUGHREACT modeling due to the high kinetics and thermodynamics involved in diagenetic reactions, favoring the polynomial equations used by simulators.

DISCUSSION
As showed in the results, the use of different codes produces different results during modeling.The convergence of results is more evident between TOUGHREACT and DISSOL-THERMAL than between DISSOL-THERMAL and GWB software.Generally at low temperatures, the simulation results show some discrepancies, because each software uses different numerical methods to solve a different set of reactions.
The choice of the software for geochemical modeling should consider the differences of each code, in a way to reproduce the reactions of the geological environment of interest with greater approximation of reality.
At low temperatures (T < 80°C), it is not possible to observe similarity among the results of the different codes due to the low kinetic energy involved in the diagenetic reactions.Each simulator has its own algorithms based on mathematical codes that are dependent on kinetic and thermodynamic parameters of each mineral, contained in database.
In the physical-chemical approach to higher temperature reactions, the kinetic energy increases to the point to generate similar results, regardless of the model.In the present study, it is shown that at elevated temperatures (100 to 140°C) there is an increase in kinetic energy, incrementing the intensity of the diagenetic reactions.

Minerals
Chemical formula Volume fraction Consequently, the GWB, TOUGHREACT and DISSOL-THERMAL systems generate similar results.However, GWB and THOUGHREACT have different approaches while working with high temperatures when compared to the DISSOL-THERMAL.This can be observed by the behavior of illite and K-feldspar.In general, above 100°C, GWB and TOUGHREACT reproduce similar behavior.The numerical method of both systems uses Figure 3. Geochemical modeling of a sandstone from North Sea using the Geochemist's Workbench and TOUGHREACT software, in comparison with the results obtained by Baccar and Fritz (1993).polynomial equations of 8 th grade to resolve the physical-chemical interactions in geological environments.These equations are dependent on kinetic and thermodynamic data of minerals that form parageneses in sandstone reservoirs.For temperatures above 140°C, these kinetic and thermodynamic data acquire an anomalous behavior for some minerals, indicated by the constant value of 500.If the software finds the 500 value, the previous value is used.This is just a reference approach to indicate that that value has not been calculated.
TOUGHREACT uses equations of state (EOS) to calculate the mass balance.The EOS perform two important functions: 1. to calculate the volume fraction of each mineral at the time of water-rock interaction; 2. to use the dissolution/precipitation of minerals identifying the change of each mineral phase.
From the variation in the minerals volume, it is possible to calculate the porosity and permeability.
The software that is most similar to DISSOL-THERMAL in the condition of this study is TOUGHREACT.This similarity is due to the mathematical treatment as discretization technique (method to solve the differential equations) and the finite difference method, providing flexibility in the description of the geometry of the reservoir.These numerical approximations increase the resolution of water-rock interaction.
Basically, this comparative study revealed two diagenetic stages, taking into account the temperature.The first stage corresponds to the temperature range of 25 -100°C, and the second stage, to the temperature range of 120 -160°C.In the lower temperature stage (25 -100°C), TOUGHREACT and DISSOL-THERMAL simulators reproduce the diagenetic reactions among minerals albite, anorthite and quartz.However, GWB and DISSOL-THERMAL software assume a similarity of geochemical reactions among a larger number of minerals such as illite, kaolinite, K-feldspar and anorthite.In the higher temperature stage (120 -160°C), we find a similar behavior in the 3 systems, due to the convergence of their numerical methods, as it is evident in Fig. 3.
Users must be careful when choosing a geochemical modeling software because, depending on the temperature range, different systems can produce very discrepant behaviors.Although GWB and DISSOL-THERMAL simulators reproduced the largest number of similar diagenetic reactions, users must be cautious in interpreting their data.

CONCLUSIONS
The present study of geochemical modeling comparing the application of GWB, TOUGHREACT and DISSOL-THERMAL codes for the simulation of diagenetic reactions observed in the Snorre Field reservoir sandstones, Norwegian North Sea, resulted in the following conclusions: 1. Diagenetic modeling at low temperatures (25 to 80°C) with different simulators reproduce quite different behaviors.In this temperature range, which involves low energy, the reaction kinetics and the thermodynamic parameters are treated according to the particularities of each code.The numerical methods used by GWB and TOUGHREACT software are able to achieve only limited convergence of results, with low precision expected for diagenetic reactions at low-temperature geological environments.2. For the range of temperature from 100 to 140°C, in which there is higher kinetic energy involved and consequently the diagenetic reactions take place with higher intensity, the simulations with all the three systems generate similar results, because they use the same mass balance and the same temperature variation for the precipitation/dissolution of minerals.The mathematical interpolations used by simulators in this temperature range are more stable and numerically accurate, giving rise to a greater convergence of results.3.For higher temperatures, from 150 to 160°C, GWB and TOUGHREACT reproduced the same behavior of dissolution/precipitation of minerals.However, when compared to DISSOL-THERMAL code, those results are discrepant.In this temperature range, water-rock interaction has a very high kinetic energy, causing fast reactions within the geological environments, making the prediction of the reaction less possible.In this way, each simulator seeks its own best set of numerical solutions to achieve the equilibrium.
Users must be cautious when choosing geochemical modeling software, as an important factor for choosing a simulator is to know the temperature range that will be used.With the temperature information, the user should make a preliminary study with the available codes to verify the convergence of results to the observed diagenetic reactions.

Figure 1 .
Figure 1.Location map (blocks 34/7 and 34/4) of the Snorre Field at the Tampen Spur in the North Sea (from Morad et al. 1990).

Figure 2 .
Figure 2. Burial history and associated stages of diagenesis in the upper Lunde Member sandstones (from Morad et al. 1990).FM: formation.

Table 2 .
Rate constant and equilibrium constant used by each software.