A coupling alternative to reactive transport simulations for long-term prediction of chemical reactions in heterogeneous CO 2 storage systems

Fully coupled, multi-phase reactive transport simulations of CO2 storage systems can be approximated by a simplified one-way coupling of hydrodynamics and reactive chemistry. The main characteristics of such systems, and hypotheses underlying the proposed alternative coupling, are (i) that the presence of CO2 is the only driving force for chemical reactions and (ii) that its migration in the reservoir is only marginally affected by immobilisation due to chemical reactions. In the simplified coupling, the exposure time to CO2 of each element of the hydrodynamic grid is estimated by non-reactive simulations and the reaction path of one single batch geochemical model is applied to each grid element during its exposure time. In heterogeneous settings, analytical scaling relationships provide the dependency of velocity and amount of reactions to porosity and gas saturation. The analysis of TOUGHREACT fully coupled reactive transport simulations of CO2 injection in saline aquifer, inspired to the Ketzin pilot site (Germany), both in homogeneous and heterogeneous settings, confirms that the reaction paths predicted by fully coupled simulations in every element of the grid show a high degree of self-similarity. A threshold value for the minimum concentration of dissolved CO2 considered chemically active is shown to mitigate the effects of the discrepancy between dissolved CO2 migration in non-reactive and fully coupled simulations. In real life, the optimal threshold value is unknown and has to be estimated, e.g. by means of 1-D or 2-D simulations, resulting in an uncertainty ultimately due to the process de-coupling. However, such uncertainty is more than acceptable given that the alternative coupling enables using grids of the order of millions of elements, profiting from much better description of heterogeneous reservoirs at a fraction of the calculation time of fully coupled models.


Introduction
Long-term, reservoir-scale, multi-phase reactive transport simulations in heterogeneous settings are computationally extremely challenging, often forcing to set up oversimplified models if compared to purely hydrodynamic simulations.Typically, 1-D or 2-D models are preferred wherever symmetry allows for it (e.g.Gaus et al., 2005;White et al., 2005;Audigane et al., 2007;Xu et al., 2010;Beyer et al., 2012); only coarse spatial discretisations are adopted for 3-D models (Nghiem et al., 2004;Kühn et al., 2006;Kühn and Günther, 2007;Zheng et al., 2009).As a result, most studies of reactive transport only consider very simple geometries and homogeneous media, thus disregarding spatial heterogeneities at reservoir scale, which in turn are routinely considered by the usually much more detailed geologic models and pure hydrodynamic reservoir simulations.These oversimplifications concern the chemistry as well, leading one to consider only a subset of the potentially reactive minerals.Moreover, the need for sensitivity and uncertainty analysis due to the large amount of uncertain parameters in geochemical simulations (De Lucia et al., 2012;Dethlefsen et al., 2012), and in particular with orders-of-magnitude uncertainties in the determination of field-scale reaction rates (Hellevang and Aagaard, 2013), cannot be met if the computational effort for one single geochemical simulation is practically too high.The resulting models can be used to highlight qualitative results or to provide rough estimations of the ongoing processes at reservoir scale (Gaus et al., 2008, and references therein).To our knowledge, there is no example of fully coupled reactive transport simulations considering complex chemistry on spatial discretisations with a resolution comparable with the one usually adopted for pure hydrodynamic simulations.
Published by Copernicus Publications on behalf of the European Geosciences Union.
However, a number of relevant questions require a rather careful description of reservoir heterogeneity, or depend on a fine resolution of reservoir features, well beyond the possibility of fully coupled simulators and computationally affordable coarse grids.The migration path in a highly heterogeneous system, for example, heavily affects the volume of reservoir which gets exposed to the injected reactant and thus, in the long term, is affected by chemical reactions (Lengler et al., 2010;De Lucia et al., 2011).Hence, heavily upscaled or geometrically simplified models may not be able to satisfactorily capture such migration, adding a decisive imprecision to the coarse coupled simulations.
Especially to solve one of these problems, namely the estimation of the long-term, reservoir-scale mineralisation at the Ketzin pilot site (Martens et al., 2013), a simplified scheme for coupling chemistry and hydrodynamics was introduced (Klein et al., 2013;Kempka et al., 2013b).The purpose of the novel method is not to upscale the simulation grid or apply multi-grid methods, but rather to take advantage of some of the characteristics of the modelled processes, by introducing an approximation in the coupling itself.

