Impacts of inorganic salts ions on the polar components desorption efficiency from tight sandstone:A molecular dynamics simulation and QCM-D study

Shun Wang ，Jing Wang ，*，Hui-Qing Liu ，Hong-Quan Zhang ，F(xiàn)ang-Na Liu ，Ronald Omara Erik

a State Key Laboratory of Petroleum Resources and Prospecting in China University of Petroleum，Beijing，102249，China

b MOE Key Laboratory of Petroleum Engineering in China University of Petroleum，Beijing，102249，China

c Research Institute of Shanxi Yanchang Petroleum (Group) Co.，Ltd.，Xian，Shanxi，710075，China

Keywords:Tight sandstone reservoir Smart water Molecular simulation Quartz crystal microbalance with dissipation

ABSTRACT Ions in brine significantly affect EOR.However，the mechanism of EOR with different brine is still controversial.By Combining Molecular Dynamics (MD) method and Quartz Crystal Microbalance with Dissipation (QCM-D) technology to analyze ions distribution and the mechanisms in detaching acidic components on the sandstone，an effective method to determine the detaching capacity was established.The results show that detaching capacity is related to ions distribution and hydration capacity.In the oil/brine/rock system，ions far from the rock are favorable for detaching，while ions near the rock are unfavorable for detaching due to ion bridging effect.The hydrogen bond between water and naphthenic acid is key to detaching.Cations strengthen the detaching by forming hydrated ions with water，and the detaching capacity is negatively correlated with hydrated ions radius and positively correlated with the water coordination number.The detaching determination coefficient was established by considering the ions distribution，ions types，and hydration strength，then verified by QCM-D.The brine detaching capacity with different Ca2+/Mg2+ ratios was predicted based on MD and detaching determination coefficient，and verified by QCM-D.The optimal Ca2+/Mg2+ ratio gave 7:3.This study provides theoretical guidance for targeted regulation of brine composition to improve the recovery of tight sandstone reservoir.

1.Introduction

The wettability of tight sandstone reservoirs is significantly different from that of conventional sandstone reservoirs.For conventional sandstone reservoirs，the formation is saturated with brine before reservoir migration occurs.According to the traditional petroleum geology theory，the oil reservoir is hydrophilic before oil migration(Alvarez et al.，2016).In this case，an oil-waterrock distribution model is formed.There is a brine layer between the polar compounds and the rock surface，and the bridging effect of charged ions in the brine layer is the main function of polar components adsorption on the rock surface.However，tight oil reservoirs mainly migrate over a short distance，and some of them show the characteristics of the source and reservoir symbiosis.With most of the tight sandstone reservoirs in China belonging to continental sedimentary reservoirs，polar compounds existing in the crude oil significantly damage the water film during the oil charging process.As a result，the crude oil in the reservoir is adsorbed on the rock surface mainly through the polar interaction between the functional groups of polar components and the rock surface(Xue et al.，2015;Sun et al.，2014).How to efficiently break the polar interaction between the polar components and the rock surface is key to improving tight sandstone reservoirs'recovery.In recent years，various EOR technologies have been greatly developed，among which smart water flooding is a new EOR method.In order to reveal how intelligent water properties affect the wettability of sandstone reservoirs and determine the potential production mechanism，several laboratory core-scale experiments have been conducted(Tang and Morrow，1999;Chavan et al.，2019;Zhang et al.，2007;Park et al.，2018;Amira et al.，2012;Ahmed and Hisham，2014;Mokhtari and Ayatollahi，2019;Li et al.，2019;Mahmud et al.，2020).For example，using core displacement experiments and Zeta potential tests，Nasralla and Nasr-El-Din(2011)demonstrated that Na+changes the electrical charge at both oil/brine and rock/brine interface to highly negative，which results in repulsion forces between the two interfaces，and hence wettability alteration and improvement in oil recovery.In addition，Na+was shown to be more effective than Ca2+and Mg2+in enhancing oil recovery.Lee et al.(Lee et al.，2010)believe that there is a significant effect on EOR when the injection water salinity of sandstone reservoir is lower than 4000 ppm.However，core-scale displacement experiments require saturation of formation water.This process assumes that the formation was initially water-wet，which is inconsistent with the oil-wet condition of the tight sandstone reservoir.Therefore，the research results can only play a guiding role in regulating the wettability of conventional sandstone reservoirs.

In recent years，QCM-D has been widely used in the petroleum industry as a fine measurement technology (Abudu and Goual，2009;Kelesoglu et al.，2012;Alagha et al.，2013;Farooq et al.，2019)，with a measurement accuracy of up to ng level.In addition，the environment in contact with the quartz crystal surface can be controlled by adjusting the injection fluid order.This particular control provides favorable conditions for studying the adsorption of crude oil in anhydrous conditions and subsequent desorption mechanisms under the action of brine in tight sandstone reservoirs.Arije et al.(Al-Khafaji et al.，2017)studied the influence of seawater with different dilution ratios on the desorption of crude oil on quartz surface via QCM-D technology and pointed out that the increase of Ca2+and Mg2+was not conducive to crude oil removal.They believed that increasing the concentration of Ca2+and Mg2+in the brine leads to the adsorption of acidic components on the negatively charged quartz surface through a strong cationic bridging effect.In addition，they found that seawater diluted 10 times was the most efficient at detaching crude oil from the quartz surface.Further dilution resulted in a decrease in detaching efficiency due to the insufficient repulsive force of diffusion double layer between oil/water interface and the water/solid interface when the salinity was below a certain value.Although significant progress has been made in studying the wettability of sandstone using QCM-D technology，scholars have not reached a consensus on optimal value of brine salinity and ion ratio due to differences in experimental conditions.In addition，the microscopic process of polar components detaching from the rock surface under inorganic salt ions conditions is still ambiguous.

