Modeling of Contact Angle versus pH and Derivation of Contact Angle versus Pressure in Deep Saline Aquifers under Geological Carbon Storage

Geological storage of anthropogenic carbon dioxide is regarded as a technically and economically viable strategy for mitigating carbon dioxide induced climate warming. 
 
Central to geological storage of anthropogenic carbon dioxide is the water rock interaction, which has a direct bearing on pH induced wettability evolution in saline aquifers. Consequently, understanding contact angle trend versus injected gas pressure is useful, considering its relationship to pH evolution in formation brine due to dissolved gas at prevailing temperatures and salinities. Several research works have published experimental data on contact angle versus pressure pertaining to geological conditions of anthropogenic carbon storage. In the present study, we have used thermodynamic theories relating to a surface charge model, contact angle and the classical Nernst equation to derive a logarithmic pH dependent contact angle equation. Considering the relationship between carbon dioxide solubility and pressure for a given temperature and salinity as well as the link between pH and the extent of solubility, we have plotted calculated contact angles versus corresponding pressures. Results of the plots obtained compare well with literature values. Therefore, given the lack of theoretical approach regarding contact angle versus pressure, our research work fills the knowledge gap considering the novelty in the derivation of the pH dependent contact angle equation.


Introduction
In recent times, researchers in the industry and academia have made concerted efforts to understand both the technical and geological aspects of carbon dioxide geological storage as well as its potential environmental impacts. In this regard, noteworthy data have been published on the effect of water rock interactions on contact angle/wettability, which affects the distribution of injected carbon dioxide and resident formation brine in the aquifer. (Jung & Wan, 2012;Jung & Wan, 2011, Kim et al., 2012. Much of contact angle data related to saline aquifer carbon storage have been obtained based on measurements on aquifer rock minerals under supercritical conditions that reflect deep saline aquifer geological conditions (Kim et al., 2012). The consensus of researchers so far regarding the water rock interaction and wettability evolution of saline aquifers under geological carbon storage is that aquifers will undergo dewetting (Saraji et al., 2013;Kim et al., 2012). Consequently, trends in contact angle versus pressure (Jung and Wan, 2011;Jung and Wan, 2013) which translates to pH decrease with pressure show increasing contact angle.
When carbon dioxide is injected into the formation above the threshold capillary pressure, a two-phase flow regime emerges where there is an interface between invaded gas and resident brine and an interface between resident brine and pore surface. Under this condition, it is possible to have a definite thickness of brine wetting film between the vapor phase and pore surface (Tokunaga, 2012) and an overlap of electric double layers and their repulsion. Following the dissolution of carbon dioxide and its subsequent hydration into carbonic acid and its dissociation, the adsorption of hydrogen ions on pore surface due to water rock interaction will reduce surface charge density, which will reduce double layer repulsion to destabilize the thin wetting film. These phenomena are what cause aquifer brine pH reduction and dewetting of aquifer rocks minerals by supercritical carbon dioxide (Saraji et al., 2013;Kim and Wan, 2012;Jung and Wan, 2011;Jung and Wan, 2012). theoretical point of view. We contend that pH induced surface charge regulation theories (Kosmulski, 2010)) can be exploited to embark on a theoretical understanding of dewetting trends versus pressure in saline aquifers under different conditions of geological carbon storage. Theoretically, for a saline aquifer under geological carbon storage, the solubility of carbon dioxide which is a precursor to carbonic acid formation and pH reduction in saline aquifers is controlled by salinity, pressure and temperature (Springer et al., 2012;Li et al., 2015. Therefore, brine pH will also depend on these physicochemical properties. Also, several experimental data on carbon dioxide-brine interfacial tension show that for a given temperature and salinity of aquifer, interfacial tension decreases with pressure and achieves a stable value at high pressures above the critical pressure (Chalbaud et al., 2010). The objective of this paper was to theoretically study dewetting trends (contact angle increase) versus pressure as revealed experimentally in the literature. Accordingly, we studied theoretically, pH induced dewetting of silica surface using our own derived equation and literature source data. We have shown a logarithmic relationship between pH and contact angle, which was vital for revealing observable contact angle versus pressure trend as found in the literature. For this reason, we chose silica as the predominant sandstone aquifer mineral for two reasons. First, sandstone formations are the most abundant with characteristic intergranular porosity (Ehrenberg & Nadeau, 2005), which make them the best candidates for geological carbon storage. Second, several sequestration projects worldwide have sandstone formations as host geologic repositories (Labus et al., 2015;Liu et al., 2012;Zhou et al., 2010).

Surface Charge Regulation Induced Wetting and Dewetting
Surfaces may often be naturally charged, or electric fields may be employed to manipulate fluids, causing electric effects to be crucial components that influence wetting or dewetting phenomena. (Nita, et al., 2017). For instance, the behavior of polyelectrolyte adsorption on a substrate is dependent upon charge density, pH, temperature, and ionic strength, which has been exploited to achieve desired wettability (Yoo et al., 1998). The mechanism of electrostatic induced wetting alteration is intimately linked to the isoelectric point (IEP-pzc) concept which is the electrokinetic equivalent of the point of zero charge pH (pzc) (Kosmulski, Isoelectric points and points of zero charge of metal (hydr)oxides: 50years after Parks' review, 2016). At the IEP, the particle surface is electrically neutral, and the magnitude of zeta potential is equal to zero (Moulin & Roques, 2003). In light of the stabilization of the thin wetting film based on the classical theory of Derjaguin-Landau-Verwey-Overbeek (DLVO) (Hall et al., 1983); electrostatic repulsion stabilization of the film vanishes, leaving a predominant van der Waal contribution (Hano et al., 2012), which causes dewetting when adsorption of potential determining ions occur. Consequently, several reasons regarding wetting transition related to low salinity water flooding have been proposed, central to the entire subject of the double layer expansion by surface charge density increase which can be theoretically and quantitatively linked to disjoining pressure forces (Ding & Rahman, 2017;(AlQuraishi & AlHussinan, 2015;Mehana & Fahes, 2018;Ashraf et al., 2010). In the context of dewetting of silica exposed to aqueous species of dissolved carbon dioxide (Kim et al., 2012;Jung and Wan, 2011;Jung and Wan, 2012) , decrease in surface charge density due to adsorption of protons from dissolved and dissociated species of carbon dioxide will be responsible for decrease in double layer repulsion, which corresponds to decrease disjoining pressure forces and eminent dewetting. In light of pH dependence of solid-liquid interfacial (Amadu and Adango, 2019), adsorption of hydrogen ions on aquifer rocks following carbon dioxide injection will increase solid-liquid interfacial tension to a maximum (Parks, 1984) at the point of zero charge pH, which in light of Young's phenomenological equation (Rusanov et al., 2004) corresponds to increase in contact angle.

Vapor-Liquid Equilibrium Theory for Carbon Dioxide Solubility in Saline Aquifers
When carbon dioxide is injected into the saline aquifer at a given temperature and salinity of formation brine, there will be solubility of brine in the vapor phase and solubility of vapor phase in brine phase such that phase mole fractions can be calculated as (Hassanzadeh et al., 2008): In which /  is the salting out effect of carbon dioxide in sodium chloride solution, 2 CO m is the solubility of carbon dioxide in brine [mol/dm 3 ] and  − 2 CO V is the partial molar volume of carbon dioxide at infinite dilution [mol/dm 3 ].

Linking pH to Formation Brine Salinity and Temperature and Pressure
The carbon dioxide system in brine is characterized by four measurable parameters, namely the total alkalinity, the total inorganic carbon consisting of the sum of dissolved carbon dioxide (carbonate, and the bicarbonate), the pH and either the fugacity of carbon dioxide or its partial in the vapor phase (ref). At a given temperature and salinity of formation brine, carbon dioxide dissolves and hydrate to form carbonic acid in accordance with the following equation (Frank, et al., 2002): The first dissociation constant of constant of carbonic acid from Eq. (11) is given as (ref In this equation, 1 2 CO K is the first dissociation constant of carbonic acid and square bracket denotes the activity of species, which is the product of species activity coefficient and molar concentration. In Eq. (5), the extent of dissociation of the mixture

HCO H
ions is related to the ratio of the sum of the molar conductance of the mixture and the sum of the molar conductance of the ions. Assuming  moles per liter is the amount of dissolved and hydrated carbon dioxide of carbonic acid that dissociates to establish the equilibrium expressed by Eq. (5) can be written based on stoichiometry as: Where 3 2 CO H con is the activity of carbonic acid initially present. Assuming 3 2 CO H con is equal to the solubility of carbon dioxide at prevailing pressure temperature and pressure, Eq. (6) can be written based on the solubility as: (7) is a quadratic equation in alfa and the solution using the quadratic theory gives  as: In line with the carbon dioxide dissociation constant being temperature and salinity dependent (Aissa et al., 2015) its value will reflect these physicochemical properties.
Using the chemical definition of pH, the pH of water with a given salinity, temperature and pressure in a saline aquifer under geological carbon storage will be given from Eq. (8) as: In which  is the activity coefficient of hydrogen ion in solution.
Activity coefficient is related to ionic strength as (Schneider et al., 2004): Henceforth Eq. (9) will be written to reflect salinity, temperature and pressure dependence of pH as: The temperature and salinity dependence of the first ionization constant of carbonic acid will be invoked. This is given in the literature as (Aissa et al., 2015) In this equation, T is the absolute temperature [K]and I is the ionic strength [moll -1 ] which corresponds to salinity S

Derivation of pH Dependent Contact Angle Equation
In line with the effect of surface charge on the stability of the thin wetting film (Ciunel et al., 2005) and the theoretical aspect of surface charge density effect on contact angle (Puah et al., 2010), which has a direct bearing on wettability, the derivation of a pH dependent contact angle equation is necessary to achieve our principal objective. In this regard, we will assume a clean silica surface where the surface is characterized by a definite number density of surface hydroxyl functional groups or silanol (Zhuravlev & Potapov, 2006 ). We will further assume that the surface of the silica is in contact with a predominantly sodium chloride brine (Hanor, 1994). Given that the pH of saline aquifers under normal conditions of geological storage is near neutral (Cooke et al., 2000) and the point of zero charge pH of silica is on the average 3 (Amadu & Miadonye, 2017), the initial pH in the captive bubble method used by Farokhpoo et al., (2013) (see Figure 1), would be above the point of zero charge pH. The surface of silica typically exhibits amphoteric behavior and becomes charged from protonation/ deprotonation reactions (Illés & Tombácz, 2006): The following equilibrium reactions will describe the ionization of surface hydroxyl groups above the point of zero charge pH (Revil, 1998;Azam et al. 2013): Figure 1. Schematic view of a sessile-drop contact angle system Farokhpoor et al., 2013) The intrinsic ionization constants for Eq. (13) is a K . The derivation of the surface charge density assumes formation brine in contact with pore surface consisting predominantly of a one-meter square of surface ionizeable groups (silanol).
At a given temperature and salinity, dissolution of carbon dioxide will provide hydrogen ions in solution. Above the point of zero charge pH of silica, ionizable surface groups will produce negatively charged sites which will favor hydrogen ion adsorption. Based on the Gibbs excess equation (Chattorage, 2001) the following can be written to link change in solid-liquid interfacial tension to hydrogen ion adsorption: In this equation, is contact angle-degrees, pH is the negative logarithm to base 10 of the hydrogen ion activity of aqueous solution in contact with solid surface, F = Faradays constant [Cmol -1 ], is surface charge density-[Cm -2 ], R is universal gas constant [Jmol -1 K -1 ],T is absolute temperature,  is interfacial tension between fluid phases [Nm -1 ].
If the charge density is predominantly due to a single class of dissociable functional group with a distinct ionization constant the at any pH of the aqueous solution the surface charge density relationship to bulk solution pH is given as (Godt, 1981)  Using the theoretical definition of pH, the following can be written (Behrens & Grier, 2001)   The relationship between surface potential, pH of bulk aqueous solution and the point of zero charge pH is given by the Nernst equation as (Fairbank et al., 1997):

Carbon Dioxide Solubility Data
Based on the thermodynamic concepts of vapor liquid equilibrium presented in Section 2.2 of the present study, Duan and Sun, (2003) (Duan & Sun, 2003) have measured excellent data on carbon dioxide solubility in brine versus pressure and temperature under varying degrees of salinity encountered under geological carbon storage, covering different temperatures and pressures that reflect anticipated supercritical conditions of geological storage (Dávila et al., 2016). Solubility data have also been published (Spycher at al.,2005).
Based on Hanor's (Hanor, 1994) classification of sedimentary basin brine, we assumed a sodium chloride dominated saline aquifer brine. Therefore, we extracted solubility data from Duan and Sun (2003) and Spycher et al., (2005) covering anticipated salinities, temperature and pressure for geological storage. We used solubility data at 333.15 K from both references.
We also assume the predominant surface ionizable groups on the surface of sandstone are silanols.
We chose 3 different concentrations of aqueous sodium chloride solution. They are 1 M, 0 M, 2 M and 4 M.
Concentrations of 1 to 4 M can be found in formation brines (Jafari & Jung, 2020 To calculate contact angle versus pH, Eq. (25) was used with all the required physical constants. Excel was used for the calculation, whereby the pH calculated versus pressure using Eq. (25) and the solubility versus pressure for a set of salinities and temperatures extracted from the cited literature were used as variables.
To calculate pH using Eq (11), the activity coefficient of hydrogen ion at every pressure at a given temperature and salinity was calculated using Eq (10), which was then multiplied by the corresponding hydrogen ion concentration.
Calculated pH values were then substituted into Eq. (25) with corresponding calculated parameters.
Experimental works so far indicate that the interfacial tension between carbon dioxide and brine decreases with pressure for a give temperature and attains a stable value at high pressures above 10 to 15 MPa (Chalbaud et al., 2010), while increasing with salinity at a given temperature carbon dioxide-brine interfacial tension was determined as a function of temperature and salinity Using Appendix A: Interfacial tensions for 0, M, 1 M, 2 M and 4 M were 0.025 Nm -1 , 0.029 Nm -1 , 0.032 Nm -1 and 0.036 Nm -1 respectively. Appendix A was used for deducing interfacial tension values.
Substitution of deduced carbon dioxide-brine interfacial tension at a corresponding temperature and salinity into Eq containing corresponding calculated parameters facilitated calculation of contact angle versus pH for different salinities of brine. Since pH is calculated as a function of carbon dioxide pressure for different salinities at a given temperature, calculated contact angles can be plotted against pressure for a given temperature and salinity.
The value of dissociation constant Ka used in this study is 5.7 molar (Dove & Craven, 2005). The value of the universal gas constant used is 3.8 Jmol -1 K -1 ( (Katmar-Software, 2020). The activity coefficient of hydrogen ion was calculated using We used an effective hydration radius of 0.03nm for hydrogen ion (Lee, 2020). The number density of surface silanols is 4.5 m -2 (Zhuravlev & Potapov, 2006). The value for the point of zero charge pH ( pzc pH ) of silica surface used in this study is 3 (Amadu and Adango, 2017).

Results and Discussion
Generally, the temperature of carbon dioxide in deep saline aquifers will be equal to or above critical pressure. At such depths, the pressure of carbon dioxide due partly to pore pressure and partly to overburden pressure will be far above critical. Consequently, the interfacial tension is likely to achieve a stable value in accordance with observed trends from experimental data (Chalbaud et al., 2010). Therefore, integration of the theoretical model ((Eq.24)), assuming constant values of interfacial tension at given salinity and temperature is justifiable. Solubility data provided by Dun and Sun and Spycher et al (2005) reflect temperature, pressure and salinity conditions of geological carbon storage. Reaction of supercritical CO 2 with brine causes a significant pH drop to 3 units in the saline pore fluid (Qafoku et al., 2015), due to carbonic acid (as dissolved CO2) in the brine (Rathnaweera et al., 2016), which is confirmed by and Table 2, showing pH, pressure and contact angles for different salinities at 333.15. Kharaka et al., 2006 have also reported pH drop of formation brine from 7 to 2 pH units. (Kharaka, et al., 2006) In all tables, calculated values of pH and contact angle are reported to 2 decimal places. In Table 1, while results for 0 molar solution show consistent increase in contact angle with pH decrease, there is a slight deviation for 4 molar solution at pressures above 100 MPa and this is probably due to the rounding off to 2 decimal places.   Figure 3 also shows a similar plot based on Dun and Sun solubility data. The figure shows that below the critical pressure of carbon dioxide, the pressure gradient of contact angle is steep for all plots at a given temperature.

Conclusion
Solubility of injected carbon dioxide in saline aquifers is the principal cause of pH evolution from neutral to values below 3 units (Kharaka et al., (2006). Strongly linked to pH evolution is the phenomenon of pH induced surface charging of amphoteric surfaces (Kosmulski, 2010), which results in sold-liquid interfacial tension increases for the case of carbon dioxide injection into saline aquifer. In this regard, the Gibb's excess adsorption equation links change in solid-liquid interfacial tension to surface coverage which together with the classical Nernst equation and the phenomenological Young's equation provide the requisite thermodynamic impetus for modeling pH versus contact angle equation. In this study, we have successfully modelled pH versus contact angle equation by exploiting a thermodynamic theory relating to surface charge density. Since carbon dioxide solubility under geological conditions of sequestration is given as a function of pressure for a given salinity and temperature, it is possible to link pH to contact angle and pressure. The following sum up the conclusion of this theoretical study: 1. Contact angle versus pH for a saline aquifer under geological carbon storage due to pH evolution in response to dissolved injected carbon dioxide can be described by a logarithmic function, 2. The relationship between pH and pressure provides the opportunity to plot contact angle versus pressure, 3. Plots of contact angle versus pressure has a step gradient for pressures below the critical pressure of carbon dioxide and a gradient approaching zero for pressures above 10 MPa, which agrees with observed experimental trends in the literature (Jung & Wan, 2012), 4. Based on the model derived in this study, the higher the salinity the bigger is the contact angle at a given pressure for a given temperature, which is also consistent with literature source data (Jung & Wan, 2012).