Targeted class of systems and formulation of hypotheses for applicability of the one-way coupling
Many potential CO 2 storage systems are located in siliciclastic sandstones, such as in Germany and more generally in northern Europe, in Rotliegend, Bundsandstein and Keuper sequences (De Lucia et al., 2012;Beyer et al., 2012;Klein et al., 2013).Such geological settings are typically characterised by high salinity of the formation brines and limited reactivity of the reservoir rock, given the high proportion of silica and limited amount of faster reacting minerals such as carbonates.Under these circumstances, relevant alteration of the rock matrix due to injection of CO 2 and, more importantly from an engineering point of view, its precipitation in new carbonate minerals are expected to occur hundreds or thousands of years after injection stops (IPCC, 2005).The slow kinetics of the involved reactive chemistry are thereby the principal limiting factor (Marini, 2006;Gaus, 2010).The coupling with hydrodynamic transport of solutes plays therefore only a secondary role, also given the high salinity of the pristine formation fluids, which makes many reactants available for reactions at any given time.In summary, a significantly large class of potential underground CO 2 storage systems is expected to display the following characteristics: 1. Separability of timescales: the timescale of relevant mineral alterations is much larger than that of the hydrodynamic processes: chemical alterations capable of significantly affecting the host rock in terms of porosity and permeability are occurring when the system has reached substantial hydrodynamic equilibrium.Conversely, relevant transport of solutes is happening only at the beginning of the simulation time, when chemistry has not yet become relevant.
2. CO 2 is the only transport-limited reactant: the major driving force for alterations is the presence of the injected CO 2 , either in a separate, dense phase or dissolved in the formation brine.
3. All other reactions are kinetically limited: as a consequence of the previous facts, the presence of other reactants is far less relevant with regard to reactions' outcome.Given the typical high salinity of the targeted formation fluids, the presence and thus the transport of other reactants is not a limiting factor for the reactions.
These three main characteristics are thus considered as assumptions for the application of the one-way coupling.It is important to fully understand the resulting implications.In particular, the first assumption implies that the feedback between chemistry and hydrodynamic transport is limited and therefore negligible, either because quantitatively marginal or because too slow, at the timescale where each of the two processes is most relevant.This encompasses both the concentration changes of the driving-force CO 2 and of the other reactants due to reactions and thus their migration, and the alteration of petrophysical properties of the medium (porosity, permeability).
The assumptions are of course a simplification and will not hold for all possible CO 2 storage sites.In particular, not only carbonate rocks but also siliciclastic reservoir with relevant carbonate cementation cannot satisfy them, since significant amounts of reactions are expected at the same (short) timescale as the hydrodynamic processes, making it impossible to un-couple them.Another consequence of the assumptions is that chemical processes, such as salt precipitation, do not fall in their applicability domain; namely this is a transport-limited and fast process which violates the first and third assumption.However, such processes are generally expected only in the the immediate proximity of the injection wells; furthermore, a hypothetical loss of injectivity, while being a problem during the operational life of the storage site, will not significantly affect the long-term chemistry at reservoir scale, which is the scale of the models which are considered in this study.For these reasons, dry-out in the vicinity of injection wells can be be safely disregarded.Nevertheless, a site-specific assessment is recommended in order to verify the adherence the investigated system to the assumptions.We can anticipate here that none of the reactive simulations for the case study of this paper showed possible salt precipitation, but this could be also due to the coarse discretisation used; however, in the specific case of the Ketzin site, there was no sign of injectivity loss throughout operations, and repeated Pulsed Neutron Gamma (PNG) logging campaigns interpreted very moderate, spatially concentrated and probably transient salt precipitation around the borehole (Baumann et al., 2014), thus confirming its negligibility.
It is also important to note that the considerations and methods proposed in this study are not strictly limited to CO 2 storage systems, but can be employed wherever the system behaves in accordance with the three assumptions, just substituting the CO 2 with a limited amount of transport-limited reactants controlling the geochemical reactions.Therefore, it will be demonstrated in the following that the three assumptions above are sufficient conditions for establishing the simplified coupling.
The paper is organised as follows: the proposed one-way coupling approach is described in detail in Sect.2.1.Section 3 introduces the case study, with the reference geochemical model, the simulation grid and the general set-up of the non-reactive and fully coupled reactive transport simulations performed for the validation.Their respective results are analysed in Sect. 4 aiming at verifying the major hypotheses.Finally, a comparison of the outcomes of simplified and full coupling concludes the study.

Implementation
The concept of the one-way coupling was firstly introduced by Klein et al. (2013).The idea is to combine, in a postprocessing approach, independent batch geochemical simulations with non-reactive flow simulations, while applying some criterion to at least ensure partial mass balance and thus compensate for the lack of feedback between chemistry and hydrodynamics.Hereby, the non-reactive simulations are intended to calculate the mutual solubility of CO 2 in formation brine and conversely of H 2 O into the CO 2 -rich phase, which are actually chemical processes and in particular dependent on the ionic strength of the brines.In particular, the proposed approach can be considered as a radical operator splitting, where the chemistry operator is applied at the end of simulation time instead of after each time step of the transport operator.
In detail, the procedure is constituted by the following steps: 1.The non-reactive hydrodynamic simulations should be rub for the whole simulation time.
2. The reference 0-D batch geochemical simulation representative should be run for the whole reservoir and for the whole simulation time.
3. The arrival time and residence time of the injected CO 2 in each grid element should be derived from the nonreactive simulations, distinguishing between gaseous or dissolved CO 2 ; we will refer to this time interval as the exposure time to CO 2 .If the element is exposed only to dissolved CO 2 , its concentration has to exceed a threshold to be considered exposed to CO 2 : the arrival time is therefore taken as the simulation time where the concentration first trespasses such a threshold.

4.
During the exposure time in each grid element, the gas saturation or the concentration of dissolved CO 2 can dynamically vary.We define the characteristic gas saturation and concentration of dissolved CO 2 as the respective maximum value reached during the exposure time.
5. The reference batch geochemical model is analytically scaled to each element of the grid for the length of its exposure time.The scaling is needed to consider the element's porosity and the characteristic gas saturation.Insights about the analytical relationship are given in Sect.2.2.
6.If the element is exposed only to dissolved CO 2 , the exposure time is truncated to prevent the CO 2 mineralisation to exceed the characteristic concentration reached in the element.
7. If the element is exposed to gaseous CO 2 , it is assumed that enough CO 2 reaches the element to achieve the maximum mineralisation without further need for a mass balance.
8. The chemical state of each element at the end of its exposure time is prolonged to the end of the simulation time.
9. The mineralisation in each element is finally summed in time to achieve global mineralisation in the reservoir.
The lack of feedback between hydrodynamics and chemistry is of course the major drawback of the simplified, volume-based one-way coupling, since it leads to lack of mass balance between chemistry and transport.A control regarding the mass balance in the elements exposed only to dissolved CO 2 has been therefore built into the algorithm.Such bounding is considered not necessary if a significant saturation of gaseous CO 2 is reached in the element, since, given the hypothesis of small amount of reactions, it is assumed that in such a case enough CO 2 is available to reach full mineralisation for the given time frame.
The simplified coupling builds on top of non-reactive hydrodynamic simulations, where indeed no CO 2 is immobilised in precipitating minerals.As a consequence, the spatial extent of the CO 2 plume calculated by those hydrodynamic simulations is always equal to or larger than the one predicted if chemistry is considered, namely by the fully coupled simulations.
It is therefore of utter importance that the modelled system satisfies the first assumption of Sect. 1.In the following, through comparison with fully coupled simulation, it will be demonstrated that the discrepancy of the CO 2 migration predicted by the non-reactive simulations is not too large with respect to the case with active chemistry.Furthermore, it will be shown how the proposed method of defining a threshold on the dissolved CO 2 concentration allowed to activate the chemistry mitigates this issue to a large extent.This constitutes the first main hypothesis that needs to be verified to validate the simplified coupling against fully coupled simulations.
For an initially geochemically homogeneous scenario, which means initial concentrations of solutes and volume fraction of minerals, the simplified coupling foresees the use of one single batch geochemical model which is practically replicated onto all elements of the simulation grid, possibly scaled following the local porosity and gas / water saturation ratio (cf.next session for more details on this subject).In other words, the simplified coupling assumes that in each element of the simulation grid the reaction path is qualitatively and quantitatively similar to one single underlying model: we will refer to this property in the following as self-similarity of reactions in the reservoir, and it descends from the combination of the first and third assumptions above: negligible feedback of transport to chemistry, and strong kinetic control of the reactions with CO 2 as only transport-limited reactant.This is a quite strong conjecture, and constitutes the second major hypothesis that needs to be verified by comparison with fully coupled reactive transport simulations.
The gas content and and the dissolved concentration vary during the simulation and also in particular during the exposure time of grid elements to the injected CO 2 .By comparison with fully coupled simulations, it can be shown that the maximum gas saturation during exposure time controls the speed of reactions; the maximum concentration of dissolved CO 2 controls the total amount of mineralisation in the elements where no gaseous CO 2 arrives.Therefore, the respective maximum value is retained as characteristic for the two variables, in the sense that it controls the outcome of the chemistry.Other tested possibilities were the average, timeweighted average and the median values during the exposure time.
Lastly, the medium considered throughout this study is homogeneous in terms of volume fractions of reactive minerals: heterogeneity is in the present study only referred to porosity and permeability.In presence of different facies or a spatially heterogeneous distribution of volume fractions, the method can still in principle be applied using as many batch geochemical simulations as there are different facies in the reservoir, assigning the corresponding chemistry to each element of the hydrodynamic grid.However this possibility has not been explored in the present study.