With the rapid development of computer technology，molecular simulation method is applied to analyze the interactions between oil/water/rock systems from a smaller spatial and temporal scale.Since molecular simulation allows complete control of components and simulation conditions in the model and provides a detailed and intuitive atomic-scale simulation process，greatly facilitating the interpretation of experimental results at the microscopic level，these advantages have aroused widespread interest among researchers(Koleimi et al.，2019，2020;Kirch et al.，2018;Santos et al.，2019;Lu et al.，2009;Huai et al.，2020).A large number of fluid-solid interactions in reservoirs have been performed using molecular simulations.The results show that the fluid-solid interactions in reservoir mainly includes the following types:(1) Non-bond interactions at oil/water interface and water/rock interface，including van der Waals force，electrostatic force and hydrogen bond.Koleini et al.investigated the effect of salinity on the adsorption of crude oil components at the calcite brine interface by using molecular simulations，and showed that a diffusion double electric layer was formed near the rock surface.This non-bond interaction hindered the adsorption of polar components on the rock，while the presence of sodium ions can form electrostatic bridges that allow negatively charged polar components to adsorb on the calcite surface(Koleini et al.，2020);(2) Physical adsorption of organic polar components with specific functional groups on rock surface by electrostatic force or van der Waals force.Ataman et al.investigated the oxygen containing functional group and the polarity of the interaction between carbonate minerals by using density functional theory，and show that carboxylic acid root adsorption ability stronger than the alcohols and aldehydes in the calcite surface (Ataman et al.，2016).Their another work showed that the amino adsorption ability stronger than nitrile and sulfonium group (Ataman et al.，2016);(3) The interactions between organic macromolecules including asphaltenes and resin(Ahmadi et al.，2020).Ahmadi et al.studied the interactions between different surfactants and asphaltene components in the process of steam flooding by molecular simulation(bib_Ahmadi_and_Chen_2020Ahmadi and Chen，2020a，b)，effects of benzene and oxygen atoms in surfactants on the interaction between asphaltene，water and anionic surfactant(Ahmadi and Chen，2021) as well as the influence mechanism of anionic and non-anionic surfactant on asphaltene aggregation(Ahmadi and Chen，2020a，b).They believe that anionic surfactants can significantly reduce π-π interactions between asphaltene molecules and weaken the effect of asphaltene accumulation.Moreover，anionic surfactants are more likely to form hydrogen bonds than other surfactants due to their four adsorption sites;The location and number of heteroatoms in asphaltenes play an important role in polarity and hydrogen bond formation.Benzene rings in surfactants can significantly increase the van der Waals interaction between surfactants and asphaltenes due to the face-toface or face-to-edge staking of the benzene ring of the surfactant and aromatic rings of asphaltene molecules;The miscibility of different surfactants and asphaltenes with different aromatics is significantly different，and anionic surfactants are more likely to be miscible with asphaltenes with high aromatics than non-anionic surfactants.The oil/water/rock three-layer model is the typical molecular model for studying the wettability alteration.The influence on the wettability is studied by changing the type and content of ions in the water layer.For instance，Koleini et al.(Koleini et al.，2019) used the three-layer model to study the influence of different brine on the wettability of carbonate reservoirs.They believe that salinity is the critical factor affecting smart water flooding.Under high salinity conditions，aggregating ions will anchor the polar components of the crude oil to the mineral surface tightly.However，there are few reports on the effect of brine on sandstone wettability using the MD method.In addition，due to the special formation mechanism of tight sandstone reservoirs，the oil/water/solid three-layer model is challenging to be used to study the effect of brine on the detaching efficiency of polar components in tight sandstone reservoirs.

In this paper，naphthenic acid was selected as the polar component of crude oil.Firstly，the adsorption characteristics of naphthenic acid on quartz surfaces without water were studied by the QCM-D technique.Secondly，a molecular model of naphthenic acid adsorption was established based on experimental results.Unlike previous oil/water/solid three-layer model，the negatively charged sandstone surface after ionization of the silanol hydroxyl group is constructed in this paper.This molecular model corresponds to the process of detaching naphthenic acid by brine in the subsequent QCM-D experiment and realizing the correspondence between micro and macro.The microscopic process of detaching naphthenic acid from quartz surface by brine was reproduced more realistically.The distribution and action mechanism of NaCl，CaCl2，MgCl2，Ca2+/Mg2+and Na2SO4in the oil-water system were studied using equilibrium dynamics simulation (EMD) and nonequilibrium dynamics simulation (NEMD).Then，the detaching determination coefficient k was established based on the distribution of ions and the mechanism of action，and the results verified by QCM-D.Finally，based on molecular simulation and detaching determination coefficient，the detaching capacity under different Ca2+/Mg2+ratio was predicted and then verified by the QCM-D experiment，and the optimal Ca2+/Mg2+ratio was obtained.By combining macro QCM-D experiment and micro MD simulation，the micro mechanism of brine detaching polar components from tight sandstone was deeply revealed，which provides a theoretical basis for regulating brine properties efficiently and improving the tight sandstone recovery.

2.Material and methods

2.1.Aqueous solutions

Brine was made by dissolving NaCl(99.5%，Aladdin，China)，CaCl2(98.5%，Aladdin，China)，MgCl2-6H2O (99.5%，Aladdin，China) and Na2SO4(98.0%，Aladdin，China) in Milli-Q water.Ensuring that all aqueous solutions with cationic concentration of 0.5 M，the composition of the solutions is listed in Table 1.

Table 1Composition of aqueous solutions.

2.2.Naphthenic acid solution

Naphthenic acid solution，50 wt%，was made by dissolving Naphthenic acid (98.5%，Aladdin，China) in n-decane which represents the most abundant alkane component in crude oil.Due to the complexity of crude oil composition，it is difficult to determine the adsorption capacity of a specific component on rock surface and the subsequent detaching capacity and detaching mechanism if crude oil is directly used in experiments.In addition，our previous work have shown that the adsorption performance of non-polar alkane components on mineral surface is almost negligible compared with that of polar components(Wang et al.，2021).On the one hand，ndecane can be used as a solvent to dilute naphthenic acid，which is conducive to the smooth progress of the experiment (If pure naphthenic acid was used directly，the quartz balance would be out of range).On the other hand，it will not affect the measurement of naphthenic acid due to n-decane is almost not tightly adsorbed on the quartz surface.

2.3.Quartz crystal microbalance (QCM) measurements

The principle and experimental procedures used in the QCM-D measurements have been reported in many articles (Nourani et al.，2016;Erzuah et al.，2018;Denise，2020;Li et al.，2018)，and are summarized here for clarity.

2.3.1.Principle

