Chemical Compositions in Salinity Waterflooding of Carbonate Reservoirs: Theory

Higher oil recovery after waterflood in carbonate reservoirs is attributed to increasing water wettability of the rock that in turn relies on complicated surface chemistry. In addition, calcite mineral reacts with aqueous solutions and can alter substantially the composition of injected water by mineral dissolution. Carefully designed chemical and/or brine flood compositions in the laboratory may not remain intact while the injected solutions pass through the reactive reservoir rock. This is especially true for a low-salinity waterflood process, where some finely tuned brine compositions can improve flood performances, whereas others cannot. We present a 1D reactive transport numerical model that captures the changes in injected compositions during water flow through porous carbonate rock. We include highly coupled bulk aqueous and surface carbonate-reaction chemistry, detailed reaction and mass transfer kinetics, 2:1 calcium ion exchange, and axial dispersion. At typical calcite reaction rates, local equilibrium is established immediately upon injection. In SI, we validate the reactive transport model against analytic solutions for rock dissolution, ion exchange, and longitudinal dispersion, each considered separately. Accordingly, using an open-source algorithm (Charlton and Parkhurst in Comput Geosci 37(10):1653–1663, 2011. https://doi.org/10.1016/j.cageo.2011.02.005), we outline a design tool to specify chemical/brine flooding formulations that correct for composition alteration by the carbonate rock. Subsequent works compare proposed theory against experiments on core plugs of Indiana limestone and give examples of how injected salinity compositions deviate from those designed in the laboratory for water-wettability improvement.


3 a V n
Mineral surface area per solids for species adsorption ( m −1 ) a V r Reactive mineral surface area per solid volume ( m −1 ) C iE ith element concentration ( mol L −1 or M) C j jth species concentration ( mol L −1 or M) CEC Cation exchange capacity ( mol L −1 ) D Column diameter ( m) Da rxn = a V k k rxn L∕u Dissolution reaction Damköhler number Da mt = a V k k mt L∕u Mass-transfer Damköhler number D g Characteristic grain size ( m), See Appendix A D iE Dispersion coefficient of element i ( m 2 s −1 ) D Lj Longitudinal dispersion coefficient of species j ( m 2 s −1 ) D m Molecular diffusion coefficient ( m 2 s −1 ) Fluid accessible porosity k Volume fraction of mineral k K Bulk equilibrium constant (see Table 1) K S Surface reaction equilibrium constant (see Table 2) K CaNa Equilibrium exchange constant ( m −1 ) K SP Calcite solubility product ( mol 2 L −2 ) k m Convection-diffusion mass transfer coefficient ( ms −1 ) k rxn Intrinsic dissolution reaction-rate constant ( ms −1 ) L Column length ( m) n CEC Cation exchange capacity ( mol m −2 ) n iE Adsorbed molar concentration of element i ( mol m −2 ) n jk Concentration of species j adsorbed on mineral k ( mol m −2 ) Pe L = uL∕D L Longitudinal Péclet number p ijk Stoichiometric coefficient of element i in species j adsorbed on mineral k Q Volumetric flow rate ( m 3 s −1 ) Re Reynolds number, see SI r j Dissolution rate of mineral j ( mol m −2 s −1 ) r iE Net production rate of element i from minerals ( mol m −2 s −1 ) r ik Net production rate of element i from mineral k ( mol m −2 s −1 ) Sc Schmidt number, see SI Sh Sherwood number, see SI s ij Stoichiometric coefficient of element i in species j -plane surface charge density ( μCcm −2 ) t Time ( s) u Superficial velocity ( ms −1 ) x Axial distance (m) x Ca 2+ = 2C Ca 2+ ∕C Cl − Normalized bulk concentration of calcium ions x Na + = C Na + ∕C Cl − Normalized bulk concentration of sodium ions y Ca 2+ = 2n Ca 2+ ∕n CEC Normalized surface concentration of calcium ions y Na + = n Na + ∕n CEC Normalized surface concentration of sodium ions Abbreviations LSW Low-salinity waterflooding

Introduction
Carbonate reservoirs contain over 50 % of current proven oil reserves Ehrenberg and Nadeau (2005), Schlumberger Ltd (2020), BP plc (2020). Improved oil recovery is challenging in carbonates, primarily because of their complex heterogeneities. Importantly, carbonate rocks also exhibit significant aqueous solubility. This means that during waterflood, contact with the reservoir rock can significantly alter the injected brine chemistry. Any water-based flooding process that requires injection of chemicals, even simple brines, is plagued by alteration of the injection composition due to mineral dissolution and precipitation, and by adsorption and ion exchange with the porous-rock minerals. Oil-displacement processes that require specifically tailored laboratory chemical formulations can be changed enough by rock interaction to be ineffective under field conditions. A pertinent example is low-salinity waterflooding (LSW), in which injected brine is typically more dilute compared to the in situ brine. LSW is appealing because of its simplicity and low cost. Although extensively studied, the oil-recovery mechanism(s) of LSW remains obscure (Sohal et al. 2016;Purswani et al. 2017;Katende and Sagala 2019;Bartels et al. 2019;Hao et al. 2019). There is some consensus that injected brine compositions should be chosen to alter rock wettability toward more water wet by controlling, for example, surface charge at the rock/brine interface and/or by detaching asphaltene molecules from the rock surfaces (Norrman et al. 2018;Hu et al. 2018;Collini et al. 2020;Shaik et al. 2020). Mineral surface chemistry is a sensitive function of in situ brine composition, in particular aqueous pH. In reactive porous media, chemical-flooding formulations devised in the laboratory do not transfer directly to field applications.
The goal of this work is to establish how carbonate rock changes brine chemistry during flow through a reactive porous medium. We provide a new tool to establish injection compositions that achieve the desired compositions in the pores of the rock. When flowing in a carbonate rock, aqueous compositions are altered by rock dissolution/precipitation kinetics and mass transfer of multiple minerals, ion exchange/adsorption on those minerals, initial reservoir CO 2 level, and solute dispersion, among others. As the solution brine composition adjusts during flow in a reactive rock, so do rock surface compositions in equilibrium with the pore solutions. Local crude oil/brine/rock wettability depends on both surface and bulk compositions, although the exact mechanism(s) for wettability alteration is not currently known.
Our computation methodology is that of transport in reactive porous media which has a rich history, e.g. (Steefel et al. 2005;Tournassat and Steefel 2019). In the specific case of aqueous flow in carbonate media, considerable effort has also been expended (Nasralla et al. 2015;Chandrasekhar et al. 2016;Korrani et al. 2016;Farajzadeh et al. 2017;Eftekhari et al. 2017;Maes and Geiger 2018;Bonto et al. 2019). With the exception of Maes and Geiger (2018) who use a lattice Boltzmann approach, all works approximate numerically the continuum convective mass-transport equations. Many researchers (Nasralla et al. 2015;Chandrasekhar et al. 2016;Korrani et al. 2016;Farajzadeh et al. 2017;Bonto et al. 2019) couple their in-house numerical transport codes with the well-documented and verified reactive-species equilibria calculated by a USGS code called PHREEQC (Parkhurst and Appelo 1999). Unfortunately, only Korrani et al. (2016) test their calculations against experiment. Few researchers evaluate surface species that are essential to understanding of surface wettability. To our knowledge, no efforts include all effects of dissolution kinetics and mass transfer, ion exchange, and dispersion for multiple aqueous-reaction species in the bulk and at the surface. There are similar efforts to account for such effects in a geochemical context (Fantle and DePaolo 2006;Zhang et al. 2020). However, these works do not consider multiple reacting species with ion exchange, dispersion, and realistic dissolution and mass-transfer kinetics. None utilizes available open-source codes.
In our description of reactive flow in carbonate rock, we include hydrodynamic dispersion, kinetic and mass-transfer effects in the dissolution of calcite together with exchange reactions on clay, and calcite surface equilibria. Our 1D single-phase reactive-transport code is built on PHREEQC. To ascertain surface-charge density profiles, we include surface speciation equilibria.
Here, we focus on the general theoretical formulation and on the validation of the numerics by benchmarking against separate analytic solutions of reactive transport, ionexchange dynamics, and tracer dispersion. In a companion paper, model comparison is made to experimental ion-concentration histories from flow through Indiana limestone core plugs (Indiana Limestone Institute 1979). Application of our findings to examine in situ compositions relevant for LSW will be addressed in a companion paper. Below we describe our theoretical approach followed by example calculations. Description and validation of the numerical approach are outlined in Supporting Information.

Theory
Prediction of aqueous compositions during frontal displacement through porous carbonate rock is a challenge because of the rock reactivity, and because of the involved speciation chemistry in the bulk aqueous and rock-surface phases. For example, Table 1 lists typical carbonate species and their equilibrium constants. In general, the convection-diffusion equation must be written for each element including aqueous speciation reactions, rock reaction kinetics, rock surface adsorption/ion exchange and surface ion complexation (Lake et al. 2003;Felder and Rousseau 2005;Al-Shalabi et al. 2014;Yutkin et al. 2018), and axial dispersion. Here, we consider a 1D reactive transport model.

Reactive Transport
Tables 1 and 2 describe the species that participate in the aqueous phase-carbonate reactions, along with their equilibrium constants at ambient temperature for bulk and surface chemistry, respectively . Units in the equilibrium constants are listed to aid the reader. Table 2 is slightly changed from that of Yutkin et al. (2018)) to include sodium surface species in reaction S 6 . The symbol > in Table 2 denotes species at the calcite surface in surface concentration units. Each of the bulk aqueous species, in addition to chloride ions, must be tracked during flow through a porous carbonate rock. This is a daunting task, especially when several reactive reservoir minerals, such as anhydrite and/ or dolomite, are present in addition to calcite. We assume insignificant reactive loss of the minerals under typical conditions of low-salinity waterflood. All aqueous species are in reaction equilibrium at ambient temperature.
To account for aqueous equilibrium reactions, we conserve element masses (Felder and Rousseau 2005;Al-Shalabi et al. 2014). During frontal displacement in a 1D reactive porous medium, conservation of element i demands that (Yeh and Tripathi 1991;Steefel and Lasaga 1994;Walter et al. 1994;Qiao et al. 2015;Eftekhari et al. 2017) where is the fluid-accessible porosity; C iE (t, x) is the ith element pore molar concentration; a V n is mineral surface area per solids volume for species adsorption; a V r is reactive mineral surface area per solids volume; n iE (t, x) is the adsorbed molar concentration of element i on solid surfaces per unit solids area ( mol/m 2 ); t is time, x is axial distance; u is superficial velocity (i.e., u∕ is the frontal advance rate); D iE is the dispersion coefficient of element i; and r iE is the net production rate of element i from minerals per unit solids surface area ( mol/m 2 ∕s ). All aqueous species are in reaction equilibrium as enumerated in Table 1, except for Reaction 6 for the carbonate mineral which is accounted for separately by dissolution reaction kinetics. Element concentrations are those of all aqueous species j in Table 1: (1)  a The hyphen appearing in the surface species distinguishes between ions at the calcite surface (species to the left of the hyphen) and those at the -plane (species to the right of the hyphen) # Reaction a log 10 K S Units where C j is the concentration of species j in the pores and s ij is the stoichiometric coefficient of element i in species j. Likewise, the adsorbed element amounts are related to adsorbed species j in Table 2 by where k is the volume fraction of mineral k; a Vk is the surface area per unit solid volume of mineral k; n jk (t, x) is the molar concentration of species j adsorbed on mineral k per unit solid area ( mol/m 2 ); and p ijk is the stoichiometric coefficient of element i in species j adsorbed on mineral k. Table 2 lists typical surface species for carbonate rock ). The element dispersion term is given by the expression where D Lj is the longitudinal dispersion coefficient of species j. Lastly, the element mineral dissolution rate in Eq. 1 is expressed for minerals k as where r ik is the net production rate of element i from mineral k per unit surface area of the solid ( mol/m 2 ∕s). The boundary conditions for Eq. 1 by Danckwerts (1953), and Pearson (1959) prescribe continuity of the flux of each species at the column inlet and outlet. When written in terms of elements i, these boundary conditions become and with where x = 0 is the medium inlet and x = L is the medium outlet.
Equations 1-8 are general. For our particular application to limestone rock, we consider only reactive calcite and combined ion exchange on calcite and clay minerals. We include clay as an adsorptive mineral because of the ion-exchange characteristics we find in our Indiana Limestone corefloods described in a companion paper. For the aqueous phase, Table 1 lists 12 species not including excess water. To this number, we add any brine species that do not participate in the carbonate-reaction or ion-exchange x equilibria. Because we consider brine as an aqueous NaCl/CaCl 2 mixture, chloride ions are added as a species. Table 2 lists another eight surface reaction species denoted by the prefix >. Following all bulk and surface species is a challenge. However, when those species are in local equilibrium, only elements need be tracked (Al-Shalabi et al. 2014;Yutkin et al. 2018). In our example calculations and experimental corefloods, there is no free gas phase, and Reaction 1 does not participate, although it is pertinent to reservoir simulations. As a result, the concentration of aqueous carbon dioxide in the water-filled pores of the medium is set by calcite equilibrium. In the remaining Reactions 3-10 in Table 1, there are four independent elements: sodium, calcium, hydrogen (sans excess water), and carbonate, i.e., carbon and oxygen appear in all species in a 3/1 ratio. There is no need to impose electroneutrality in Reactions 3-10 because each is charge-neutral (Al-Shalabi et al. 2014;Yutkin et al. 2018). Sodium ions do not originate from calcite dissolution. Consequently, Reactions 9 and 10 only occur when sodium ions are present in the aqueous brine. Chloride ion becomes the fifth independent element. Overall electroneutrality is imposed on the brine. Once the transient profiles of the elements and corresponding aqueous species are established, surface-species concentration profiles follow from the known surface-reaction equilibria in Table 2. Implicit in the surface calculations is the assumption that reactive minerals remain in excess and are not depleted to alter rock porosity.
Not all element balances contain all terms in Eq. 1. For example, chloride ions do not participate in ion exchange and in carbonate-species reactions. Consequently, for chloride ions, ion-exchange and reaction-rate terms are absent, and Eq. 1 reads The dispersion coefficient in Eq. 9 does not have a species-specific subscript because all aqueous species in Table 1 are dilute and have about the same aqueous molecular diffusion coefficients, D m ∼ 10 −9 m 2 ∕s . Equation 9 indicates that chloride ions behave as an ideal tracer. Without dispersion, chloride ions transport along the core as a salinity wave with no holdup (Pope et al. 1978).
Conversely, calcium retains all terms in Eq. 1. Thus, we have that where i designates calcium. To understand the role of dispersion, we vary the longitudinal dispersion Péclet number, Pe L = uL∕D L . Table 1 indicates that calcium species j = 1 − 3 correspond to Ca 2+ , CaOH + , CaHCO + 3 aqueous ions, respectively. We label mineral 1 as calcite. Mineral 2 represents adsorptive minerals including calcite and clay. On calcite and clays, sodium/calcium ion exchange obeys local equilibrium or where the ion-exchange isotherm n Ca (C Ca 2+ ) is specified by classical mass action demonstrated in Appendix B In Eq. 12, y Ca 2+ = 2n Ca 2+ ∕n CEC , x Ca 2+ = 2C Ca 2+ ∕C Cl − , K CaNa , is the ion-exchange equilibrium constant, C Cl − = 2C Ca 2+ + C Na + and n CEC = 2n Ca 2+ + n Na + is the cation exchange capacity expressed in units of mole of total exchange sites per unit surface area. In adopting Eq. 12, we assume that ion-exchange rates are fast and that local equilibrium in attained. However, typical clay minerals in reservoir rock occur in book stacks, and the assumption of local equilibrium may be questioned (Mansa et al. 2017;Khoury 2019;Zhang et al. 2001;Wilson et al. 2014;Khalifa and Morad 2015). It is traditional to express exchange capacity in equivalents per liter of pore volume, or CEC = 10 −3 n CEC a V 2 2 ∕ (Pope et al. 1978;Hill and Lake 1978). From Eq. 12, ion-exchange amounts depend directly on ionexchange capacity and inversely on chloride-ion concentration.
The remaining term to specify in Eq. 10 is the dissolution rate of calcite. We utilize the expression derived by Yutkin et al. (2018, Appendix A): where k rxn is the intrinsic dissolution reaction-rate constant (m/s) and K SP is the solubility product of calcium carbonate (i.e., K 6 in Table 1). Based on the extensive review of Morse and Arvidson (2002) and the references therein, we find that Eq. 13 is experimentally verified. At low pH, k rxn is a strongly increasing function of hydrogen-ion concentration. At pH values above about 6, however, k rxn is constant at approximately 10 −5 m/s under ambient conditions (Morse and Arvidson 2002).
Equation 13, however, does not account for the convective mass-transfer resistance encountered by the dissolving calcium and carbonate ions (Bird et al. 2007). For pH values less than about 10, Appendix A in SI demonstrates that where k m is the convection-diffusion mass transfer coefficient appropriate to flow through porous media (Bird et al. 2007). Expressions for the mass transfer coefficient in packed beds are available in Rexwinkel et al. (1997) and Bird et al. (2007, page 486). We adopt Eq. 14 to quantify the dissolution kinetics of calcite.