The analytical scaling relationship
One major feature integrated in the proposed simplified coupling approach is the treatment of spatial heterogeneity of porosity and of gas saturation in the porous medium, the latter also has a dynamical variation along the simulation run, but ultimately taken into account as a constant, characteristic value for each element during its exposure time.Since all equations in geochemistry are expressed per unit mass of solution or solvent, if one wants to actually calculate the total reaction in a given volume of reservoir rock, a scaling is needed considering the total water present in the rock volume, the total reactive surface of each mineral available for reactions and so on.Under assumption of homogeneous mineral volume fractions, all these variables are ultimately a function of porosity and gas saturation, which hence control speed and total amount of the chemical reactions in a rock volume.Klein et al. (2013) introduced analytical relationships relating the outcome of a single geochemical simulation, scaled for a given porosity and gas saturation, to any other given porosity and saturation, enabling the treatment of heterogeneity without running many different batch geochemical models.In the following it will be demonstrated that the equations are a direct consequence of the particular form of the kinetic law implemented in the models, and of the choices made representing the minerals' reactive surface and possibly its variation along the simulation.
Considering the reference rock volume V r , a common rate expression per mass unit of solvent (with dimensions of [mol 1+γ s −1 kg H 2 O −1−γ ]) takes the form (Lasaga, 1998) where is the specific (per mass unit of solvent) reactive surface of the mineral in contact with the solution, the saturation ratio of the mineral and α, β two parameters reflecting deviations from the ideal transition state theory, to be empirically determined for the mineral; a is the activity of a solute species acting as catalyser and its exponent γ defines the order of the kinetic law.The specific reactive surface can be written as where V m [m 3 mineral kg H 2 O −1 ] is the volume of the mineral in contact with the unit mass of solvent in the considered system and A m [m 2 m −3 mineral] the specific reactive surface per unit of mineral volume.A m is a parameter intrinsic to the mineral -e.g.related to grain size.The expression of V m in terms of porosity and water saturation for the reference rock volume V r reads where ϕ the porosity, f m the volume fraction of mineral m referred to all minerals and ρ H 2 O the density (mass of water per unit volume of solution).The presence of S w at the denominator of this expression means that care must be given to very low values of water saturation.In practice, the validity of this parametrisation is limited to water saturation being larger than a given threshold (few percents).For lower values of water saturation, chemistry cannot be calculated; however, its disregard would be irrelevant given the tiny amount of reactions happening in presence of tiny amount of solvent.The rate r, including its dependency on porosity and water saturation, reads as Consider now two distinct rock volumes, V 1 and V 2 , with initially equal mineral volume fractions and chemical composition (concentrations) of the reactive solution and differing only per porosity and saturation respectively ϕ 1 , S 1 and ϕ 2 , S 2 .The rates as seen by the volume unit of solution are (5) It is not necessary to explicitly integrate these equations to compute their respective timescales.It is namely straightforward to recognize that the chemical reaction will proceed in the two volumes with different velocities, and thus different characteristic timescales, but going through a succession of equal states.For this purpose it is illustrative to proceed by discretising Eq. ( 6) considering the amount of reaction per mass unit of solvent ξ [mol kg H 2 O −1 ], which is often called progress of reaction.We have Sought are the times t 1 and t 2 at which the reactions in V 1 and V 2 reach the same progress, or ξ 1 (t 1 ) = ξ 2 (t 2 ).Consider time steps t 1 = t 1 − t 0 and t 2 = t 2 − t 0 small enough that the rates r 1 and r 2 can be considered constant and the changes in water saturation and porosity become negligible.One can write Substituting Eq. ( 5) into Eq.( 7) and simplifying for the f m , ρ H 2 O and the terms, which by definition are equal for both V 1 and V 2 at time t 0 , we finally obtain Upon reaching t 1 and t 2 respectively, the solutions in V 1 and V 2 are again completely equivalent, and so are the (1 − ) terms, which only depend on activity coefficients and aqueous concentrations.Hence, the procedure can be repeated, meaning that for any given time Furthermore, if ξ is the amount of reaction per unit mass of solvent, the total amount of reaction N in the given rock volume is Dividing the two equations and rearranging, it finally reads Equations ( 9) and ( 12) are the scaling equations introduced by Klein et al. (2013), which are given here in a slightly more general form.
Notably, this result holds also for kinetic laws where the total rate is a linear combination of terms of the same form of Eq. ( 1), and for mixed kinetics/equilibrium simulation as well, since equilibrium is a special case of kinetics where the kinetic constant k is very large.One major assumption still has to be respected, i.e. that the reactions do not significantly affect the porosity ϕ and the water saturation S. Furthermore, in this parametrisation porosity and water content control the timescale of chemical reactions.Analogous expressions can be obtained changing the form of kinetic laws and the parametrisation of the reactive surfaces.