Adsorption/desorption measurements were carried out using a dissipative QCM-Z500 from KSV(Helsinki，F(xiàn)inland).The dissipative quartz crystal microbalance(QCM-D)simultaneously measures the changes in resonance frequency(f)and dissipation(D).When an AC voltage is applied across the electrodes，the crystal oscillates with a characteristic frequency.The frequency changes induced upon contact with solutions can be related to mass loading，liquid loading，and liquid trapping.Liquid trapping by interfacial cavities and pores were neglected due to smooth crystal surfaces.The equations describing mass and liquid loading are given as follows:

Mass loading (Sauerbrey Equation) (Sauerbrey，1959):

When ΔD/(Δfads/3) ＞2 × 10-6，Sauerbrey equation is no longer applicable，and Voigt model is used to fit the viscoelastic film；For the loose viscoelastic film，the data analysis is very complex，and the frequency dissipation changes are caused by the different densities and viscosities of the solution.Therefore，the Kanazawa equation should be used to calculate the fluid loading effect：

Liquid loading(Kanazawa Equation)(Liu et al.，2021;Kanazawa and Gordon，1985;Liu et al.，2021):

Where，n is the harmonic number，f0is the fundamental resonant frequency(5 MHz)，Δm is the adsorbed mass，A is the active area of the crystal (0.785 × 10-4m2)，ρqis the specific density of quartz(2650 kg/m3)，vqis the shear wave velocity in quartz(3340 m/s)，C is the constant of the quartz crystal (0.177 mg/Hz m2)，ρ is the density of solution (kg/m3)，and η is the viscosity of solution (kg m-1s-1).Subscriptions 1 and 2 refer to the solutions of different densities and viscosities.

Assuming that the adsorbed layer is thin，rigid，uniformly distributed，and does not slip at the interface(Knag et al.，2004)，the mass of the adhered layer can be calculated by Eq.(1)，when a Newtonian liquid is present，the liquid also becomes coupled to the crystal oscillation，causing frequency shifts due to the change in the density and viscosity of the fluid.This frequency change can be estimated by using Eq.(2) (Sauerbrey，1959).In this paper，the naphthenic acid solution was diluted and prepared by n-decane，and the properties of the naphthenic acid solution were significantly different from that of the brine prepared by ultra-pure water.Therefore，the influence caused by the properties of the solution could not be ignored，and the following formula was used to eliminate this influence:

where，Δfadsis the frequency change caused by naphthenic acid detaching，Δfmeasuredis the frequency change caused by salt water loading，and Δfliqloadis the frequency change caused by solvation effect.

2.3.2.Procedure

Prior to each experiment，the crystals were rinsed with ethanol，followed by Milli-Q water.Next，the crystals were dried with N2and placed in 20，000 ppm aqueous sodium dodecyl sulfate solution for at least 30 min.Afterward，the crystals were rinsed with Milli-Q water and dried with air again.The QCM chamber and connections were washed with ethanol and Milli-Q water before each measurement.Before starting an experiment，the instrument was calibrated in air (Farooq et al.，2011).

Adsorption procedures:1) Base line:The chamber was flushed with pure n-decane to obtain a stable baseline.2)Naphthenic acid solution injection:2-2.5 mL of the naphthenic acid solution was injected at 60 μL/min through the temperature loop into the measurement chamber.The naphthenic acid solution was kept in the temperature loop for 5 min to equilibrate with the temperature.The naphthenic acid solution was injected twice to ensure saturation of the crystal surfaces.3) N-decane injection:2 mL n-decane was flushed over the crystal with a relatively high flow rate to rinse weakly bound naphthenic acid from the surface.

Desorption procedures：The first 2 steps are the same as the adsorption measurement.3) Different aqueous solutions are injected into the crystal to remove the naphthenic acid directly.The temperature was kept at 25°C in all experiments.The analysis is performed using the third overtone of the fundamental frequency.

3.Molecular dynamics simulation

3.1.Model construction

In this paper，molecular simulation was used to study the microdetaching mechanism of naphthenic acid on quartz surfaces by different ions.In order to reproduce the actual process more realistically，a molecular model of negatively charged sandstone surface due to ionization of silanol hydroxyl group was constructed.

The construction process of the three-dimensional periodic molecular model is as follows:(1) the quartz crystal was derived from the Materials Studio software model library and cut along the(0 0 1) surface.Then，the crystal was extended to the quartz wall surface of 3.93 nm × 4.25 nm × 1.49 nm (2) Double layer naphthenic acid model，in which the upper layer contains 15 naphthenic acid molecules，and the non-polar part faces the lower layer，while the lower layer contains 25 naphthenic acid molecules with the polar part faces the quartz wall.The length and width of the naphthenic acid layer were consistent with the quartz wall during the construction.The naphthenic acid layer is then spliced onto the quartz surface.(3)Brine layer，ensure that the length and width of the brine layer are consistent with the quartz surface and splice with the previous model.(4) Ionize the hydroxyl groups of silicols on the quartz wall and add the corresponding amount of sodium ion balance charge in the third step in advance.(5) Under the naphthenic acid layer and quartz wall fixation condition，the model was optimized by 50，000 steps.Then the dynamics simulation was carried out by 10ps (the simulation here is only for mixing the water layer with the naphthenic acid layer and is not included in the formal simulation calculation).The initial configuration is shown in Fig.1 then the data file is exported.The formal simulation is carried out using LAMMPS，a large-scale atomic/molecular parallel simulator of Sandia National Laboratory (Plimpton，1995).

3.2.Potential models