Numerics
Numerical solution of Eqs. 10-14 along with Eq. 1, 6, and 7 written for carbonate, sodium, and hydrogen elements is challenging. All species in Table 1 must be accounted for through Eq. 2. We use the 1D explicit central finite-difference implementation of the publicly available PHREEQC code of the US Geological Survey (Charlton and Parkhurst 2011;Yutkin et al. 2018) along with the well-documented extensive aqueous chemistry speciation calculations embodied in PHREEQC. Implementation details are described in Appendix C of SI. Appendix C also benchmarks our numerics against the derived analytic solutions. Distances and times in the species profiles and histories below are dimensionless and reported � as x/L and ut∕ L = PV , respectively. Table 3 2 lists the pertinent parameters of the example calculations below unless otherwise noted.
In the calculations that follow, reaction and mass-transfer rates are varied through reaction and mass-transfer Damköhler numbers: Da rxn = 1 a V 1 k rxn L∕u and Da mt = 1 a V 1 k m L∕u , respectively. Ion-exchange strength is gauged through the cation exchange capacity ( CEC ), and dispersion is varied through the Péclet number: Pe L = uL∕D L .

Results and Discussion
We outline the theoretical results in the order of ion exchange alone, mineral reaction alone, and combined ion exchange with mineral reaction. The roles of reaction rates, mass transport rates, cation exchange capacity, and dispersion are discussed as appropriate. We consider ion exchange alone first.