Scaling relationship for the TOUGHREACT simulator
The TOUGHREACT simulator (Xu et al., 2011) implements the calculation of reactive surfaces canceling out completely its dependence on water saturation, in order to prevent diverging rates when water saturation is small (Xu et al., 2008, Appendix G, p. 175).Hence, the scaling equations for this simulator read 3 Validation of the simplified scheme The procedure for validating the one-way coupling involves the comparison and analysis of a case study for which both the simplified coupling approach and the fully coupled reactive transport simulation are applied.
The case study has been chosen in order to ensure typical conditions for the application range of the simplified coupling, and is inspired by published models for the Ketzin pilot site (Kempka and Kühn, 2013;Klein et al., 2013).For all simulations the established TOUGHREACT/ECO2N (Pruess and Spycher, 2007;Xu et al., 2011) simulator was adopted.It implements a solubility model for CO 2 in brines which is accurate until 100 • C and 600 bar pressure for a Na / Cl dominated solution up to 6 molal (Spycher and Pruess, 2005).The non-reactive simulations have been obtained running the same model and switching off chemistry, and the reactive in sequential non-iterative mode.
It is worth mentioning that the simplified coupling is a general method and is independent from the simulator used for the validation.
The influence of spatial heterogeneity, here limited to porosity and permeability, and of the characteristic timescale of reactions relative to the characteristic timescale of CO 2 migration were also investigated.

Geochemical model
The chemical system chosen for the validation reflects the analyses published by Förster et al. (2010), Norden et al. (2010) and Würdemann et al. (2010) and is based on the model for the Ketzin test site published by Klein et al. (2013).The chosen model displays a typical complexity for CO 2 storage problems; there is no need and no claim of being an accurate quantitative reference for the site, nor to investigate its plausibility or the quality of the thermodynamic database adopted for the simulations.
The Stuttgart formation is the target reservoir of the Ketzin pilot site.It is mainly constituted by sand channels and floodplain facies.Reservoir simulations (Kempka and Kühn, 2013) show that the sandy facies receives the major amount of the injected CO 2 due to their higher conductivity, and can therefore be retained as the reference facies for investigating the chemical processes.
The rock-forming phases considered reactive are anhydrite, albite, hematite, illite and chlorite, together amounting to 24 % in volume of the rock.The rest of the rock -principally quartz, up to 40 % in volume -is considered inert.Secondary phases included in the model are kaolinite and the carbonates siderite, magnesite, dolomite and calcite.The choice of endmembers for the carbonates is a proxy for a Fe-Mg-Ca-carbonate solid solution which is more likely to form under these conditions; dolomite has also been included as proxy for a Mg-Ca solid solution.These simplifications are principally dictated by numerical reasons and constitute in our opinion an acceptable approximation for the complexity of the chemical system, also given the purpose of this model.The cement phase anhydrite and the secondary minerals calcite and dolomite are considered at local equilibrium in the model, whereas all other minerals follow the parametrisation for kinetic laws given by Palandri and Kharaka (2004) and the estimation of specific surface areas by Klein et al. (2013).Notably, 1 order of magnitude is preserved between the specific reactive surface areas of clay minerals (kaolinite, chlorite and illite) and the non-clays.For the initial brine composition, the available analyses of Würdemann et al. (2010) have been considered without further refinement other than the concentrations of SiO 2 and AlO − 2 which were taken from the re-evaluation performed by Klein et al. (2013).The initial composition is listed in Table 1.The batch model, also simulated with TOUGHREACT, is defined as a porous medium of 1 m 3 volume, represented as one single cubic element in the program.The minerals' volume fractions listed in Table 2 are assigned to the medium.The rest of pore space not occupied by the brine is assumed filled with CO 2 at an initial pressure of 70 bar, which is an average value derived by the analysis of the hydrodynamic simulations of Kempka and Kühn (2013).The impact of pore pressure was investigated by Klein et al. (2013) by means of batch models with constant pressure, and found negligible with respect to the speed of reactions at least in the range of pressure expected in Ketzin.In the present case the system is closed, and thus the simulation is performed at constant volume for the gas phase.The change in porosity and the consumption of water due to reactions, as well as the solubility of H 2 O in the gas phase, are negligible in terms of volume change.Hence, the initial pressure governs how much CO 2 is contained in the reference volume and therefore available for reactions.The reference simulation is for this reason chosen with an initial porosity of ϕ = 0.5 and a water saturation of S w = 0.5.
Figure 1 displays the outcome of the reference batch simulation.Anhydrite (not shown in the picture), is undersaturated at the beginning of the simulation and, being included as equilibrium phase in the model, is from the very first step of the simulation dissolved, staying substantially constant afterwards.Chlorite is dissolved and initially also albite.However, the latter inverts its trend and starts re-precipitating after around 3000 simulation years.Illite starts precipitating from the beginning, reaching a plateau after around 3000 years and staying constant afterwards.At this point also kaolinite starts precipitating.The most relevant feature of the model is the sequence of precipitating carbonate minerals.The first, starting from ca. 1500 simulation years, is the iron-bearing siderite, which is followed by a temporary appearance of calcite.Calcite starts redissolving as dolomite starts precipitating.Precipitation of siderite and dolomite continues until the end of simulation time.There is no appearance of magnesite in the model and a negligible dissolution of hematite.
Prolonging the simulation to around 200 000 years sees the model approaching but not reaching an equilibrium state.At this point the precipitated volume of siderite reaches around 3000 cm 3 m −3 of rock, or 0.3 % of the total rock volume.The precipitation of carbonates and illite, however, is partially compensated by dissolution of chlorite and anhydrite, and thus the calculated change in porosity after 200 000 years is −0.00056 porosity units, which represents a relative change of around −0.11 % from the initial value of 0.5.It can be safely considered negligible in terms of impact on the hydrodynamic properties of the medium.