The accuracy of molecular dynamics rests on the correct choice of mathematical equations and interaction parameters for the potential energy.In this study，water molecules were described by the TIP4P model(Abascal et al.，2005).Many classical force fields，such as LFF (Lopes et al.，2006;Chilukoti et al.，2014) CLAYFF (Skelton et al.2011，2011，2011;Bourg and Steefel，2012;Cygan et al.，2004)，CWCA (Cruz-Chu et al.，2006)，and LR (Lee and Rossky.1994)，are often used to describe the interaction between water and silica.In these models，the CLAYFF force field is in good agreement with first-principles simulation results because it is developed for hydrates and multi-component minerals and their interactions with fluid interfaces(Skelton et al.2011，2011，2011).In addition，the Clayff force field can reproduce the orientation characteristics of organic matter adsorption on the silica surface(Ledyastuti et al.，2012) and predict the adsorption isotherm of propane in the quartz pore accurately (Le et al.，2015).Therefore，the potential energy model will be used here to describe the quartz.Na+/Cl-，Ca2+/Mg2+，andare described by the same field parameters as used by Smith(Smith and Dang.1994)，Aqvist (Aqvist，1990)，and Williams (Williams et al.，2014).The OPLS-AA model(Jorgensen et al.，1996) of the all-atomic force field is adopted to characterize naphthenic acid molecules，and the LJ potential energy between different atoms (with a cut off 1.2 nm) was calculated using the geometric mean mixing criterion.For over 20 years，the OPLS-A A force field has been efficiently used in molecular modelling studies of proteins，carbohydrates，nucleic acids，and applied in computer modelling of many organic compounds，including decane and shorter alkanes (Bratek et al.，2020).The outstanding feature of the OPLS-AA force field is that the potential energy parameters between different atoms are obtained by fitting the laboratory experiment results，which is highly reliable.This force field model has been widely used to reproduce the thermodynamic properties and configuration characteristics of alkanes，polymers，and even biological macromolecules (Wang et al.，2015;Soetens et al.，2001;Michalec and Lísa，2017;Cerar et al.，2017;Gontrani et al.，2009)，and the research results have been widely recognized.The general form of its mathematical model is as follows (Jorgensen et al.，1996):

Where，E is potential energy，kcal/mol;Kris the harmonic weight of bond stretching，kcal/(mol·2);r is the distance between two bonding atoms，?;r0is the equilibrium bond length，?;Kθ is the harmonic weight of bond angle bending，kcal/(mol·radian2);θ is the bond angle，radian;θ0is the equilibrium bond angle，radian;α is dihedral angle，radian;C1，C2，C3，and C4are the coefficients of the Fourier series，kcal/mol;α0is the phase angle，radian;rijis the distance between atoms i and j，?;σijis the collision radius，?;εijis the potential well depth，kcal/mol;qiand qjare the charges carried by atoms i and j;rcutis the truncation radius，?.

Fig.1.Initial configuration.

Eqs.(6)-(9) correspond to bond stretching potential energy，bond bending potential energy，dihedral angular potential energy，and non-bonding potential energy (including Lennard-Jones potential energy and electrostatic potential energy)，respectively.Eq.(9) can also be used to describe the interatomic interaction force with more than 3 chemical bonds spaced within molecules，where fijbetween 1 and 4 atoms is 0.5，and fijis 1 in all other cases(Jorgensen et al.，1996).The truncation radius of van der Waals force is 1.20 nm.The geometric mean mixing criterion calculates the non-bonding potential energy between different atoms.Meanwhile，to improve the simulation efficiency，the Particle-Particle-Particle -Mesh (PPPM) (Luty et al.，1994) algorithm calculates the long-range static force.

3.3.Simulation details

The Large-scale Atomic/Molecular Massively Parallel Simulator(LAMMPS)，distributed by Sandia National Laboratories，was utilized to perform MD simulations (Plimpton，1995).First，the potential energy of the initial configuration(Fig.1)was minimized by iterative adjustment of these atomic coordinates using the conjugate gradient algorithm with the quartz and naphthenic acid in the initial configuration fixed.The system was then relaxed for 500 ps under NVT ensemble (constant number of atoms，isovolumetric，and isothermic conditions)，and the temperature of the system is controlled by Berendsen algorithm at 300 K.After that，keep the quartz wall fixed and cancel the fixation of naphthenic acid，the Nosé-Hoover algorithm was used to control the system's temperature (coupling constant:0.1 ps)，2000ps simulation was carried out under the NVT ensemble，and data were collected at 1ps intervals in the last 500ps for statistical analysis(Wang，2016).When the system reached equilibrium，acceleration along the x-direction was applied to the brine，1000ps NEMD simulation was carried out under the NVT ensemble to simulate the detaching process.The visualization results in this paper are realized by Ovito Basic 3.4.4 software (Stukowski，2010).

4.Results and discussions

4.1.Adsorption characteristics of naphthenic acid based on QCM-D

Fig.2 shows the frequency and dissipation shifts during the naphthenic acid adsorption on the quartz surface.The resonant frequency of the quartz immediately decreased when the naphthenic acid solution was injected into the quartz surface，indicating that naphthenic acid molecules were quickly adsorbed on the quartz surface.The resonant frequency did not change until 20 min，and the adsorption amount on the quartz surface also reached saturation.After that，n-decane was injected at a high flow rate to rinses the weakly bound components，and the resonance frequency increased but did not return to the initial value，indicating that some naphthenic acid molecules got adsorbed on the quartz surface irreversibly.Moreover，the low dissipation value after the ndecane solution rinses the weakly bound component meanings that naphthenic acid forms a rigid film on the quartz surface.Since there is no aqueous solution involved in the whole adsorption process，the quartz surface is electrically neutral and covered by silanol groups，and naphthenic acids is adsorbed on the silanol site through polar interaction.

Fig.2.Frequency and dissipation shifts caused by naphthenic acid adsorption on quartz

Since the viscosity of the naphthenic acid solution is higher than that of the n-decane solvent，when calculating the adsorption amount of naphthenic acid the frequency shifts caused by density and viscosity cannot be ignored.The viscosity and density values of n-decane and naphthenic acid solutions at the experimental temperature are shown in Table 2，and the frequency shifts caused by the difference in density and viscosity can be calculated by Eq.(2).The total mass density of naphthenic acid adsorbed on the quartz was calculated by Eq.(1).

Table 3 records the adsorption mass of naphthenic acid on quartz surface measured by repeated experiments.Experiments 1-8 are all the systems with n-decane rinsing the weakly bound components (as shown in Fig.2).This process aims to distinguish the solid and weak bound adsorption layers of naphthenic acid on quartz.It can be seen from Table 3 that there are solid and weak binding adsorption layers for naphthenic acid adsorption on quartz.The repeated experiments found that the adsorption mass density of solid binding layer was 106.20-141.60 ng/cm2，with an average of about 128.33 ng/cm2.This is because the solid binding adsorption layer mainly depends on the strong polar interaction between naphthenic acid and quartz surface，and the number of adsorption sites on quartz surface is constant，representing a slight fluctuation on adsorption mass.On the contrary，the weak binding component mainly depends on its weak non-polar interaction with the solid binding layer，and the adsorption mass has certain randomness，resulting in a wide range of total adsorption mass in experiments 1-8.We conducted another 10 groups of experiments(Experiment 9-18)to clarify on the total adsorption mass and took the average adsorption mass 210.98 ng/cm2as the total adsorption mass for subsequent research.