Ion Exchange
Figures 1 and 2 illustrate the role of ion exchange for sodium/calcium exchange on calcite/clay minerals with no dispersion. Logarithmic concentration profiles at 0.5 PV injected are shown on the left, and histories are shown on the right. Calcite mineral is inert (i.e., Da rxn = 0 ). Cation exchange capacity is 0.345 mol/L pore fluid in Fig. 1 and  3.45 × 10 −3 mol/L in Fig. 2, respectively. These values cover a range characteristic of Berea sandstone (Hill and Lake 1978;Pope et al. 1978). The equilibrium exchange constant K CaNa is 6600 m −1 , estimated from Hill and Lake (1978), Pope et al. (1978).  Examination of the calcium ion-exchange profiles in Fig. 1a reveals a salinity shock dropping the calcium concentration from the initial porous-medium value of 0.1 M down to 4 × 10 −3 M . The calcium concentration behind the salinity shock is set according to equal rock adsorption of calcium at the initial and injected chloride concentrations outlined in Appendix C of SI. A zone of constant calcium concentration (at 4 × 10 −3 M ) follows the salinity shock. This zone precedes a spreading-wave zone of calcium desorption from the rock surface (Hill and Lake 1978). The fastest wave of this zone corresponds to the calcium concentration of 4 × 10 −3 M and penetrates out to about x∕L = 0.35 where it intersects the salinity-induced zone. The slowest calcium concentration wave corresponds to zero injected calcium and penetrates out to about x∕L = 0.01 . There are four calcium-concentration zones in the column: the initial calcium concentration, the salinity-induced concentration zone behind the indifferent salinity wave, a spreading calcium desorption zone, and a final injected-concentration zone. Because hydronium ions are not involved in classical sodium/calcium mass-action ion exchange, there is no effect of pH on ion exchange. With no mineral reaction, injected hydronium ions simply flush those initially present. Chloride behaves as an inert tracer; sodium follows from charge balance. At a CEC of 0.345 mol/L, washout of the ion-exchangeable calcium in Fig. 1a takes many pore volumes. Even without rock or aqueous-phase reactions, ion exchange significantly alters injection chemistry. Detailed information on ion-exchange wave dynamics is available in Appendix C of SI. Figure 2 illustrates the importance of the cation exchange capacity on concentration profiles (a) and histories (b) under the same conditions as those in Fig. 1. However, the cation exchange capacity is lowered by a factor of 100 to 3.45 × 10 −3 mol/L. With the considerably lower ion-exchange capacity, much less calcium resides on the rock surface. As observed in Fig. 2a, the abrupt calcium-concentration decline remains at the chloride tracer front but the salinity zone collapses. Spreading calcium desorption waves emanate directly following the tracer chloride wave. Washout of the exchangeable calcium in Fig. 2b is much more extensive than that in Fig. 1. Figure 3 investigates the role of axial dispersion on classical sodium/calcium ion exchange, where the amount of dispersion is characterized by a longitudinal Péclet number of Pe L = 100 . The remaining conditions are identical to those in Fig. 1. As anticipated, all shock fronts smooth with the calcium zone spreading more at the rear than at the front. The injected concentrations do not reach the column outlet for many pore volumes. Hydrodynamic dispersion completely smears the constant salinity zone. The seemingly large(er) dispersion is due to the ordinate logarithmic scale. The 1D advective algorithm in PHREEQC well captures ion exchange in porous media (see Appendix C for validation against analytic solutions).