Coupled simulations: spatial discretisation and CO 2 injection
The simulation grid chosen for the coupled simulations in this study displays the typical complexity -and limitationsof reactive transport simulations in the domain of CO 2 storage, being coarse enough to allow for fully coupled reactive transport simulations in an affordable CPU time.
The spatial discretisation of about 2.5 km × 3.5 km nearfield subset of the Stuttgart formation of the Ketzin anticline is comprised of 2950 hexahedral elements arranged in one single layer with a 3-D structure following the topography of the target formation (Fig. 2).The thickness of the layer is around 75 m and thus is comprised of both main stratigraphic units of the Stuttgart formation.In particular, the sloping of the anticline structure drives a significant upwards migration of the CO 2 and was therefore retained as a major distinctive feature the simplified coupling needs to be benchmarked against.Faults and discontinuities defined in the geological model have been geometrically smoothed in the one-layer model.
Pore volume multiplicators at the border elements provide constant boundary conditions throughout the simula- tions.The construction model ensures a realistic pressure build-up and is large enough that the injected CO 2 does not arrive at the borders in the investigated time span.Two distinct cases were investigated: an homogeneous medium with constant porosity ϕ = 0.2 and permeability of 1 × 10 −13 m 2 , and a spatially heterogeneous medium with porosity ranging from 0.08 to 0.22 and permeability from 2.8 × 10 −14 to 3.0 × 10 −12 m 2 .The heterogeneous distributions of porosity and permeability have been obtained by upscaling to the simulation grid respectively with arithmetic and geometric average the reservoir model described by Norden and Frykman (2013) and Kempka et al. (2013a) and history matched by Kempka and Kühn (2013).No further refinement of the upscaling was done, i.e. to obtain the same total pore volume in the simulation grids between the homogeneous and heterogeneous cases, or more sophisticated permeability upscaling.
The initial pressure for the models is set after equilibration runs assuming a hydrostatic gradient.All simulations are isothermal at a temperature of 35 • C. 67 000 t of CO 2 are injected at a constant rate in 5 years in the element highlighted in Fig. 2.This represents the amount of CO 2 injected at the Ketzin site.After injection, the system is let free to evolve until a total simulation time of 2200 years is reached.
It has to be stressed that the many adopted simplifications do not allow one to consider the resulting model as an attempt to investigate the Ketzin site in a way coherent with the available monitoring data and previous modelling work.Again, the purpose of these numerical experiments is to validate the simplified coupling approach while mimicking the complexity (in this case: heterogeneity of porosity and permeability) of a real-life problem.

CO 2 migration in reactive and non-reactive simulations
For the purpose of this study it is important to evaluate the discrepancy between the reactive and fully coupled nonreactive simulations concerning the migration pattern ofand hence the reservoir volume exposed to -the injected CO 2 .The simplified coupling is applied on the outcome of non-reactive simulations, and it must therefore be ensured that the bias with respect to the fully coupled simulations is quantitatively manageable.The reservoir volume exposed to the injected CO 2 , either gaseous or dissolved, is therefore one of the crucial parameters in view of the application of the simplified coupling.After injection, the CO 2 would rapidly start its upward migration towards the anticline top, spreading and progressively dissolving in the formation brine along its way.Due to the coarseness of the simulation grid, the gaseous phase tends to progressively disappear, so that at the end of simulation time the majority of elements in which gas is still present display only a small residual gas saturation and almost the whole injected CO 2 is dissolved in the formation brine.The spatial distribution of dissolved CO 2 is displayed in Fig. 3, which collates snapshots after 2100 simulation years for both the homogenous and heterogeneous cases.The significant differences in migration in the homogeneous and heterogeneous cases stem on one hand from the difference in total pore volume around the injection between the two simulations, and on the other hand from the preferential flow paths which enhance the total spreading of the CO 2 (Lengler et al., 2010).
The non-reactive simulations show a quite enhanced spreading of the dissolved CO 2 compared to the reactive case.However, this mainly consists of elements with a small concentration located at the border of the CO 2 plume.Its central part, where the highest concentrations are, shows a similar spatial pattern.This is particularly evident in the homogeneous case, where the cloud has a much lower expansion than in the heterogeneous case.
Figures 4 and 5 show the quantification of the exposed volumes respectively for the the homogeneous and the heterogeneous cases.In the diagrams the reservoir volume exposed to gaseous CO 2 (the black curves) is compared with the reservoir volume exposed only to dissolved CO 2 in both fully coupled and non-reactive simulations (blue and red curves).It appears evident that soon after injection the reservoir elements showing a gaseous phase reach a maximum and start decreasing immediately afterwards.However, the discrepancy between reactive and non-reactive simulations is quite small in both cases: the presence of gas phase in an element is not significantly overestimated by the non-reactive simulations.
Much more significant is the difference between the volumes exposed to dissolved CO 2 .Hereby two threshold values used to include the elements in the statistic are highlighted in blue and red, pointing out the large sensitivity of the outcomes on this choice.This is due to the fact that many elements are actually exposed only to low CO 2 concentrations.Furthermore, a major discrepancy is found, both in the homogeneous and heterogeneous case, between the fully coupled and the non-reactive simulations: the latter curves display a constant increase of exposed volume, whereas the reactive reach a maximum after around 500 years for the homogeneous case and around 1000 for the heterogeneous, afterwards substantially decreasing or hitting a plateau.
This behaviour actually means that in the non-reactive simulations the CO 2 plume continues increasing its extent during the whole simulation time.Such discrepancy is actually a measure of the overestimation of the migration pattern caused by not accounting for chemical reactions.
However, the importance of the threshold is spatial and has an impact on the definition of the considered exposure time.The larger the threshold is, the later the arrival time of dissolved CO 2 considered.In the simplified coupling the total mineralisation of the elements exposed only to dissolved CO 2 is limited on one hand to the actual maximum amount of CO 2 reaching the element, but on the other it also gets limited by the arrival time of CO 2 , which is retained as the starting point of the chemical reactions.Therefore, the significance of the threshold is 2-fold, and it is the combination of its spatial At this simulation time there is a significant difference between the reservoir volume affected by dissolved CO 2 depending on whether or not chemistry is included in the simulation.Since the simplified coupling makes use of non-reactive simulations to estimate the exposure time of each element, this discrepancy will impact its outcome.However, this effect is mitigated by the introduction of a threshold filtering out small concentrations of aqueous CO 2 as geochemically inactive (see text for further details).and temporal effects which helps mitigate the overestimation of CO 2 migration returned by the non-reactive simulations.In summary, the reservoir volume exposed to gaseous CO 2 is accurately predicted by the non-reactive simulations and does not need a particular treatment.On the contrary, for the elements exposed only to dissolved CO 2 , the discrepancy between the non-reactive and fully coupled reactive simulations is significant and imposes the choice of a threshold value for the CO 2 concentration considered geochemically active to mitigate it.