Table 2Density and viscosity for solutions.

Table 3Mass density of naphthenic acid adsorption on quartz.

Where，Δfmeasuredis frequency shifts caused by the naphthenic acid solution injection，measured value;Δfliqloadis solvent loading effect，calculated value;Δfadsis Δfmeasured-Δfliqload，representing the frequency shifts caused by naphthenic acid adsorption;Δmadsis naphthenic acid adsorption mass density;Δfsis frequency shifts caused by naphthenic acid solid binding layer，measured value;Δmsis mass density of naphthenic acid solid binding layer.

4.2.Molecular simulation

The molecular model was constructed based on the conclusion of the adsorption experiment，as shown in Fig.1.The mass rationing of naphthenic acid into solid adsorption and weak adsorption is 5:3 in this model，consistent with the experimental results.After 2.5ns EMD and 1.0ns NEMD simulation，the final configurations are shown in Fig.3，with the water molecules hidden for a clear view.Due to the complex microscopic interactions，the distribution of ions with different characteristics can be found.The mechanismbehind this distribution on the detaching of naphthenic acid were analyzed from static，dynamic，and interaction aspects.

Fig.3a shows，there are two states for naphthenic acid molecules in the system generally after EMD.A continuous layered structure distribution present in H2O，MgCl2，and 50%CaCl2systems more obviously.However，the naphthenic acid molecules are aggregated and distributed in clusters，which are discontinuous in the x-direction.No obvious layered structures are pronounced in zdirection in NaCl，CaCl2，and Na2SO4systems.As shown in Fig.3b，after NEMD，naphthenic acid molecules in all systems moved away from the quartz，and the distribution became looser in the z-direction.However，there was still an apparent cluster behavior in CaCl2and Na2SO4systems.

4.2.1.Distribution of ions and carboxylate

Fig.4 shows the number density distribution of different ions and carboxylic acid in the system after the equilibrium (left column) and non-equilibrium (right column) systems reached a steady state.By analyzing the distribution positions of ions and carboxylate under equilibrium dynamics conditions，the mechanism of distribution of ions in the oil-water system can be deeply understood.Compared with the equilibrium state，the dynamic change of ion distribution in detaching naphthenic acid can be reflected in the non-equilibrium dynamics，which reveals the law of movement of ions and the detaching mechanism of naphthenic acid.

After EMD，the negatively charged quartz surface is more inclined to adsorb sodium ions by electrostatic interaction，and Na+presents a very high number density distribution in the vicinity of the quartz surface in all systems as shown in Fig.4.In the system with NaCl (Fig.4b)，besides adsorbing on the quartz surface，other Na+coupled with Cl-through electrostatic interaction，forming another electric layer about 3 nm away from the quartz surface.The electric layer is thick and overlaps with the weakly bound naphthenic acid layer in the z-direction.In the system with CaCl2(Fig.4c)，part of Ca2+also adsorbed on the quartz surface besides a large amount of Na+.Since the system with CaCl2has twice as many Cl-at the same cation concentration as the system with NaCl.The Cl-generally stays away from the quartz surface under electrostatic interaction and competes with the quartz surface to adsorb cations.Therefore，some Na+can be adsorbed by Cl-resulting in quartz surface charge imbalance，and some Ca2+moving to the quartz surface to balance the charge.Similarly，the remaining Ca2+and Clformed an electric layer about 3.5 nm away from the quartz surface，which is thinner than that in the NaCl system.The distribution of ions is relatively concentrated，and there is almost no overlap with the weak binding layer of naphthenic acid;There is almost no electric double layer in the system with MgCl2(Fig.4 d).Except for the higher Na+density on the quartz surface，the other Na+and Mg2+distributed dispersedly in the z-direction.Due to the smaller molecular weight，Mg2+moves faster than Ca2+in the same charge condition，so the competitive adsorption between Mg2+and Na+on the quartz surface becomes intense.Because of competition between Mg2+and Na+，Cl-moves towards the quartz surface and finally reached equilibrium under the electrostatic interaction.As a result，the ions are relatively dispersed in the z-direction，and the cation density decreases with the increased distance from the quartz surface.The distribution region of ions almost covers the naphthenic acid region completely;For the system with 50%CaCl2(Fig.4 e)，since the concentration of Cl-remains constant and the concentration of Mg2+is halved，Cl-is less likely to move towards the quartz surface.As shown in Fig.4e，Cl-concentrates about 3.5 nm away from the quartz surface and absorbs a certain amount of Na+by electrostatic interaction，which creates conditions for Mg2+to be adsorbed on the quartz surface.Therefore，the numberdensity of Mg2+ions near the quartz surface is higher in the 50%CaCl2system.Due to the competition between cations on the quartz surface，some ions overlap with naphthenic acid between the two electrical layers;Theand Na+formed two electrical layers at 1 nm and 3 nm from the quartz surface in the system with Na2SO4(Fig.4 f).Since Na+ions of the Na2SO4system is twice that of the NaCl system at the same concentration of anions.However，the quartz surface tends to adsorb Na+，causing part ofto be carried close to the quartz surface.In spite of this，mostleaves the quartz surface by electrostatic repulsion.

Fig.3.Final configurations a:Final configuration for EMD b:Steady configuration for NEMD From left to right:H2O，NaCl，CaCl2，MgCl2，50% CaCl2，Na2SO4.

In the process of detaching naphthenic acid (i.e.，NEMD)，the distribution mechanism of the ions and carboxylic acid groups are consistent with the EMD on the whole.However，there are apparent changes in the system with MgCl2and 50%CaCl2systems.During the detaching process，the Mg2+in both systems tend to move closer towards the quartz surface.

4.2.2.Mechanism of detaching naphthenic acid by ions

Fig.4.Carboxylate and ions distribution EMD (left column) and NEMD (right column) a:H2O，b:NaCl，c:CaCl2，d:MgCl2，e:50% CaCl2，f:Na2SO4.