Mineral Reaction
Figures 4, 5, and 6 report the effects of calcite dissolution on concentration profiles and histories with no ion exchange or dispersion. Injected and initial column conditions are the same as those in Figs. 1, 2, and 3. Figure 4 illustrates the role of the reaction Damköhler numbers between 10 −3 and 10 2 for a fixed mass-transfer Damköhler number of 100, large enough for mass transfer to have no effect. The reaction front proceeds through the medium at the frontal advance rate. For all Da rxn , there is a reaction zone at the beginning of the reactive medium where calcite dissolves. Beyond that reaction zone, calcium concentration is constant: dissolution products are convected away as fast as they dissolve. At Da rxn ≥ 10 , however, the injected solution beyond the narrow inlet zone reaches equilibrium with the carbonate rock. No further dissolution occurs. Because reaction Damköhler numbers are very large in reservoir applications, this result confirms our previous order-of-magnitude estimates that rock dissolution occurs at the injection wellbore ). Consequently, once past the reaction zone, any rock-dissolution-based mechanisms in low salinity waterfloods are impossible (Chen et al. 2017). From the corresponding concentration histories in Fig. 4, rock reaction neither speeds nor slows species transport advance. Reaction kinetics does not cause tailing of the effluent concentrations as sometimes reported (Austad et al. 2015). Rather, effluent concentrations are constant after breakthrough and approach rock equilibrium for high reaction Damköhler numbers.  Fig. 4, except that the mass-transfer Damköhler number varies between 10 −3 and 10 2 for a fixed Da rxn = 100 , characteristic of rock equilibrium beyond the narrow inlet reaction zone. Concentration profiles and histories are similar to those in Fig. 4 except for a weaker dependence on Da mt . Mass transfer slows dissolution kinetics dictated by Eq. 14. Beyond the reaction zone, the dissolution product concentrations are constant. Figure 6 is identical to Figs. 4 and 5, except now sodium, hydronium (i.e., pH), and chloride concentration profiles and histories are also shown in addition to that for calcium.  . Injected pH is 7 (not shown on the scale of the graph). However, after the reaction zone, calcite dissolution raises the pH from 7 to almost 10. Thus, rock reaction controls aqueous chemistry occurring during a displacement process. In our calculations, CO 2 is confined to the aqueous phase. In reservoirs, however, CO 2 may be present as free gas and/or dissolved in oil phases. In which case, we find pH is controlled by CO 2 concentration in the reservoir. Therefore, composition design of advanced waterflooding processes should account for conditions expected in a reservoir, including CO 2 sources or sinks.  Fig. 7, except those for pH, are virtually indistinguishable from those in Fig. 1. At a CEC of 0.345 mol/L, ion exchange dominates reaction chemistry. Classical ion exchange is independent of pH. Thus, the pH profile and history in Fig. 7 arise because of calcite dissolution. Injected pH (= 7) rises at the column inlet to approximately 9.6, similar to Fig. 6, but then falls to 9.25 near the salinity front. This is due to the ion-exchange salinity zone influencing aqueous chemistry in addition to calcite dissolution. With reaction chemistry alone in Fig. 6, the pH remains at 9.9 up to the salinity shock. In the history of Fig. 7b, breakthrough of the salinity-shock wave occurs at one pore volume followed by the salinity zone, by spreading ion-exchange waves, and, finally, by rock-reaction equilibrium concentrations. Figure 8 is identical to Fig. 7, except that more species are portrayed. Particularly important are the profile and history of the -plane surface charge density, . Since the calcite surface contains equal numbers of negative and positive sites , the -plane charge governs the subsequent adsorption of charged aqueous species. For the conditions of Fig. 8, is positive at 0.95 μC/cm 2 initially in the core and falls to essentially neutral as the waterflood proceeds. The charge density in the -plane depends on aqueous sodium and CO 2 concentrations. Sodium ions in Reaction S 6 disfavor Reaction S 5 and expose negative charge obtained though Reaction S 3 . Therefore, although positively charged, adsorbed sodium ions promote negative charge density in the -plane. Additional aqueous CO 2 creates more aqueous bicarbonate ions that also makes the -plane more negative via reaction S 3 in Table 2. Figure 8 highlights the power of our reactive-transport model, because surface charge density is often thought to play an important role in reservoir wettability.

