A Comprehensive Equilibrium Model for the Phosphonate Scale Inhibitor-Carbonate System including Coupled Adsorption/Precipitation (Γ/Π)

Inorganic scale deposition is one of the main flow assurance issues in hydrocarbon production, and injection systems, leading to significant loss in production and subsequent expense to mitigate. Scale Inhibitors (SI) are chemicals that are commonly used to prevent the inorganic scale from building up in the system. SIs are usually injected into the porous media, where they can react and be retained, and released back into the produced brine as the production commences. This in situ scale prevention method is known as “ squeeze treatment ” , and the lifetime of this process is determined when the SI concentration in the return fluid falls below the MIC (Minimum Inhibitor Concentration), where the SI is sufficiently effective at preventing or efficiently reducing scale formation. In this study, a general chemical equilibrium model has been developed to simulate the chemical reaction of phosphonate scale inhibitors and their retention in a carbonate system (calcite), in the presence of an aqueous phase containing free divalent cations. The model couples together the following processes (i) the carbonate system, (ii) speciation of the SI, modelled as a weak polyacid, H n A, (iii) the metal (Ca 2 + , Mg 2 + ) binding – SI chelant interactions, (iv) phosphonate SI acid impurity reactions, (v) adsorption of the free and complex species to the rock surface and (vi) and finally, the precipitation of complex species (SI-Ca-Mg). To find the equilibrium conditions, charge balance and mass balances for key components in the system were considered,

Inorganic scale deposition is one of the main flow assurance issues in hydrocarbon production, and injection systems, leading to significant loss in production and subsequent expense to mitigate.Scale Inhibitors (SI) are chemicals that are commonly used to prevent the inorganic scale from building up in the system.SIs are usually injected into the porous media, where they can react and be retained, and released back into the produced brine as the production commences.This in situ scale prevention method is known as "squeeze treatment", and the lifetime of this process is determined when the SI concentration in the return fluid falls below the MIC (Minimum Inhibitor Concentration), where the SI is sufficiently effective at preventing or efficiently reducing scale formation.In this study, a general chemical equilibrium model has been developed to simulate the chemical reaction of phosphonate scale inhibitors and their retention in a carbonate system (calcite), in the presence of an aqueous phase containing free divalent cations.The model couples together the following processes (i) the carbonate system, (ii) speciation of the SI, modelled as a weak polyacid, H n A, (iii) the metal (Ca 2+ , Mg 2+ ) binding -SI chelant interactions, (iv) phosphonate SI acid impurity reactions, (v) adsorption of the free and complex species to the rock surface and (vi) and finally, the precipitation of complex species (SI-Ca-Mg).To find the equilibrium conditions, charge balance and mass balances for key components in the system were considered,

Introduction
Inorganic scales are one of the main flow assurance problems in the water and oil and gas industry or everywhere incompatible fluids can mix up.Inorganic scales can reduce the flow area of pipelines leading to increased pressure drop, or even blocking of the pipeline.Scale deposition can occur in the near wellbore area, resulting in a reduction of permeability and consequently reduction in production.Any remedial treatments involve well shut-in, increasing the loss of production, resulting in the nonprofit time increase of the pipeline.Severe scale deposition may result in losing production wells, which may result in loss of revenue in millions of dollars.So, mitigation of the scale problem is much more than favorable for operator companies to prevent the scales from building up rather than remediating them.Generally, as the water cut increases, the likelihood of scale deposition increases [33].As most of the oil and gas reservoirs are mature ones, their water cut is increasing and at the same time water injection processes are causing higher water cuts more feasible to occur resulting in inorganic scales more likely to occur.
One of the main inorganic scale mitigation methods is squeeze treatments.In squeeze treatments, Scale Inhibitors (SI) are injected into the formation in a concentrated batch followed by an over-flush stage to push the SI slug deeper into the reservoir formation [8].Scale inhibitors are usually weak poly acids, denoted by H n A for an n-protic system, that can be speciated [23], and can create complexes with the cations that are the building block of the inorganic scales and in this way remove these ions from the system preventing the scales to build-up [13,20].Scale deposition will be prevented, as long as the concentration of SI in solution is above a certain threshold concentration, known as minimum inhibitor concentration (MIC).The squeeze treatment lifetime is determined by the time or produced water volume when the SI concentration falls below the MIC which results in another treatment requirement.Thus, the treatment lifetime is governed by the retention level of the SI in the rock formation which depends on the rock, fluid, and operating condition of the system.So, the performance of the SI retention in the system is maintaining the SI concentration above the MIC for longer periods [18,28,34].
Retention of scale inhibitors is governed by two main mechanisms of adsorption, Γ(C), and precipitation, Π(C) [5,24,33].The overall efficiency of the retention process is the result of the coherent action of these two mechanisms.A large number of previous studies have shown that in most cases, both mechanisms can actively contribute to the retention and a coupled adsorption/precipitation mechanism determines the retention process efficiency [5,19,24,26,32,33].
In the adsorption mechanism, the SI molecules are retained on the rock surface through electrostatic interactions (adsorption mechanism), which is the main mechanism in cases where a low concentration of divalent cations exists in the system [5,6,14,25,35].Adsorption can be described by an adsorption isotherm, described by various empirical expressions such as Langmuir, Freundlich, and BET, denoted by Γ(C), where C is the [SI] that exists in the solution phase.Based on the adsorption models, as the [SI] in the solution decreases, the adsorption amount is decreased, and in this way, which means some of the adsorbed SI molecules can be released in the solution and support the [SI] in the system.
The other mechanism is precipitation.SI-dissociated species can create complexes with divalent cations such as Ca 2+ and Mg 2+ , which usually exist in the formation and injected waters and also can be released from the reservoir rock in carbonate formations [22].A portion of the complexes will precipitate based on their solubilities.This part can be redissolved again as the [SI] in the solution decreases and participates in providing [SI] in the solution.Several studies in the literature reported that precipitation may significantly increase the squeeze treatment lifetime by increasing the mass of SI retained in the porous media.Precipitation mechanism is the dominant mechanism in systems that are rich in divalent cations such as Ca 2+ and Mg 2+ , for example in carbonate reservoirs [2,12,16,28].In systems where divalent cations concentration is significant, adsorption and precipitation determine the retention of SI in the porous media [16,27].
In terms of the SI interactions in porous media, lots of experimental investigations have been conducted in the literature.Most of these investigations focus on the effect of different parameters on the retention of SI in different mineralogizes [4,7,9− 11,13,14,16,18,22,27,31,33].Some of these studies are mainly focused on the SI interactions in reactive formations, especially carbonates.The interaction of SI in reactive formations is much more complex than in sandstones, since they are actively involved in SI dissociation and complexation, by pH changes, and supplying divalent cations to the system.In this regard, the effect of different operational parameters in this system such as pH, temperature, water composition, and solid surface area on the retention of SIs in carbonates have been experimentally investigated [7,9,11,36].With particular focus on the reactive formations, some investigations have researched the stoichiometry of the SI in these systems [23,29,30].
In terms of mathematical modeling, most of the studies have focused on the squeeze treatment design and transport model in the porous media.A few studies that have focused on the prediction of SI concentration considering the different mechanisms of the retention based on the adsorption isotherms (Γ) and precipitation (Π) [3,4,17,24,26,32].
Despite these extensive experimental and modeling investigations, there is not any detailed reaction-based modeling of these interactions in carbonate systems.It can be due to this fact that lots of different reactions can occur in these systems with a particular huge set of reaction constants for each system.However, there are some experimental studies that have investigated and modeled some of the possible reactions and their associated constants in particular systems in this category which consider only some part of the full system (i.e.speciation of DETPMP or Complexation of Ca-DETPMP) [21,23,29,30].On other side, there are not any accurate reported investigations for the precipitation characterization of SI complexes in these systems which directly controls the precipitation part of the process.This can be due to the fact that precipitation nearly always occurs with adsorption in such systems and identifying them is not an easy and straightforward task and usually are considered together with adsorption; So, SI/carbonate/divalent brine systems are not very well characterized for all the possible complexes which can occur [23,31,35].
In this study, a comprehensive model is developed by considering all the reactions that might occur in a general phosphonate SI-carbonatebrine system containing divalent cations.Also, the impurity reactions are considered along with the main reactions and the whole system will be solved by considering the mass and charge balances in such a system which will result in the characterization of the system for solution equilibrium.Following that, a coupled adsorption/precipitation isotherm is introduced to model the amount of adsorption and precipitation in different concentrations and mass to volume ratios.The coupled adsorption/precipitation isotherm is then considered along with the solution equilibrium model to characterize the full equilibrium of phosphonate SI-carbonate-brine system with the adsorption and precipitation which ends up with determining the concentration of all species in the system in any form of solution, adsorption, and precipitation states.