The visualization results show that ions and water molecules can coordinate to form hydrated ions，and the formation of hydrated ions directly affects the detaching of naphthenic acid.We characterized the microstructure of different systems by analyzing the RDF (Radial Distribution Function) among inorganic ions，carboxylic acids，and water molecules.By combining RDF，ions and carboxylate distribution characteristics，the role and mechanism of inorganic ions in the detaching naphthenic acid were revealed.RDF refers to the possibility of finding a target particle within a certain distance from the central particle to reflect the binding strength of the particle and the target particle(Hsissou et al.，2020).RDF plots provide essential information about macromolecular structures around specific atoms.Combining a RDF plot and visualization of molecular trajectories helps to grasp interactions between desired molecules (Ahmadi and Chen.2021).

Where，g(r) refers to the prob capacity of finding y atom at r away from x atom;δr is the thickness of spherical shell;n(r) is the number of particles in spherical shell，which is used to reflect the microstructure between fluids.

As both naphthenic acid and water molecules are polar molecules，the carboxyl group on naphthenic acid is likely to have polar interactions with water molecules.RDF between carboxyl group H of naphthenic acid and O of water molecule was calculated as shown in Fig.5a.Due to the high electronegativity of the oxygen，the hydrogen of the hydroxyl group (-OH) can make strong hydrogen bond with the oxygen of water molecules.It is also worth mentioning that the double-bonded oxygen in the carboxyl(-COOH)can also make a hydrogen bond;however，this hydrogen bond is much weaker than that taking place with hydroxyl(Ahmadi and Chen.2021).Hydrogen bonds bind water molecules to naphthenic acid molecules tightly and play a significant role in the detaching process.

RDF between ions and water molecules is shown in Fig.5b.For a pure water system，although hydrogen bonds can be formed between H and O atoms in water molecules，while the low peak value of RDF indicates that the interactions between water molecules themselves is lower than that between ions and water molecules，which means it is difficult to detach naphthenic acid continuously and effectively;For cations，the water coordination number for Ca2+and Mg2+is greater than Na+，while the hydrated radius for Mg2+is smaller than that of Ca2+and Na+.This means that the hydrogen bond between Mg2+and water molecules is strongest among the three cations.According to the distribution characteristics of ions in Section 4.2.1，the number density of cations in the system with MgCl2decreased with the increase of the distance from the quartz surface，so it was not conducive to detaching naphthenic acid;Ca2+shares similar water coordination number with Mg2+.However，the interaction distance between Ca2+and water molecules is larger than that between Mg2+and water molecules.Even though some Ca2+are closely adsorbed on the quartz surface，the “anchoring” effect of naphthenic acid is weak with Mg2+.In addition，most of the Ca2+are distributed far away from the quartz surface.It can promote the detaching of naphthenic acid due to the above distribution feature combined with high water molecular coordination number.The radius of hydrated Na+is roughly the same as Ca2+，that means Na+and Ca2+interact with water molecules at nearly the same distance，but the peak value of RDF between Na+and water molecule is lower than that between Ca2+and water molecule，which means that the interaction between Na+and water molecule is the weakest among the three cations.Despite，in the sodium chloride system，there are no other Na+on the quartz surface that can“anchor”naphthenic acid，other than the sodium ions that balance the charge.Therefore，the system with NaCl can promote the detaching of naphthenic acid.The peak value of RDF betweenand water molecule is relatively low，which means that the hydrogen bond betweenand water molecule is weak.However，there are more Na+in the system with Na2SO4.Despite a certain number ofmoving around the quartz surface，mostexist far from the quartz surface，and Na+coordinated with this part ofhave strong interaction with water molecules，and have strong promotion effect on the detaching of naphthenic acid;For the 50% CaCl2system，Mg2+is more likely to be adsorbed near the quartz surface，but the adsorption strength is lower than that of Ca2+(At the quartz surface，the adsorption position of Ca2+z=1.75 nm，and that of Mg2+z=2 nm).In addition，the remaining Mg2+and Ca2+were distributed far from the quartz surface to promote the detaching of naphthenic acid molecules synergistically.The detaching mechanism of naphthenic acid by water molecules and ions is shown in Fig.6.

The above analysis identified the position of ions in aqueous solution and the mechanism of detaching naphthenic acid from a qualitative perspective of interaction.In order to quantitatively explore the detaching capacity of ions to naphthenic acid，we coupled the ion position distribution with a detaching mechanism to analyze the detaching capacity of ions to naphthenic acid under same concentration condition.By parameterizing the distribution of ions positions and detaching mechanism，the detaching determination coefficient k is proposed.

4.2.3.Detaching determination coefficient

According to analysis，the influence rule of distribution position on detaching is as follows:For the same kind of cations，the closer to the quartz surface，the less likely to detach naphthenic acid.In Fig.4，the area enclosed by ions number density curves represents the number of ions.The x-axis is divided into three parts from the quartz surface to where the ion density equals zero.The anchoring region(1.5nm-3.0 nm)，the transition region(3.0nm-4.0 nm)，and the detaching region(4.0nm-6.0 nm)were named from bottom to top，and the areas of each part were S1，S2，and S3，respectively，as shown in Fig.7.

Fig.5.Radial distribution function among particles pair a:RDF between naphthenic acid and water b:RDF between ions and water.

Fig.6.Diagram of detaching mechanism.

The interaction between ions in the anchoring region and quartz is very sensitive to distance due to the anchoring region being close to the quartz surface.Even for the same ions，the closer to the quartz surface，the stronger the anchoring strength.Therefore，instead of simply corresponding the integral curve area to the number of ions，the corresponding adjustment should be made according to the curve peak position in the anchoring region.In order to simplify the calculation，if the curve peak position appears to the left of the middle line of the anchoring region，it means that the average binding strength of cations to quartz surface in the anchoring area is higher than that of cations to quartz wall in the midline position，and the integral area S1is multiplied by 1.2 as the effective area S1e.If the curve peak position appears to the right of the middle line of the anchoring region，it means that the average binding strength of cations and quartz surface in the anchoring area is lower than that of cations and quartz wall in the midline position，the integral area S1is multiplied by 0.8 as the effective area S1e，and the calculation formula is as follows:

Where，S1eis the effective area of the anchoring area;S1is the actual integral area of the anchoring region;rmaxis the peak position of the anchoring region;rcis the middle line of the anchoring area，rc=2.25 here.

Fig.7.Diagram of regional division (taking the system with MgCl2 as an example).

Since the detaching region is far from the quartz surface，the same ions located in the detaching region can be considered to have the same detaching capacity.Therefore，the area of the detaching region can be directly regarded as the number of ions.In addition，with the increased distance between the ions in the transition region and the quartz surface，the anchoring capacity decreases，and the detaching capacity increases.Therefore，according to the ratio of the detaching region normalized area to the anchor region normalized areafor all cationic，the transition region is divided into detaching transition region and anchoring transition region as following rules:

Where，S1eis the normalized effective area of the anchoring region;S2is the normalized area of the transition region;S3is the normalized area of the detaching region;S21is the normalized area of the anchoring transition region;S23is the normalized area of the detaching transition region;0.75 is the correction coefficient，which is used to correct the deviation by human caused in the division of regions.After the equilibrium of each system，the normalized area occupied by different ions in the three regions was calculated by integral，and the calculated results are shown in Table 4.The sum of the normalized effective area of the anchoring region and the anchoring transition region area was defined as the total anchoring effective area，and the sum of the normalized area of the detaching region and the detaching transition region area was defined as the total detaching effective area for subsequent calculation.

Firstly，the interaction between different ions in the anchoring region and quartz surface is calculated to determine the adsorption capacity of quartz surface to the ions in the anchoring region.The intensity of interaction between ions and the quartz surface isdetermined by the non-bonding interaction between ions and oxygen atoms on the quartz surface.The middle position of the anchoring region is about 0.75 nm away from the quartz surface，and the interaction intensity between ions and the quartz surface at this location is defined as the average interaction intensity between ions in the anchoring region and the quartz surface.The interaction parameters of ions and oxygen atoms are shown in Table 5.The interaction between ions and oxygen atoms adopts geometric mixing rules.

Table 4Normalized area of each region in different systems.

Secondly，the strength of interaction between different particles and water molecules is calculated to determine the binding strength of various particles to water molecules in the process of detaching and anchoring.Hydrated ions formed by different ions with different radius and water molecular coordination numbers，which is reflected in the peak position and the peak value of the g(r)，respectively.The hydration strength is inversely proportional to the radius and positively proportional to the coordination number.So，the ratio of the peak value to the x-coordinate is defined as the hydration strength h.The normalized hydration strengthh for different particles are shown in Table 6.

Table 6Normalized hydration strength.

The product of the detaching region normalized area (include detaching region and detaching transition region) of each particle in the system and the corresponding hydration strength is defined as the detaching strength k1:

Similarly，the anchoring strength k2is as follows:

k1-k2is the detaching determination coefficient k：

Where，h is normalized hydration strength;E is the normalized interaction intensity between particles and quartz surface;c is the binding coefficient;Define k for pure water system as 0，and c=1.637 was obtained.

The detaching determination coefficient k was established by considering the number of ions in three regions，types，hydrationstrength，and the distance between the quartz surface and ions in the anchoring region.When k ＞0，it is beneficial to the detaching of naphthenic acid，and the higher k is，the stronger the detaching capacity is，vice versa.The above formula is used to calculate the k for each system，as shown in Fig.8.

Compared with the pure water system，all the other systems promoted the naphthenic acid detaching except MgCl2system，and the detaching capacity in descending order was as follows:50%CaCl2＞NaCl ＞Na2SO4＞CaCl2＞H2O ＞MgCl2.

4.3.Detaching efficiency based on QCM-D

The efficiency of detaching naphthenic acid with different brine was studied based on the adsorption experiments 9-14 in Table 3.Since the viscosity and density of the solid binding component remaining on the quartz surface cannot be determined after the weak binding component is rinsed by n-decane，the frequency shifts caused by the viscosity and density cannot be calculated if the brine is injected at this time.For this reason，after naphthenic acid adsorption equilibrium is followed by injection of brine，frequency loading can be calculated by using the viscosity and density difference between naphthenic acid and brine.In addition，according to the analysis in Section 4.1，the weakly bound components in each experiment differ significantly.In order to facilitate the comparison of the desorption efficiency of brine，the average value of 18 groups of data in Table 2 was taken as the standard sample，and Eq.(21)was used to calculate the desorption efficiency.The diagram of detaching efficiency calculation is shown in Fig.9.

Table 7Detaching efficiency for different brine.

Fig.8.Detaching determination coefficient for different brine.

Table 5Normalized interaction intensity between particles and quartz surface.

Fig.9.The diagram of detaching efficiency calculation.

Where，R is detaching efficiency，%;Δm is the total mass of naphthenic acid adsorbed on the quartz，ng/cm2;Δm1is the mass change caused by saline rinsing，ng/cm2;Δm2is the mass change caused by ions adsorption on the quartz surface in the blank saline experiment，ng/cm2;Δmadsis the average total mass of naphthenic acid adsorption on the quartz surface，ng/cm2.

The frequency shifts during the brine detaching naphthenic acid from the quartz surface are shown in Fig.10a.It can be seen that the resonance frequency recovered quickly with the injection of brine.On the one hand，there is a significant viscosity and density difference between brine and naphthenic acid solution.On the other hand，the material exchange happens when brine is in contact with the quartz surface.Part of the naphthenic acid molecules are detached resulting in an increase in resonance frequency，while ions adsorption resulting in the decrease of resonance frequency at the same time.The viscosity and density of brine are treated as pure water approximately because of the low concentration of brine.The frequency shifts caused by the difference in viscosity and density between oil and water are calculated by Eq.(2).In addition，the influence of ion adsorption on the quartz surface was excluded through blank experiments，as shown in Fig.10b.

Fig.10.Detaching efficiency for different brine based on the QCM-D a:The frequency shifts during the brine detaching naphthenic acid b:The frequency shifts for blank experiments.

Fig.11.Detaching efficiency for different brine.