Spatial self-similarity of reactions
The second major hypothesis is the self-similarity of chemical reactions in the fully coupled simulations.In fact, the hypothesis underlying the simplified coupling is that the same reaction path is replicated in every element of the simulation rock − or pH units q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q Max S g : 0.07 − ϕ : 0.2 q scaled 0D simulation element of fully coupled simulation pH 6.Comparison between the mineral changes predicted by the fully coupled simulations in one single element as opposed to the outcome of the scaled 0-D reference geochemical model.Depicted, with the same colours of Fig. 1, are siderite (red), albite (blue), chlorite (magenta) and calcite (black).Also, for completeness the pH in the two cases is represented (dashed line for the scaled batch and the cross symbols for the fully coupled simulation).The calculated similarity function for this element is 4.5.grid, translated in time by the corresponding arrival time of CO 2 and scaled according to the actual porosity and characteristic gas saturation of the element.Such a replicated reaction path happens to be exactly the reference 0-D batch geochemical model depicted in Fig. 1.The discrepancies among the elements of the coupled simulations can be imputed to the hydrodynamic transport of solute species, which is the only physical process responsible for perturbing the reactions, since by definition they all start from the same initial state in the whole domain.
The comparison is visually represented in Figs. 6 and 7, which show two exemplary elements taken respectively from the homogeneous and heterogeneous simulation.The first is an element exposed to gaseous CO 2 , while the second is exposed only to dissolved CO 2 .The dotted lines represent the outcome of the fully coupled simulation and the solid lines the output of the (scaled) 0-D model for one element of the simulation grid.The time axes of the scaled 0-D simulation have been translated to match the actual arrival time of CO 2 in that element.
We can anticipate that from a quantitative point of view, the agreement between the fully coupled and the 0-D scaled simulations is generally more than acceptable in the majority of grid elements, which supports the hypothesis the simplified coupling is based on.Common discrepancies are represented by a slight delay in precipitation of some minerals (i.e.calcite in Fig. 6, the black line).rock − or pH units q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q q Max S g : 0 − ϕ : 0.16 q scaled 0D simulation element of fully coupled simulation pH CO 2 arrival Figure 7.Comparison between the mineral changes in one element of the heterogeneous case with weak similarity (52.26).This element has only been exposed to dissolved CO 2 .Even with such a high value of the similarity function, the difference between the precipitated carbonates at the end of simulation time is below 10 %.
A numerical measure of the self-similarity is needed in order to control the initial hypothesis.This is achieved comparing the mineral changes in each element to the 0-D reference simulation.Noted with M i , the changes of a mineral M (in mass unit per m 3 rock) at the time steps i = 1, . .., N from the fully coupled simulations and with M * i the changes from the reference model (after translation of its time axis to match the arrival time of CO 2 in the element), the similarity is defined as the average quadratic relative error over the time samples: From such a definition it follows that the similarity is a positive function, equal to zero if the compared reactions are exactly the same, and taking higher values for increasing discrepancy.Similarities can be computed for each mineral separately or for one of their possible linear combinations.In the following we will concentrate on the most relevant linear combination for the purpose of this study, which is the total amount of precipitating carbonates.
The results confirm that indeed the hypothesis of selfsimilarity holds: in Fig. 8 are shown the histograms of the calculated similarities for each element of the fully coupled simulations reached by CO 2 in the homogeneous and heterogeneous case respectively.For readability the histograms have been truncated to a value of 20 for the similarity function, which represents the 95th percentile in the homogeneous and the 99th percentile in the heterogeneous case.In both cases the vast majority of values are concentrated around 5 or less, with particularly small values -and thus high a degree of similarity -for the heterogeneous simulations.Figures 6 and 7 display similarities of 4.5 and 52.26 respectively.Even for such a high value as in the latter case, the discrepancy of total amount of mineralised CO 2 in the element is around 10 % at the end of the simulation time: in that particular case the delay in the mineralisation of siderite and the temporary appearance of calcite explain this difference.
Figures 9 and 10 display the spatial distribution of the similarity function in the reservoir, for the homogeneous and heterogenous cases.It appears that the best (lowest values of the similarity function, in blue) similarities are located, in both cases, in the extern part of the CO 2 cloud, whereas the core part of the cloud shows lower degree of similarity.On the one hand this is due to the longer exposure time and hence duration of reactive chemistry in such portion of reservoir, which amplifies and propagates the discrepancies; on the other hand, this is an indirect measure of the influence of hydrodynamic transport: it is in the central portion of the cloud that the brine flow is the highest and thus most perturbs the chemistry.
In summary, for at least 95 % of the elements in which CO 2 mineralisation occurs, the fully coupled simulations predict reactive chemistry which is qualitatively and quantitatively nearly identical to the reference 0-D simulation.The second major hypothesis for the application of the simplified coupling can be thus considered as verified.

Total mineralisation
The next step is to evaluate the actual performance of the simplified coupling with respect to the outcome this method has been built for: the prediction of the total CO 2 mineralisation in the reservoir.
This result is shown in Fig. 11, highlighting the influence of the choice of the threshold on the result.
The diagram confirms that indeed the results of the simplified coupling (the dotted and the dashed lines) are quite on par with the outcome of the full coupling (the solid lines) for both the homogeneous and heterogeneous cases.It is also apparent that the chosen threshold value has a significant impact on the predictions.For the simulations of this study the optimal concentration threshold, found by trial and error, is comprised between 0.0002 and 0.0005 in mass fraction.The evolution of the mineralisation is nicely reproduced throughout the simulation time, with the exception of a slightly more significant deviation (in form of delay) in the first 500 simulation years for the homogeneous simulation and around 1000 years in the heterogeneous.After that, the simplified coupling approximates with excellent agreement the prediction of fully coupled simulations: under 3 % in terms of the total injected CO 2 .
Based on Fig. 11, a further consideration is possible, in particular from the observation of the heterogeneous case (the blue curves).As result from the enhanced spreading of the injected CO 2 and the heterogeneity in porosity, contributing to the speed of reactions, the mineralisation is quite quick if compared with the homogeneous case, reaching a value of around 60 % of the injected CO 2 after 1700 years.Under these conditions, a further prolongation of the simulation cannot be justified: the discrepancy in spreading and migration of CO 2 as estimated by the non-reactive simulations, which are of course not accounting for in the feedback of chemistry, will be obviously overestimated to an excessive degree, violating one the assumptions the simplified coupling is based on.Therefore, Fig. 11 is bounded on the y axis: the outcome of the comparison is not credible if the mineralisation is too high.
The influence of the kinetics relative to the transport scenario was furthermore tested by replacing the chemical system considered until now with one with kinetic constants 1 order of magnitude smaller.The comparison procedure was repeated in an homogeneous case, running the new 0-D batch geochemical model and the fully coupled reactive transport, while the non-reactive simulations are the same as in the previous case.The considered simulation time for this case is 10 000 years.The outcome in terms of total mineralisation is shown in Fig. 12.
Here one would expect that the slower kinetics would produce lower mineralisation, and thus that the difference between reactive and non-reactive simulations would be reduced; overall, a better agreement between full coupling and simplified coupling would be anticipated.However, surprisingly, the discrepancy between the two approaches proved to be larger.Only a very small concentration threshold allows one to match the predictions of the full coupling with the same precision of the case of faster kinetics.The explanation of this result is the decrease in self-similarity of reactions: in other words, for slower kinetics the relative importance of hydrodynamic transport of solutes becomes more relevant, and thus one single geochemical model is less capable of representing the reactions in all elements of the spatial discretisation.A faster kinetics would reinforce the hypothesis that only the presence of CO 2 in each element drives its chemistry, reducing the relative importance of hydrodynamic transport.In this case it is more likely to experience more severe errors in the mass balance for the injected CO 2 .Continuing with this reasoning, an even faster kinetic would eventually lead to violating the third underlying assumption of Sect. 1 in which the reactions are largely kinetically limited.These ideas suggest the existence of a window of applicability for the simplified coupling and ultimately a more rigorous way to define (and compute a priori) if and when transport and reactions can be un-coupled for a given system, and when fully coupled reactive transport is the only possible option.
In summary, the simplified coupling achieves an estimation of the total amount of mineralisation which is in excel-lent agreement with the fully coupled simulations if an optimal value for the threshold is employed.