Ion Exchange and Mineral Reaction
If the cation exchange capacity is lowered to CEC = 3.45 × 10 −3 mol/L in Fig. 9, the resulting concentration profiles and histories are not exactly the same as those in Fig. 6 for reaction alone. Even at a low CEC , the fingerprints of ion exchange are visible. Figure 9 resembles more closely Fig. 2 for ion exchange alone. The spreading waves of calcium desorption surprisingly overpower rock dissolution. However, in Fig. 9b, the calcium concentration is set by calcite equilibrium (at ∼ 1.2 × 10 −4 M ), contrary to zero calcium concentration in Fig. 2b at large throughput (not shown in the figure). Finally, Fig. 10 repeats Fig. 7 but now with dispersion at a longitudinal Péclet number of Pe L = 100 . Remnants of Fig. 9 are seen, but are smoothed, as expected.

Conclusions
This work studies numerically the 1D simultaneous reactive transport, mass-transfer, ion exchange, and dispersion during flow through a calcium-carbonate porous medium. Dissolved aqueous carbonate species are in equilibrium as is ion-exchange kinetics. Carbonate dissolution kinetics is quantified by a rigorous rate expression. The key to the modeling effort is the continuum mass conservation of elements. All aqueous carbonate and brine species are followed as is pH and CO 2 . For the first time, equilibrium surface species are also followed including the surface charge density and zeta potentials during frontal displacement provided that pore-solution concentrations are not strongly influenced by surface equilibria. Because reaction kinetics are fast, as characterized by the reaction and mass transfer Damköhler numbers, injected pore solutions attain equilibrium with the carbonate rock except for a narrow zone near the medium entrance.
All results are obtained using the open-source software PHREEQC, including aqueous speciation and numerical solution to the 1D transport equations. Because of the highly nonlinear coupled behavior of species during mineral dissolution and surface ion exchange, we have performed careful verification of the PHREEQC transport code against analytic solutions over a wide range of the input parameter values, and, where applicable, we have reported the reliability ranges. Currently, our model implementation is restricted to a single mineral but is readily extended to multiple minerals. We assert that 1D transport in the reactive porous media with the simultaneous and coupled bulk and surface speciation, ion exchange, and dispersion is well captured by the open-source code PHREEQC. Our calculations provide a new tool to assess how injected aqueous brine solutions alter upon contact with reactive rock. Accordingly, effective in situ brine formulation can be designed to improve oil recovery from reactive-mineral reservoirs.