The efficiency of naphthenic acid detaching is shown in Table 7 and Fig.11.It can be seen from Fig.11 that the efficiency of the other systems in the naphthenic acid desorption increased compared with the H2O system，excluding the MgCl2system.The detaching efficiency of the 50% CaCl2system was the highest，which was about 3 times that of the pure water system.Using core displacement experiments，Rezaeidoust et al.，demonstrated that Mg2+can replace Ca2+from the rock surface，and Ca2+can carry the adsorbed acidic components away from the rock surface to improve the acidic components recovery (Rezaeidoust et al.，2009).This mechanism is in essence consistent with the mechanism in this paper that Ca2+and Mg2+combination can improve the detaching efficiency.The detaching efficiency of the other systems in the order from high to low was as follows:NaCl ＞Na2SO4＞CaCl2.Likewise，traditional laboratory experiments have also reached a similar conclusion.Using core displacement experiment，Huang demonstrated that the effect of Ca2+to improve oil recovery for sandstone reservoir is better than that of Ca2+(Huang，2019).By measuring contact angle，interfacial tension，and core displacement，Chai et al.concluded that reducing Ca2+concentration and increasingconcentration had a good effect on enhancing oil recovery (Chai et al.，2021).Through laboratory experiments，Lager et al.proved that divalent cations are more capable of bridging and absorbing acid components in crude oil than monovalent cations.Therefore，Na+have a better capability of naphthenic acid desorption than Ca2+and Mg2+(Lager et al.，2008)，which is consistent with the conclusion of this paper.which verifies the correctness of the detaching determination coefficient in Section 4.2.3.

Table 8Normalized area of each region in different systems.

4.4.Prediction and application

The detaching efficiency of naphthenic acid by brine mixed with CaCl2and MgCl2in different ratio was studied by equilibrium dynamics simulation.After the simulation reached equilibrium，the number density distribution of ions was shown in Fig.12，and the areas of ions in different regions were calculated by integration，as shown in Table 8.By substituting the area data into Eq.(20)，the relationship between k and CaCl2percentage is shown in Fig.14.As can be seen，with the increase of CaCl2percentage，k increases firstly and then decrease.When the CaCl2percentage increases to 70%，k reaches the maximum.After that，k decreased with the CaCl2percentage increase，and the detaching capacity decreased.The detaching experiment was carried out using QCM-D，as shown in Table 9 and Fig.13，and the results are in good agreement with the predicted results.When the CaCl2percentage increases to 70%，the detaching efficiency reaches the maximum，which further verifies the reliability of this model.

Table 9Detaching efficiency for brine with different Ca2+/Mg2+ ratios.

Fig.12.Carboxylate and ions distribution of the systems with change CaCl2 percentage a:10% b:30% c:70% d:90%.

Fig.13.Detaching efficiency for brine with different Ca2+/Mg2+ratios based on the QCM-D a:The frequency shifts during oil film stripping with different calcium chloride content b:The frequency shifts for blank experiments.

Fig.14.Variation of detaching efficiency and k with CaCl2 percentage.

According to Section 4.2，when MgCl2is 100%，the cation distribution in the z-axis direction is relatively dispersed due to the competitive adsorption of Na+and Mg2+，and the closer to the quartz surface，the higher the cation concentration is，which is not conducive to the detaching naphthenic acid.With the decrease of Mg2+，the competition between the divalent cation and Na+on the quartz surface weakened，leading to the gradual departure of the cation from the quartz surface and the weaker ability of“anchoring” naphthenic acid.When the CaCl2percentage exceeds 70%，part of Ca2+begins to move to the quartz surface to balance the charge.Since the adsorption strength of Ca2+on the quartz surface is greater than that of Mg2+，the “anchoring” capacity will be strengthened to reduce the detaching efficiency.

In this paper，the detaching mechanism of small polar compounds in crude oil by ions and the method of determination detaching ability are proposed.However，there is still a lack of understanding of the detaching mechanism and detaching ability determination method of macromolecular polar compounds such as asphaltenes and resins which are widely existing in crude oil.The detaching of these compounds will be studied systematically and in detail in our coming work to better guide the formulation of injection water for enhanced oil recovery.

5.Conclusions

In this paper，the detaching behavior of petroleum acid from sandstone surface assisted by inorganic salt ions was studied by combining the MD and QCM-D techniques.MD method was used to analyze distribution of inorganic salt ions and their action mechanism，and an effective method to determine their detaching capacity was proposed.The prediction results of the model are verified by QCM-D technology，and the main conclusions are as follows:

1.The capacity of inorganic salt ions to detach naphthenic acid from sandstone surface is related to its distribution and hydration capacity.

2.Under the condition of single inorganic salt，Na+ and Cl- can form an electric layer far from the quartz surface.Part of Ca2+can be closely adsorbed on the quartz surface，while the rest Ca2+ and Cl- form an electric layer far from the quartz surface.Mg2+is distributed along the z-axis，and the closer to the quartz surface，the more it accumulates;Part of Na+carries close to the quartz surface，while most of Na+ and form an electric layer far away from the quartz surface;Under the condition of equal proportion，Mg2+ can replace Ca2+ and adsorb near quartz，and Ca2+/Mg2+/Cl- can form an electrical layer far away from the quartz surface.

3.Hydrogen bond between naphthenic acid and water molecules is key to detaching naphthenic acid;Cations enhance detaching of naphthenic acid by forming hydrated ions with water molecules.The detaching capacity is negatively correlated to the radius of hydrated ions and positively correlated to the coordination number of water molecules.

4.The detaching determination coefficient k was established by considering ion distribution，ion type，and hydration strength.When k ＞0，it is beneficial to the detaching of naphthenic acid，and the higher k is，the stronger the detaching capacity is，and vice versa.

5.Under a condition of same cation concentration，the detaching efficiency was as follows: 50%CaCl2 ＞NaCl ＞Na2SO4 ＞CaCl2 ＞H2O ＞MgCl2;In the Ca2+/Mg2+system，the detaching efficiency first increases and then decreases with the increase of CaCl2 percentage，and the detaching efficiency reaches the maximum when the CaCl2 percentage is 70%.

Acknowledgements

The authors appreciate the financial supports by the National Natural Science Foundation of China (No.52074316)，Science Foundation of China University of Petroleum，Beijing (No.2462018QNXZ01)，Open Fund (No.SXCU-201905) of Shaanxi Cooperative Innovation Center of Unconventional Oil and Gas Exploration and Development(Xi’an Shiyou University)，and Major Science and Technology Project of China National Petroleum Corporation (No.2019E-2608).