The optimal choice for the threshold
A threshold value for the minimum dissolved CO 2 concentration considered active has been defined to filter out the overestimation of reservoir volume exposed to CO 2 (in both a spatial and a temporal sense) which comes from the nonreactive simulations.If too large, its choice can lead to severe underestimation, if too small, to severe overestimation of the total mineralisation.
In a real application for the simplified coupling such optimal value is of course a priori unknown, since no fully coupled simulations are available.The strategy to determine realistic values -or at least to define a realistic bandwidthis always dependent on the investigated problem, and only guidelines can be outlined here.
A first consideration regards the proportion of elements exposed to dissolved CO 2 in the simulations.In a typical storage system in saline aquifers, the reservoir volume exposed exclusively to dissolved CO 2 can be much larger, in the time frame for which chemical reactions are relevant, than that exposed to gaseous CO 2 (Kempka and Kühn, 2013).However, the simulated concentration of dissolved CO 2 is often quite small, particularly at the front of the CO 2 plume.
Coarse spatial discretisations resulting in elements with very large volumes may induce an overestimation of the spatial extent of dissolved CO 2 .In fact, most simulators allow the existence of gaseous CO 2 in a separate, dense phase only if the dissolved CO 2 excesses the saturation limit in the formation water in that element.Therefore, in particular in presence of large volume elements, and particularly at the border of the CO 2 plume, the presence (and hence propagation) of dissolved CO 2 is favoured with respect to the gaseous phase, ending up in overestimating the total dissolved CO 2 .Very small calculated concentrations of dissolved CO 2 propagate to adjacent elements to an extent which is not completely physical.Numerical dispersion may contribute to this effect.Notably, all these issues become more relevant for increasing coarseness of the simulation grid.Therefore, it can be expected that the fine hydrodynamic grids, the simplified coupling is designed to be applied on, should be far less sensitive to the threshold than the coarse model used in the present study for validation.
One-dimensional tests not showed in the present work support this hypothesis: through grid refinement not only the choice of threshold became less significant, but also the optimal threshold value was increased in the finer grids.Future work is however needed to thoroughly assess this aspect.Another interesting option to explore in future work is the def-inition of a time-dependent threshold, and in particular increasing along the progress of simulations.
A computational way to determine the optimal threshold is represented by the set-up of 1-D or 2-D models mimicking the typical element dimension found in the full 3-D grid, the typical exposure time to CO 2 , and of course the same chemistry as in the original model.

Applications and future work
The simplified coupling has been developed to apply reactive chemistry on structurally complex, fine resolved, millionelement grids routinely used in pure hydrodynamic reservoir simulations, with no need to meet compromises such as the number of elements and the spatial heterogeneity described in the models.
This approach frees the reactive transport modelling from the limitations which hamper many of its real-life applications.It applies to hydrodynamic reservoir models, which are usually well matched, without grid upscaling and therefore profiting of the best possible description of heterogeneity.This also includes boundary conditions and transitory regimes (think about injection rate of CO 2 and its consequences on reservoir pressure), which are problematic to depict using oversimplified grids.The price for this abilityif all the hypotheses are met, which is usually the case in CO 2 storage systems -is not a computational or numerical burden, but an approximation of the total outcome of reactive chemistry, determined as bandwidth following the uncertainty of only one parameter: the concentration threshold.Kempka et al. (2014) presents the application of the simplified coupling in the form described in the present work to a real-life scenario concerning the Ketzin pilot site.The simulation grid used for the hydrodynamic model amounts to 648 420 elements and the simulation time reaches 16 000 years; it would have been impossible to achieve these results with the currently available fully coupled reactive transport simulators.
The storage in saline aquifer is actually an unfavourable case for the one-way coupling, given the complexity of its hydrodynamics; the case of CO 2 storage in depleted gas reservoirs would be much easier, where the water phase is typically present only in limited amount or residual saturations (Audigane et al., 2009;De Lucia et al., 2012).In such cases the hypothesis of self-similarity of reactions will be met to an even better extent, since the solute transport in the formation brine will be quite reduced, if present at all.It will be interesting, as future work, to test the one-way coupling in such settings, also given the higher temperatures which could lead to faster reactions than those investigated in this study.Also, the scale of the injection project, considering millions of tons instead of few tens of thousands as in the Ketzin example, could significantly affect the applicability and the precision of the one-way coupling.In this regard we expect that the larger the CO 2 injection relative to the avail- able pore volume, the better the precision of the simplified coupling, since even more reservoir volume would be exposed to gaseous CO 2 for longer time than in the case study of this paper.This is however a speculative assessment, since the hydrodynamics of the reservoir could correspondingly be much more complex than anticipated here.The simplified approach is not limited to CO 2 or gas storage, but can be employed wherever its underlying assumptions are respected, and notably where the reactions are initiated or controlled by the arrival time or presence of a limited number of solute species.In this sense the presented method can be truly considered as an alternative to fully coupled reactive transport simulations and not the simplification of one particular CO 2 storage model.As conceptual proof for this flexibility, one could have considered the pH value as control variable for the initiation of chemical reactions in the hereby presented CO 2 case study, and led the whole study in terms of pH, which could lead to a much more physically meaningful definition of the threshold value.However, since the explicit evaluation of gaseous CO 2 saturation (for the scaling) and dissolved concentration (for the mass balance) would still be needed in our case, we did not follow that route.
A crucial application of the method is furthermore the ability to perform sensitivity and uncertainty analysis of the underlying components: both the hydrodynamic and the geochemical model.The advantage is of course represented by the separation of the processes, which can be therefore efficiently simulated on their own, and quickly reassembled in the post-processing approach, without the computational burden of fully coupled simulations.
Finally, a combination of full coupling and simplified oneway coupling can be used for speeding up reactive transport simulations.Namely, in the framework of sequential coupling through operator splitting, in which at each time step first hydrodynamics is calculated on the simulation grid, and afterwards chemistry on each grid element, one could imagine substituting such an expensive calculation of chemistry by scaling one single batch model as proposed in the simplified one-way approach.This would ensure consistency between CO 2 available for transport and mineralised, eliminating the major uncertainty connected to the application of the simplified approach in post-processing mode.Also some other ideas presented in the present work could lead to improvements or further computational speedups.In particular the check for self-similarity of reactions could be applied in heterogeneous settings in order to identify a reduced number of reaction paths that need to be actually included in the reactive transport simulations, much in a sense of the reduction of complexity of the chemical system discussed by De Lucia and Kühn (2013) or, with focus on more numerically efficient kinetic laws modelled as first-order decay processes, Hellevang et al. (2013).