Model overview
In this study, a mathematical chemical equilibrium model is presented to predict the reactions of phosphonate SI with carbonate formations (specifically CaCO 3 in this work), along with the accompanying SI retention mechanisms in carbonate formations.A schematic view of this overall process is shown in Fig. 1 in which a volume (V) of SI solution of initial concentration C 0 = [SI] 0 at a given initial (adjusted) pH = pH 0 with initial divalent concentrations of [Ca 2+ ] 0 and [Mg 2+ ] 0 are added to a given mass (m) of calcium carbonate (CaCO 3 ), as shown in Fig. 1(a) at t = 0. If, after some time, this system goes to equilibrium in the absence of either adsorption (Γ) or precipitation (Π), then the situation is reached in Fig. 1 This is referred to below as Modelling -Part 1.
As discussed in detail above, scale inhibitor retention in the porous media is governed by both adsorption and (Γ) and precipitation (Π) processes.Therefore, in a total model of the system, the next processes/ reactions that must be included are shown in Fig. 1(c) and Fig. 1(d); that is, either there may be only pure adsorption of SI, denoted Γ (as in Fig. ((c)), or coupled adsorption and precipitation may occur (coupled Γ/Π) as shown in Fig. 1(d)).This is referred to as Modelling -Part 2.
In solving the entire equilibrium system numerically, Parts 1 and are solved sequentially, where the initial equilibrium of the soluble system (Part 1, Fig. 1(b)) is established.Then, a procedure is implemented to solve for the adsorption/precipitation (Γ/Π) processes (Part 2, Fig. 1(c) and Fig. 1(d)), and this is repeated iterative between Parts and 2 until convergence is achieved.

Modelpart 1: the carbonate, Ca 2+ /Mg 2+ brine, SI interactions (Γ=Π=0)
To model the first part of the process, it is noted above that we need to couple together the equilibrium equations of (i) the aqueous carbonate system, (ii) the SI dissociation, modelled as a weak polyacid (H n A), and (iii) the Ca 2+ and Mg 2+ binding (chelation) to the various SI dissociated species.As a practical point for phosphate SI inhibitors, we also model the presence of "impurities" that may occur in commercial phosphonate products, as occurs in most industrial chemicals.A brief outline of the 3 sets of equations is given below, but a technical description of the derivation of the final equations is given in Appendix A.

The carbonate system
In this equilibrium system, the first set of chemical reactions that exist is the carbonate which contributes to the final pH, and the concentration of Ca.The pH in turn determines the speciation of SI and the level of Ca 2+ (and Mg 2+ ) in solution determines the level of complexation.The reactions associated are given below.The important parameters for this system are the equilibrium constant of the reactions and the solubility product of the calcite which all of them are well characterized in the literature [1,29].
Fig. 1.SI-Ca-Mg system reaction sequences with subsequent adsorption and precipitation.Units: Γ and Π in mg/g; C eq = [all SI species in solution] in mg/L; Change in concentration of SI due to Γ and Π is C Γ = m Γ/V and C Π = m Π/V which will be used below.This is essentially an "apparent adsorption" experiment; i.e.Γ app vs.
[SI T ] f.where the SI T refers to all species containing SI in solution, including free SI dissociated species and those chelated to Ca and Mg. (3)