Conclusions
The simplified one-way coupling introduced by Klein et al. (2013) has been validated by means of comparison with fully coupled reactive transport simulations in a typical CO 2 underground storage setting in a saline aquifer, exploring one homogeneous and one heterogeneous case in terms of porosity and permeability.It was demonstrated that for such case studies the chemical reactions in each element of the simulation grid are with good approximation self-similar, which means on one hand that one single 0-D geochemical model can be used as proxy for the reactions occurring in all elements of the simulation grid, and on the other hand that the hydrodynamic transport of solutes plays a secondary role in comparison to the presence or not of the injected CO 2 , which is the true driving force of the chemical reactions.The choice of a threshold value for dissolved CO 2 considered geochemically active governs the convergence of the one-way coupling with the fully coupled simulations.This is particularly true for the given case study in which the migration of and exposure to dissolved CO 2 represents the single most significant discrepancy between the non-reactive and the fully coupled simulations.For the considered case studies an optimal mass fraction concentration of around 0.0005 was found to ensure the best matching of the outcome of the fully coupled simulations.Since in a real application this parameter is a priori unknown, the outcome of the simplified coupling has to be determined rather as bandwidth, estimating the optimal threshold for example by means of 1-D or 2-D simulations.However, given the advantage of performing coupled simulations on finely discretied grids with no simplifications and upscaling of heterogeneous features of the reservoirs, the uncertainty due to the simplified coupling appears justified.Furthermore, removing the computational burden for reactive transport simulations makes the simplified approach particularly adapted to sensitivity analyses, which are much needed given the uncertainty inherent to geochemical modelling.Edited by: L. Gross

Figure 1 .
Figure1.The reference geochemical model, scaled for initial porosity ϕ = 0.5 and initial water saturation S w = 0.5.Siderite (red) and dolomite (orange) are the main precipitating carbonate phases contributing to mineral trapping; calcite (black) is predicted to appear only transiently.Kaolinite (green), illite (grey) and albite (magenta) also precipitates.The newly formed phases are compensated by dissolution of anhydrite cement (not shown) and of chlorite (blue).In terms of porosity changes these reactions can be safely considered negligible in the considered time frame.

Figure 2 .
Figure 2. The spatial discretisation of the near-field Ketzin anticline used for the reactive transport model (northing direction is along the y axis).It comprises 2950 elements in one single layer, but with a 3-D structure following the topography of the target formation (magnified 5 times in the z direction).Depicted here is the heterogeneous distribution of porosity, derived from the upscaled Ketzin geological model.

Figure 3 .
Figure 3. Spreading of dissolved CO 2 in terms of aqueous mass fraction (XCO2aq) after 2100 years.Compared are the homogeneous (left column) and heterogeneous (right column) cases; first row shows the non-reactive and second row the reactive simulations (full coupling).At this simulation time there is a significant difference between the reservoir volume affected by dissolved CO 2 depending on whether or not chemistry is included in the simulation.Since the simplified coupling makes use of non-reactive simulations to estimate the exposure time of each element, this discrepancy will impact its outcome.However, this effect is mitigated by the introduction of a threshold filtering out small concentrations of aqueous CO 2 as geochemically inactive (see text for further details).

Figure 4 .
Figure 4. Homogeneous simulations: reservoir volume exposed to CO 2 distinguishing between elements exposed to gaseous phase (black lines) and to dissolved CO 2 only (red; blue following the threshold for the minimum considered concentration in terms of mass fraction).

Figure 5 .
Figure 5. Heterogeneous simulations: reservoir volume exposed to gaseous and dissolved CO 2 .

Figure 8 .
Figure 8. Histogram of similarities of total carbonate mineralisation in the homogeneous and heterogeneous cases.

Figure 9 .
Figure 9. Spatial distribution of the self-similarity (cf.text for definition) of chemical reactions relative to total carbonate mineralisation for the fully coupled reactive transport simulations in the homogeneous case.

Figure 10 .
Figure 10.Spatial distribution of the self-similarity for the heterogeneous case.

Figure 11 .
Figure 11.Total carbonate mineralisation in reservoir: comparison between full and simplified coupling.Influence of the choice of threshold value for dissolved CO 2 considered active for the coupling in the homogeneous and heterogeneous cases.

Figure 12 .
Figure 12.Influence of the threshold in case of slower kinetics, homogeneous case.
ongoing R&D activities are funded within the project COMPLETE by the Federal Ministry of Education and Research within the GEOTECHNOLOGIEN program.Further funding is received by VGS, RWE, Vattenfall, Statoil, OMV and the Norwegian CLIMIT program.The service charges for this open-access publication have been covered by a Research Centre of the Helmholtz Association.

Table 2 .
Minerals included in the model.