SI speciation/dissociation
Phosphonate SIs are considered as weak poly acids that can be speciated, releasing H + protons, which will participate in the acidity development of the system, as well as the negative species which can be bound with divalent cations [23,29].This is a similar mechanism to that by which SIs prevent scale precipitation, binding to the scale-forming cations in nucleation, or in the growing crystal.It should be noted that SIs can only be speciated through their weakened H + in OH -groups in their structures; therefore, only a limited number of H + in their formula can be speciated, and because of that, in this perspective, they are represented as the H n A, where the n is the number of H + protons that can be speciated.In this study, DETPMP is considered which is denoted as H 10 A (n = 10), to represent the chemical structure, with 5 phosphonate groups with two hydroxyl groups that up to 10 consecutive speciation may occur, with the corresponding reported speciation constant [29]; see Fig. 2.However, the model was developed for a general case that SI can be speciated n times with corresponding speciation constants as follows.Each reaction in the successive dissociation has an associated equilibrium constant, in this case, a dissociation constant and there are n such constants (see Appendix A).

Complexation of SI with divalent cations, Ca 2+ and Mg 2+
SI dissociated species form complexes with divalent cations in the brine, most commonly Ca 2+ and Mg 2+ .These cations usually exist in the formation brine, but can also be released, due to the dissolution of rocks containing these minerals, such as carbonates; e.g.calcite CaCO 3 or dolomite CaMg(CO 3 ) 2 .These chelation reactions are part of the first steps in the precipitation mechanism, which contributes to the SI retention in the formation.Detailed data on these chelation interactions is quite scarce in the literature [16,29] and these reactions are not very well characterized.The main uncertainties relate to the values of the individual stability constants for the Metal (M)-Phosphonate (Ph) groups binding.However, previous studies do provide some approximate relative values for these reactions, particularly the relative strengths of the M-Ph binding [23] particularly for the relative strengths of the Ca-Ph and Mg-Ph links.The binding of SI negatively charged species and positive cations (Ca 2+ and Mg 2+ ) is described by the following complexation reactions: (9) The detailed way these binding constants are treated is shown in Appendix A.

Treatment of SI impurities
Commonly, commercial scale inhibitors have some degree of impurity, due to the industrial processes involved in their manufacture.As previously reported in the literature, for example, DETPMP may contain about 10% of phosphorus-containing compounds that are not DETPMP (e.g.H 3 PO 3 /H 3 PO 4 ) [20].However, these compounds contribute to the "P" signal that is measured by ICP in the "SI assay".This phenomenon would also be the true case in experimental works if the SIs have been used in "as supplied" status.Ideally, equilibrium experiments that are performed should use purified phosphonates to remove any effects from these impurities.However, to perform accurate modeling of the process being studied here using commercial products, these impurities should also be accounted for in our equilibrium equations.
For DETPMP, the impurities are mostly phosphorus and phosphoric acids (H 3 PO 3 /H 3 PO 4 ) [20], which contain phosphorus in their structure and have an effect on the experimental results in two ways.First, they contribute to the phosphor concentration in the solution, which would be wrongly attributed to SI in ICP analysis, thus overestimating the real concentration of SI-associated P. In addition, as the H 3 PO 3 /H 3 PO 4 are acids, and therefore they can influence the acidity/alkalinity of the system through the speciation process and can therefore disturb the equilibrium.
The first effect is much easier to account for since it only changes the amount of the SI concentration in the initial and final state and in the mass balances of the SI.The other effect that should be considered is the contribution of the P-acid impurities through the acidity/alkalinity of the system, which alters the system equilibrium.To consider the acidity/ alkalinity effects of these impurities, their speciation reactions must be included in our calculations, and these are given below.The speciation process of these acids is very well known, and accurate speciation constants for these reactions are available.In the literature, it has been reported that phosphorous acid is an H 2 A polyacid and can release two H + into the solution as follows: Phosphoric acid is an H 3 A polyacid and releases 3 H + as follows: The equilibrium constants for all these dissociations are given in Appendix A.

Coupling of the complete reaction set -Part 1
We then couple the 4 sets of chemical equilibrium equations discussed above; via the carbonate system, SI speciation/dissociation, SI complexation with Ca, competing SI complexation with Mg, and the SI impurity equations.The objective is to find the equilibrium state of this entire set of equations subject to any given set of initial conditions.As in all correctly posed large systems of chemical equilibrium there are a large number of independent species, N s , but rather fewer equilibrium equations.However, when the set of equilibrium equations is combined with the full set of mass balance equations and 1 charge balance equation, then the final number of equations is guaranteed to be N s .Specifically, the mass balances are of the C, SI, Mg, Ca, phosphorous species, and phosphoric species.This formulation yields a set of N s nonlinear equations which fully characterize this equilibrium condition.The full formulation is very lengthy, but it may be reduced systematically such that the entire system depends only on 3 primary variables, z ].After some algebra operation and substitutions, three final dependent equations (F 1 , F 2 , F 3 ), each a function only of the 3 main variables (z, m 1, and m 2 ) can be derived in conjunction with 5 ancillary equations.Solving this set of equations requires the solution of F 1 = 0, F 2 = 0 and F 3 = 0.These equations are as follows (see Appendix A for derivation): From the equations above (and Appendix A), it is evident that the entire system of equations is dependent on the 3 key concentrations, z = [H + ], m 1 = [Ca 2+ ] and m 2 = [Mg 2+ ].The 3 master equations (F 1 = 0, F 2 = 0 and F 3 = 0) can be solved by Newton Raphson, in conjunction with some boundary checks to prevent the divergence and unrealistic solutions.Once the z, m 1 and m 2 are calculated at equilibrium, all other species may be calculated by back substitution, yielding a full characterization of the equilibrium system.
It should be noted in this system, the following assumptions have been considered: -The volume of the solution is not changed due to the reactions (dissolution, precipitation, complexation, and speciation) -The conditions, such as temperature and pressure are assumed constant and hence the equilibrium and solubility constants do not change.It is straightforward to relax this condition but the various equilibrium constants, dissociation constants, solubility products, etc., would be required as functions of temperature (mainly) and pressure.

Model -Part 2: coupled SI adsorption/precipitation (Γ/Π)
The second part of the fully coupled model development, the adsorption (Γ) and precipitation (Π) of the (now complexed) SI species must now be included in the model.Since both mechanisms act toward removing some of the chemical species from the solution, then this will generally lead to a new equilibrium that could recursively change the equilibrium in Part 1 described above, thus some iteration will generally be required.

The coupled adsorption/precipitation (Γ/Π) isotherm
To address coupled adsorption/precipitation (Γ/Π) in this work, we use the concept of "apparent adsorption", Γ app [10,14].In an apparent adsorption experiment, we imagine the system in Fig. 1-a going to equilibrium as in Fig. 1-d.If we were to calculate all of the material lost from the solution in this experiment as if it was adsorption (although it is part Γ and part Π), then we would calculate an "apparent adsorption", Γ app , as follows: Many experimental results have been published on scale inhibitor/ carbonate systems calculating apparent adoption [10,14,2,27,[5][6][7].A schematic example of an apparent adsorption experiment is shown in Fig. 3-a, and many experimental results on scale inhibitor /carbonate have been published [5,7,10,14,27].In this figure, we note that there are two distinct regimes where at lower threshold concentrations (C < C th ) we observe pure adsorption (Γ), and then at higher concentrations (C >C th ) we observe a region of coupled adsorption/precipitation (Γ/Π).In this figure, we note that results are presented for 2 mass to volume ratios (m/V).Note that for pure adsorption, apparent adsorption data for both of these (m/V) ratios track along the same line, which is the M. Kalantari Meybodi et al.
adsorption isotherm, Γ(C).However, above the threshold concentration, then the two lines deviate for different (m/V) ratios, and this is the regime of coupled adsorption/precipitation (Γ/Π).The reasons for this behavior are fully explained and modeled in previous literature [9,11].
The analysis of a large set of apparent adsorption data shows clearly that this quantity is comprised of 2 additive functions (shown in Fig. 3); the adsorption isotherm (Γ) along with the "precipitation function" (C Π vs. C eq ) describing the purely precipitate region, which the data indicates is a linear function.This is shown in Fig. 3-c and we note that the relationship between adsorption and precipitation depends on the specific case of mass to volume ratio (m/V).Thus, to apply coupled adsorption/precipitation calculations in any actual case requires the underlying functions for Γ and Π: From the definitions in Fig. 1, for any initial [SI] = C i , when adsorption/precipitation occurs, then: However, we indicated above that C Π is the linear function of the C eq as follows: and this is demonstrated from the experimental data below.Therefore, Eq.24 may be reformed as: It should be noted that, in the pure adsorption region (Γ∕ =0, Π=0) for 0 ≤ C ≤ C th , then C Π =0 as shown in Fig. 3(c).For any initial C i (total SI concentration), the value of C eq is determined based on the solution of the following equation for any given m/V ratio, which can be solved quite straightforwardly by any root finding or optimization method.
Thus, the subcomponent functions of the Γ/Π isotherm are (i) the adsorption isotherm, Γ(C), and (ii) the precipitation function, C Π vs. C eq , and for any given (m/V) ratio and initial SI concentration, C i , these can be used as described above in Eq. 27, to solve the system for C eq .This construction is actually shown in Fig. 4 where the math shows that this leads directly to Eq. 27.
The reader may ask: for a given (m/V) ratio why not use the measured Γ app directly, since Γ app = Γ + Π = Γ + (V/m) C Π ?The reason is that this would give a solution for that specific system, but without the 2 independent subcomponent functions, Γ(C) and C Π , then we could not calculate the equilibrium system for any other arbitrary (m/V) ratios.In fact, practical laboratory Γ app bottle experiments tend to work at fairly "low" values of (m/V) ratio; in units of g/L; typically, m ≈ 10 g, V ≈ 0.04 L, and m/V ≈ 250.Inside porous media, or the grid blocks simulating a flowing system, then typically this ratio can be greater by an order of magnitude or more, say m/V ≈ 3500.This is important when considering a scale inhibitor/carbonate system flowing within a porous medium, but this will be dealt with in a later paper (Meybodi et al., in preparation).

Determining the coupled Γ/Π subcomponent models experimentally
In this section, it will be shown how the various subcomponents (Γ(C) and C Π ) of the coupled Γ/Π isotherm are determined from a previously published experimental data set [9].From these previous experiments, the apparent adsorption data are given in Fig. 5. From these experiments, there is a clear dependency of the apparent adsorption on the mass to volume ratio (m/V).To have a unified model for the DETPMP-Calcite-NSSW system that can be used for all ranges of the concentrations at the given temperature, then we must derive these subcomponents of the coupled Γ/Π isotherm model.
To carry out this task, first, the adsorption isotherm, Γ(C), should be determined from the experimental data.As shown in the literature [9], in the adsorption-precipitation of SI, there is a region of SI equilibrium concentration where adsorption is the sole mechanism of SI retention (at very low concentrations usually below a threshold level of C th < ~ 100 mg/l of equilibrium [SI]).This region is characterized as the region in which the apparent adsorption (Γ app ) graph is not a function of solid mass to the solution volume(m/V).Using the experimental data, the Langmuir adsorption isotherm was matched for the pure adsorption region, as shown in Fig. 6, which for this example resulted in Γ max =0.1321 and K=0.0773 for the Langmuir isotherm model and the full adsorption isotherm based on these values is given in Fig. 7.
The adsorption isotherm based on the experimental data, Γ(C), may then be used to find the precipitation function (C Π vs. C eq ) part of the coupled Γ/Π isotherm.This is done quite simply by subtracting the adsorption part of the apparent adsorption (i.e., Γ(C)) from the measured Γ app ; since Γ app = Γ + Π, then this leaves the Π (precipitation) part in units of mg/g, which can immediately be translated to C Π = mΠ/ V.When this procedure is carried out using the experimental data available to us, the results are shown in Fig. 8; this clearly shows straight line relationship between C Π and the equilibrium SI concentration, C eq .For the data in Fig. 8, this resulted in the fit, a=2.1399 and b=-234.1626as the precipitation matching parameters.
To check the reliability of the experimentally determined subfunctions of the coupled isotherm, the results for the recalculated apparent adsorption experiment are compared for the model with the reported experimental values in Fig. 9.These show good agreement with the unified coupled Γ/Π isotherm in predicting the apparent adsorption behavior of the coupled adsorption/precipitation system for different mass to volume ratios.

Coupling the part 1 (Equilibrium system) and part 2 (Γ/Π submodels) into the total model
Firstly, the Part 1 set of Eqs.20 to 21 are solved using the Newton-Raphson numerical method.If there were no adsorption or precipitation, then this would solve the coupled scale inhibitor, carbonate, and Ca/Mg chelant system.However, when both adsorption and precipitation occur, we must then resolve Part 2, the adsorption/precipitation (Γ/Π) part of the problem.
We then use the coupled Γ/Π isotherm subfunctions (Γ(C) and C Π ) at the given (m/V) ratio and initial SI concentration, C i , to calculate the C eq , C Γ , and C Π through an iterative approach.Then the concentration of the species that should be adsorbed or precipitated is identified, as follows.
(i) Adsorption of all species containing any dissociated or complexed (with Ca and Mg) SI may take part in the adsorption process.The amount of this allowed adsorption is the total molar amount calculated by the adsorption isotherm Γ(C).We assume the proportion of each adsorbed species in the total adsorbed amount is exactly the same as the calculated soluble fraction from the Part 1 calculation.(ii) Precipitation of only species that are bound to Ca and/or Mg are assumed to precipitate; no free SI species unbound to Ca or Mg may precipitate.As in the adsorption case, the proportion of species precipitating is taken to be identical to those calculated in solution.
After the calculation of the concentration of each species in the adsorption and precipitation phases, the remaining SI concentration in the solution is recalculated considering the reduced SI concentration due to the adsorption and precipitation processes.Since this solution has been altered by adoption-precipitation, it may be slightly out of equilibrium and thus must be recalculated.So, the equilibrium will be performed once again for the remaining species in the system, and the adsorption-precipitation is calculated based on the new equilibrium, and adjustments are conducted for the previous adsorptionprecipitation.This iterative methodology continues to satisfy the following 2 equilibrium conditions simultaneously which guarantees the whole system's equilibrium: -Equilibrium of C eq , C Γ , and C Π as the main parameters of the equilibrium system  -Equilibrium of all remaining species in the solution phase This model can be considered as the full general model for considering the coupled adsorption/ precipitation of SI in carbonate formations.Slight changes in the above procedure may be required to model particular ways of carrying out these experiments in the laboratory, and these will be explained for specific cases.

Results and discussion
In this section, the developed model based on the methodology was investigated for its validity and reliability by reproducing certain simple cases and by direct comparison with experimental results.

Model input data
The following data were used in the model for the various calculations which follow.In most of the calculations and all of the experimental conditions, North Sea sea water (NSSW) was used as the water composition which has the composition given in Table 1.Also, the initial pH of NSSW was considered to be pH = 5.7.
The scale inhibitor (SI), Diethylene Triamine Penta (Methylene Phosphonic Acid) (DETPMP) was used in all experiments and calculations, since this is a widely applied commercial product.DETPMP can be protonated by up to 10 H + and is denoted as H 10 A. Based on the previous studies in the literature [29,30], the values in Table 2 were used as the dissociation constants for each of the speciation steps.All equilibrium speciation, dissociation, and solubility constants are taken at standard T and P conditions.This is a data limitation rather than a limitation of our model and few of these constants are known at higher temperatures.
In model calculations, pure calcite (CaCO 3 ) was considered as the carbonate substrate.The set of equilibrium constants in Table 3 was used for the calcite dissolution and carbonic system reactions.
Regarding the DETPMP impurities, H 3 PO 3 and H 3 PO 4 have 2 and 3 acid speciation constants respectively with the values given in Table 4.The other parameters that are required for the model utilization are the stability constants for the equilibrium formation of the SI species with divalent cations of Ca 2+ and Mg 2+ .However, there are not any accurate and reliable values reported for these quantities, but some approximate values have been given for them in the literature [23,31].The values of 3×10 8 L/mol and 1×10 8 L/mol were considered for the formation of Ca-SI and Mg-SI complexes respectively.
In the first instance, the model validity has been checked previously [15] regarding some of the well-known cases (carbonic system, DETPMP dissociation) in the literature.These checks can confirm the generality of the model for both cases of just equilibrium and equilibrium with adsorption and precipitation.

Technical details of the modelled experiments 3.2.1. Overview of experimental procedures
For the validation of the proposed model, two sets of experimental data were used.Both sets of experiments were conducted in a very similar manner with slight modifications in procedures that make them appropriate for validating the model in two different situations.The first set of experiments that has been used for model validation in the case of no-adsorption, and no-precipitation has been conducted recently and is not yet published; in Section 3.3 below, this provides the data for Experiments 1 and 2. The other set that was used for the adsorptionprecipitation equilibrium model validation has been published previously [9]; this provides the data for Experiment 3 in Section 3.3.Before utilizing the model to simulate these experimental data, it is more favorable to have a brief overview of the experimental procedure that was used in both experiments to understand experimental results and the way that they will be modeled.
Both sets of experiments have been conducted following the procedure given in the study of Jarrahian et al. [9].This study investigated the precipitation behavior of the DETPMP-NSSW in the presence of calcite.North Sea Seawater (NSSW) was prepared based on the given composition in Table 1 by dissolving appropriate salts and then filtered through a 0.45 μm filter paper.Then a 10,000 mg/L SI stock of DETPMP (active part) was made up in the NSSW and a series of the SI stock solutions were created based on the mass balance for different SI solutions of ([SI] = 0, 50, 100, 300, 500, 1000, 2000 mg/L.Afterwards, the pH of these solutions was adjusted by introducing HCl and NaOH to the desired initial solution pH of 2, 4, and 6 in different stock solutions.Then, the initial [SI] in these solutions were measured using the Inductively Coupled Plasma Optical Emission Spectroscopy (ICP-EOS) method.Then small portions of these solutions (40 ml) were used for the main tests with accurately weighted 5 and 10 gr of calcite.The test containers were placed in the oven at 95 • C.After 24 hr, they were removed and the content was filtered through a 0.22 μm filter paper and the supernatant was placed at the ambient condition to be cooled down to 25 • C for 24 hr, and then the final pH of the supernatant was measured as well as [SI], [Ca 2+ ] and [Mg 2+ ] using ICP-OES.Then the amount of [SI] "missing" from the solution was used to calculate the apparent adsorption (Γ app ) which accounts for the sum of the adsorption (Γ) and precipitation (Π), as described above [7,9,11].This set of experiments was used for model validation in the case of precipitation-adsorption.
The other set of experiments, which was conducted recently, followed the same procedure with two small modifications.The first and more important one was that the experiments were carried out at temperature of 25 • C. In these experiments, no-precipitation was observed which makes this set of experiments more appropriate for the noadsorption, no-precipitation model verifications.Also, in the new experiments, the pH of the DETPMP-NSSW was measured before any pH adjustment.It made these results fit the model verification process where no calcite exists under the scheme of no-adsorption, and noprecipitation.

Using experimental data for model validation
To apply the proposed model to the experimental data, the initial scale inhibitor concentration, [SI] initial is used as an input of the model along with the NSSW composition given in Table 1 along with the various equilibrium and reaction constants given in Table 2, Table 3 and Table 4. Then following the procedure outlined in the methodology section, the [SI] final is used to match the final SI concentration in the model and the calculated values of [Ca 2+ ], [Mg 2+ ] and [H + ] are found; these are then compared with the corresponding values from experiments.

Modifications of model for pH adjustment
In the experimental procedure, the pH of the solution is adjusted to a fixed value (pH 2, 4, or 6) after the dissolution of SI in NSSW just before adding calcite to the system.Exactly this procedure must be applied in the model; i.e. the SI is speciated within the aqueous phase thus finding its own natural pH, and this is "numerically titrated" to the fixed starting pH.This is necessary since the titration to the starting pH changes the total acidity/ alkalinity of the system.Indeed, the intermediate pH when the DETPMP at a particular concentration is first added to the brine, is an important check on the model for that stage of the experiment (and this is shown below).Then using the new initial pH after adjustment, the full equilibrium including the calcite is established to find the final equilibrium condition of the system; these stages are carried out both experimentally and in the numerical model.

Direct comparison of the coupled model with experiment
In this section, we compare the model predictions with 3 experiments of increasing complexity as follows: (i) Experiment 1: when DETPMP is made up in brine (NSSW here) at different SI concentrations, then different pH values result.The first test of the model is to predict these pH values.(ii) Experiment 2: the solution in experiment 1 is pH adjusted and then solid CaCO3 is added to it, but this is done at room temperature such that no precipitation occurs, and adsorption is (practically) Γ = 0.The test here is to predict the final solution pH values and also the solution Ca 2+ and Mg 2+ levels.
(iii) Experiment 3: The experiment described in (ii) above is then repeated but at higher temperatures where precipitation and adsorption also occur.This involves the entire model, and the test here is to predict the final solution pH values, the final Ca 2+ and Mg 2+ levels, and also the amount of precipitate (and adsorption).

Model validation for the SI-NSSW pure solution case -Experiment 1
Firstly, the model is tested for the DETPMP-NSSW system for the solution only (no CaCO 3 ) by comparing our predictions with the pH experimental measurements before adjustment.For this purpose, various concentrations of the DETPMP were dissolved in the NSSW and the pH of the equilibrated solutions was measured.As shown in Fig. 10, a very good agreement has been achieved between model and experimental results; the pH values agreed with an absolute average relative error of 2.95%.This confirms the model validity for the SI-Brine equilibrium system determination without any carbonate interaction or adsorption/precipitation present.

Model validation for the SI, NSSW-CaCO 3 solution case (no adsorption or precipitation) -Experiment 2
A comparison of the model predictions and experimental results was then carried out for the case where pH was adjusted to pH = 4, with the addition of m = 5 g and 10 g of calcite.Note that in this process, the model had to carry out a numerical titration of the native pH at each concentration (as in Fig. 10) to the test pH = 4. Fig. 11 shows the pH values predicted by the model compared with the experimental values.Good quantitative agreement between model and experimental is found, with an absolute average error of ~5%.It should be noted that the final solution pH is the result of competition between many reactions that are occurring simultaneously, as described in Section 2.2 above.Some of these interactions behave in the direction of decreasing pH, such as SI dissociation or the presence of SI impurities, and some of them behave in the opposite manner, such as calcite dissolution.Thus, the final predictions in Fig. 11 are actually the result of several competing reactions as listed in the modelling section.
Result in Fig. 12 for these calculations, where there is no adsorption or precipitation, a large excess amount of calcite exists, it is predicted that quite similar amounts of calcite are dissolved during the equilibrium process for the m= 5 g and 10 g cases.The actual quantity of CaCO 3 dissolution, ~25 mg, would only be ~0.5% or 0.25% of the m=5 g and 10 g cases respectively, and this would be difficult to measure experimentally.
Fig. 13 shows predicted and experimental final calcium levels for this same case where again good agreement is observed between the model and experiment; with an absolute average error of ~5.5%.The calcium level in the investigated system is controlled by calcium-releasing reactions such as calcite dissolution, and calcium consumption reactions like complexation reactions of SI-free species with calcium cations.Fig. 14, shows the agreement between model and experimental values of the final solution magnesium levels for the same case; these agree with an absolute error of ~5.9%.Qualitatively, Mg behaves in the same way as Ca except there is no Mg source from the rock which is CaCO 3 .(However, if dolomite, CaMg(CO 3 ) 2 , were used then this can easily be accounted for in the model).In this case, the magnesium level does not change in the system and is simply redistributed from the free state at the start of the process to the combination of free and SI-Mg complexes at the equilibrium conditions.

Predicted speciation of the Si-Ca-Mg complexes -Experiment 2
An important aspect of SI reactions with divalent cations (Ca 2+ and    Mg 2+ ), is the distribution of the various species in equilibrium since these determine the overall stoichiometry of the formed complexes (e.g.SI-Ca N1 _Mg N2 ) and hence the consumption of divalent cations.For example, if the dominant species in the system is H 8 A − 2 , then only 1 divalent cation will be engaged per SI molecule in the complex, but if the dominant species are H 2 A − 8 , then 4 divalent cations will be involved per SI molecule in the complex [23].
Recently, Meybodi et al. ( 2023), demonstrated using a preliminary version of the model presented here, that at any final equilibrium pH, then only 2 or 3 groups comprise the dominant range of dissociated species (from H n A).Hence, the final set of SI_M 2+ complexes will be formed related to this same distribution of species.To illustrate this point, Fig. 15 shows the distribution of the free SI species for different equilibrium pH values resulting from the case for m = 10 g calcite.From the SI dissociation curve [9,15], we expect the species of H 3 A − 7 , H 4 A − 6 , and H 5 A − 5 with the dominant species being H 4 A − 6 for a final pH ~ 7; this is confirmed by the results in Fig. 15.The slight variations of this distribution as the final pH is lowered slightly are shown in Fig. 15 for the values of pH = 6.46 and 6.31.The free unbound species of the SI also show the same distribution of charged species as in Fig. 15, but as the stability constants are very high (~10 8 ) for the complexation reactions, then the concentrations of free SI species are nearly 6 orders of magnitude lower than the corresponding SI_M 2+ complex concentration.
The overall stoichiometry of the SI_M 2+ complexes was calculated based on the model results as shown in Fig. 16; i.e. the N1 and N2 in SI-Ca N1 _Mg N2 .From the literature [22], for the pH~7 and [Mg 2+ ]/[Ca 2+ ] ~3, it is expected to have nearly 2 magnesium and 1.2 calcium cations on average for each SI which confirms the results given in Fig. 16.So, at a pH of 7, it is expected to have the stoichiometry of SI_ Ca 1.2 _Mg 2 , and this is the range predicted from the model in Fig. 16.

Model validation for the SI-Ca/Mg/CaCO 3 /adsorption, precipitation (Γ/Π) case -Experiment 3
The results for the modeling of the experiment, where coupled adsorption-precipitation was observed (Experiment 3), are now presented.To model this system, we must invoke every part of the entire coupled model described as Part 1 and Part 2 above, i.e. the coupled carbonate/ SI dissociation / Ca-Mg chelating system along with the adsorption/ precipitation (Γ/Π) model.Using the derived coupled isotherm in the methodology, the coupled model was applied to model the experimental data with initial pH=4.Fig. 17 shows the model predictions for pH for Experiment 3 which includes adsorption and precipitation.The results showed that the model broadly catches the true trend of experimental data.The average model accuracy for the pH was ~ 8.2%.
Fig. 18 and Fig. 19 present the predicted and experimental final solution calcium and magnesium levels for this system.It shows that the model predicts the calcium and magnesium levels with the average absolute error of ~5.3% and ~4.72%, respectively.The increasing trend in the calcium level in Fig. 18 is due to the changing solution pH and availability of a high amount of calcite that can provide calcium to the solution whereas for magnesium, there are no sources to release this ion.Therefore, just some of the magnesium is consumed in the adsorption and precipitation case and the level of the magnesium decreases slightly (by incorporation into the SI-Ca-Mg precipitate).
The model predictions of the complex stoichiometries are presented in Fig. 20.For the same conditions of the previous case with equilibrium pH ~7 and [Mg 2+ ]/[Ca 2+ ] ~3, the average stoichiometry of Ca 1.16 -Mg 2 -SI is predicted with the model which is consistent with the previous case.This is due to the fact that equilibrium pH controls the speciation of the SI and the dominant species which can remain in the solution or be removed from the solution through precipitation and adsorption.
An important consideration in the adsorption-precipitation  equilibrium case is the distribution of the SI species between adsorption, precipitation, and solution.Fig. 21 shows this distribution for one of the investigated cases.As was shown previously, for the lower concentration of SI, the adsorption is the only mechanism that is active to take the SI out of the solution.However, after a certain concentration of SI, the precipitation starts to act along with the adsorption taking the SI from the solution or in other words, SI retention in the medium.The reason for the decreasing effectiveness of adsorption in SI retention with SI concentration, is that the Langmuir adsorption model was considered to model the adsorption which has an ultimate adsorption limitation, and from a concentration, the amount of adsorbed SI remains nearly constant and further removal is accounted by precipitation.
Finally, the distribution of the different species in 2 different cases from two different equilibrium SI concentrations and equilibrium pH is given in Fig. 22.This is completely aligned with the distribution of the species in different states based on the equilibrium SI concentration and associated mechanisms for SI retention as well as the dominant species that occur in different conditions.
Based on the results presented in this section, the model shows the good agreement with the experimental data in different situations, which confirms the validity of the model for the solution equilibrium and also the coupled equilibrium processes.This model geochemically characterizes the reactions and processes that are occurring in the SIcarbonate-brine system, which is the case in a wide range of squeeze treatments that are deployed in the real field operations.Therefore, it can be used in conjunction with a transport model for accurate modelling and simulation of these squeeze treatments, and this development is currently in progress.These simulations of squeeze treatments, employing this coupled model, will be equipped with some new features such as the prediction of pH, and stoichiometry profiles of the return stream in addition to the return SI concentration profile.Prediction of pH and stoichiometry of produced stream is of particular importance for engineering designs and treatments like corrosion control considerations.
In addition, the coupled Γ/Π isotherm, which was presented in this study, can describe the retention of SI in systems where both adsorption and precipitation are participating in SI retention.This model can be used to scale the results of a system to different m/V ratios.For example, this allows us to scale up the bottle test results to the results in porous media which have an order of magnitude difference in m/V ratio.This enables this model to be used as part of the coupled equilibrium model, for field squeeze treatment design and simulation.

Summary and conclusions
The main developments in this paper are summarised along with the detailed technical conclusions as follows: 1.This paper presents the first complete geochemical equilibrium model for the phosphonate SI-carbonate-divalent brine system coupled to the subsequent SI adsorption on the carbonate substrate and precipitation of the SI-Ca/Mg complexes that form.This model may be used to model all of the processes that occur in a wide range of SI squeeze treatments in carbonate formations.2. To develop this geochemical equilibrium model, a base solution equilibrium model was first developed to characterize the solution.
The equilibrium model first formulated all engaging reactions in the system, including the carbonic system, SI speciation, SI complexation, and reactions associated with impurities.These were considered and coupled together using mass and charge balance equations to derive the full equations system.This system was ultimately    ] which could be solved numerically.Other equilibrium concentrations in the solution could then be calculated from these key concentrations yielding a full characterization of the equilibrium system.3.Then, in order to incorporate the adsorption and precipitation of SI in carbonates, a new methodology was introduced to characterize the coupled adsorption/precipitation (Γ/Π) action of SI in such systems.The coupled Γ/Π process is described by 2 additive subfunctions, the adsorption isotherm (Γ(C)) and the (linear) precipitation function (C Π vs. C eq ).These can be solved for any given system's initial SI concentration (C i ) and mass to volume ratio (m/V).Coupling this Γ/Π model with the developed solution equilibrium model can be used to fully characterize the coupled equilibrium system by determining the concentration of all the species in any state of solution, adsorption, and precipitation.4. The model was then validated using a series of experimental datasets of increasing complexity where: (i) the solution pH was predicted and compared with a series of DETPMP concentrations made up in brine (NSSW here) at different SI concentrations (Experiment 1); (ii) the SI solution is pH adjusted and then solid CaCO 3 is added to it at room temperature where no precipitation or adsorption occurs.The model then predicts the final solution pH values and also the solution Ca 2+ and Mg 2+ levels (Experiment 2), and finally (iii) the previous solution experiment is then repeated but at higher temperatures where precipitation and adsorption also occur.The model then predicts the final solution pH value, the final Ca 2+ and Mg 2+ levels, and also the amount of precipitate (and adsorption).These model/ experiment comparisons successively test each part of the coupled model.Good agreement was found between the model developed in this work and experiments for all 3 levels of testing described above.Furthermore, the solution chemistry including the distribution of the species in different states and stoichiometry of complexes were also predicted using the full model and compared with the literature.Again, good agreement was found and this also confirms the validity of the model in predicting these important, but difficult to measure, parameters.
To date, we believe this is the most complete model of the Carbonate/SI/Ca-Mg Brine/ΓΠ system available in the literature.It is being applied to model the existing data that we have accessed over the last 3 decades.In addition, the model establishes a foundational framework for flow models in simulating coupled squeeze treatments in carbonate formations.The model has recently been incorporated into a flow/transport model and this is currently being evaluated and will be the subject of a forthcoming paper (Meybodi et al., in preparation).We aim to apply this full equilibrium/transport code to improve the efficiency of SI squeeze operations in carbonate reservoirs by improving our predictions of reactive SI/Carbonate treatments.In turn, this will reduce operational uncertainties, and lead to optimized treatments and reduced costs by offering a comprehensive characterization of the process.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
(b); i.e.Γ = 0 and Π = 0. Modelling this initial reaction part of the process (to Fig. 1(b)) involves coupling together 3 sets of equations for: (i) The aqueous carbonate system of equations, including the solid calcium carbonate (CaCO 3 ); (ii) The phosphonate scale inhibitor (SI) dissociation equations, treating the phosphonate as a weak polyacid.H n A, which for the main species studied here (DETPMP) would be H 10 A; (iii) The Ca 2+ and Mg 2+ binding equations with the (various) dissociated SI species.

Fig. 9 .
Fig. 9. Apparent adsorption comparison for the coupled isotherm and experimental data.

Fig. 11 .
Fig. 11.Comparison of model pH with experimental data for 5 gr and 10 gr calcite at initial adjusted pH=4.

Fig. 13 .
Fig. 13.Comparison of model calcium level with experimental data for 5gr and 10gr calcite at initial adjusted pH=4.

Fig. 14 .
Fig. 14.Comparison of model magnesium level with experimental data for 5gr and 10gr calcite at initial adjusted pH=4.

Table 3
Equilibrium and Solubility Constants.