Wormholing Influenced by Injection Temperature in Carbonate Rocks

Heng Xue^1,2,3*, Zuxi Huang¹, Liqiang Zhao², Hehua Wang^1,3, Pingli Liu¹

¹China ZhenHua Oil Co., Ltd., Beijing, China.
²State Key Laboratory of Oil & Gas Reservoir Geology and Exploitation Engineering (Southwest Petroleum University), Chengdu, China.
³Chengdu North Petroleum Exploration and Development Technology Co. Ltd., Chengdu, China.

DOI: 10.4236/ojogas.2019.41002                

Abstract

This paper presented empirical models of describing reaction rate vs. hydrochloric acid temperature and concentration by regressing experimental data. And this paper also introduced the dependent reaction heat model into the thermal non-equilibrium models and coupled with two-scale continuum model to obtain governing equations for describing wormholing under non-isothermal conditions. The governing equations were discretized by implicit difference method and solved by programing. The effects of temperature on wormholing have been investigated based on the simulation results on 2-D vision. A significant difference of the effluent temperature between the dependent reaction heat model in present work and the constant reaction heat model in available literatures was observed, especially in high injection rate and strong acid concentration. In addition, the tendencies of optimum injection capacity vs. injection temperature under isothermal and non-isothermal conditions were almost different. Finally, an optimum injection temperature was found by changing the injection temperature under non-isothermal conditions.

Keywords

Carbonate Reservoir, Matrix Acidizing, Temperature, Wormhole Propagation, Numerical Simulation

1. Introduction

Matrix acidizing is a most important method to stimulate carbonate reservoir through bypassing the damage zone with wormholes using a small volume of acid. Therefore, nearly all of the studies focusing on recovering or stimulating the production of carbonate reservoirs were based on investigations of the initiation and propagation of wormholes or unstable dissolution fronts. Many mathematical models were developed to describe the wormhole patterns observed in experiment results [1] [2] [3] [4] : the capillary tube model [5] [6] [7] ; the network model [8] ; the averaged model [4] [9] ; the fractal model [1] [2] . Panga et al. (2005) coupled Darcy-scale dissolution model with the pore-scale model to summarize the governing equations which is called the two-scale continuum model [10] . In this model, the macroscopic properties such as mass transfer coefficients and dispersion properties were properly considered. Since the results calculated from the two-scale continuum model matched experiment results very well, most of the extension simulation studies were based on them. Kalia and Balakotaiah (2007) studied dissolution patterns in radial flow by extending the linear model to radial model [11] . Maheshwari et al. (2013) analyzed the 3-D dissolution patterns [12] . Smith et al. (2013) simulated core CO₂ flooding of low permeability carbonates coupled with X-ray computed microtomography (XCMT) [13] .

The wormholing process is very complicated since there are many known factors influencing convection, dispersion and reaction properties, such as acid system, concentration of acid, porosity or permeability distribution, fracture or vug medium, injection conditions, acid flow types and temperature. In the physical heterogeneities, Kalia and Balakotaiah (2009) used a unique heterogeneity parameter to characterize porosity heterogeneity, and found that the dissolution depended on the orientation of fractures [14] ; Liu et al. (2012) introduced a normal distribution instead of a uniform distribution of porosity to study the wormholing in radial flow, and they also claimed that only a large number of continuous vugs could affect wormholing [15] ; Izgec et al. (2010) further addressed the effects of vugs on matrix acidizing based on experiments [16] ; Cohen et al. (2008) motivated a thorough study of geometry effects on wormholing from pore scale to wellbore scale in 3-D vision [17] . And in the acid system, most of the studies focus on diverting acid because of extensive applying in both sandstone and carbonate heterogeneous reservoirs. The simulation results of diverting effects achieved a high agreement with experimental results [18] [19] [20] [21] [22] . The acid types (HAc, CDTA, EDTA and HCl acid system) affecting wormhole patterns, pressure drops, PVBT and wormhole density also achieved an extensively discussion through both laboratory experiments and numerical simulations [15] [18] [22] [23] .

In addition to the work mentioned above, almost all of them ignored the effects of temperature on wormholing in matrix acidizing simulation. However, temperature is also an important factor influencing wormhole patterns formation. Kalia et al. (2009) have introduced the temperature model into two-scale continuum model and compared the wormhole patterns under isothermal and non-isothermal conditions by changing the viscosity of acid according to temperature under liner flow [24] . Based on the work of Kalia et al., Bousri et al. (2012) have studied the similar problem with local non-equilibrium conditions [25] and the adiabatic boundary and isothermal boundary in wormholing were discussed under radial flow [23] [26] . It demonstrated that the temperature also played a leading role in governing fluid flow as dispersion, convection and reaction did. The reaction heat is a constant as treated in the previous work while it is confirmed that it will be affected by temperature and pressure in matrix and fracture in carbonate reservoirs [27] . In order to simulate temperature profile of core more accurately, the dependent reaction heat model was introduced into two-scale continuum model to simulate wormhole initiation and propagation under different injection conditions. The simulation results calculated from both reaction heat models were compared and discussed. To understand effects of injection temperatures on wormholing, we have studied acid volume of breakthrough, optimum injection capacity and overall pressure drops under different acid injection temperatures under isothermal and non-isothermal conditions using this simulation method.

2. Mathematical Models

2.1. Darcy Flow Model

The fluid flow in porous medium conforms to the equivalent Darcy’s law given by Ratnakar et al. (2013) for linear flow [22] :

(1)

where $\Vert .\Vert$ represents the norm of a vector/matrix. U is the velocity (m/s). P is the formation pressure (MPa). K is the permeability of the reservoir (10⁻³ × μm²). μ_eff is the effective viscosity (mPa・s).

Equation (2) shows the effective viscosity model of acid fluid:

(2)

where $\varphi$ is the porosity (Dimensionless). n is the power index. K_ten is the permeability tensor. μ₀ is the apparent viscosity (mPa・s). It should be noted that Darcy law is adapted to the flow of Newtonian acid when the parameter, n, in the effective viscosity term equals to 1. While, it will be adapted to non-Newtonian acid flow when n is less than 1. In this paper, we assume that HCl acid system is the Newtonian fluid. And the Mass balance or continuity model of acid phase is given by:

(3)

2.2. Darcy-Scale Model

In order to keep track of the concentration of hydrogen ion, the model to describe the hydrogen ion using species balance is shown as follows [10] :

(4)

where D_e is the effective dispersion coefficient in the acid phase (m²/s). $a_{\text{v}}$ is surface area per unit volume in the core available for reaction (m²/m³). C_f is the cup-mixing mass concentration of the acid in the fluid phase (kmol/m³). R(C_s) is the rate of the dissolution reaction, for single step irreversible reaction $R\left(C_{\text{s}}\right)=k_{\text{s}}{C}_{\text{s}}$ (m・kmol/s m³). And the Mass balance model for calcium carbonate gives:

(5)

α is the dissolving power of the acid, defined as mass of solid dissolved by unit mass of reacted acid (kg/kg). ρ_s is the density of the solid phase (kg/m³).

2.3. Partial Liquid Equilibrium

The reaction between carbonate rock and aqueous of HCl acid system involves the procedures, the transport of the acid molecules from the bulk fluid to the rock surface, and the chemical reaction at the rock surface. And the overall reaction is controlled by the slower procedure. Almost all the previous work of simulating wormhole propagation in carbonate rocks uses hydrogen ion concentration particularly to describe chemical reaction from the equilibrium given by:

(6)

While in this high TDS (Total Dissolved Solids) system, we have to use hydrogen ion activity in equilibrium conditions. Alkattan et al. (1998) gave the developed function of overall calcite dissolution rate of hydrogen ion activity as follows [28] .

(7)

where k_c is the mass transfer coefficient (m/s). k_s is the reaction rate constant (m/s). ${\gamma }_{{\text{H}}^{+},\text{s}}$ is the activity coefficient of hydrogen ion.

When $k_{\text{c}}\gg k_{\text{s}}{\gamma }_{{\text{H}}^{+},\text{s}}$ , the overall calcite dissolution rate is kinetically controlled, then $R\left(C_{\text{s}}\right)\approx k_{\text{s}}{\gamma }_{{\text{H}}^{+},\text{s}}{C}_{\text{f}}$ . But when $k_{\text{s}}{\gamma }_{{\text{H}}^{+},\text{s}}\gg k_{\text{c}}$ , the overall calcite dissolution rate is mass transfer controlled, then $R\left(C_{\text{s}}\right)\approx k_{\text{c}}{C}_{\text{f}}$ . The activity coefficient of hydrogen ion ( ${\gamma }_{{\text{H}}^{+},\text{s}}$ ) in the process of calcite and limestone dissolution can be computed from:

(8)

The designate parameters f_γ, m_Cl, B_HCl, E and C_HCl in Equation (8) are given by Alkattan et al. (1998) [28] .

2.4. Pore Scale Model

The dissolution process determined local quantities in porous media, therefore the modified Garman-Kozeny correlation models are introduced to describe the relationship between porosity, permeability, interfacial area and pore radius.

(9)

where a₀ is the initial interfacial area of the medium (m²/m³). K₀ is the initial permeability of the reservoir (10⁻³ × μm²). r_p is the pore radius of the medium (m). r_p0 is the initial pore radius of the medium (m). β is the exponent determined by the experiment (Dimensionless). ${\varphi }_0$ is the porosity of the reservoir (Dimensionless).

Panga et al. (2005) proposed the formulas to describe the relationship between Sherwood number and local mass transfer coefficient [10] , k_c. and the diffusion coefficient determined by molecular diffusion and convection diffusion in x direction and y direction respectively as are listed as Equation (11) and Equation (12).

(10)

(11)

(12)

where Sh is the Sherwood number, which is defined as the ratio of convective to diffusive mass transport (Dimensionless). Sh_∞ is the asymptotic Sherwood number (Dimensionless). m is the ratio of pore length to pore diameter (Dimensionless). Re_p is the pore scale Reynold’s number (Dimensionless). Sc is the Schmidt number (Dimensionless). |U| is the magnitude of v and u (m/s). λ_y is the constant depending on pore geometry (λ_y = 0.1 for packed-bed of spheres). λ_x is the constant depending on pore geometry (λ_x = 0.5 for packed-bed of spheres). a_os is the constant depending on pore geometry (a_os = 1.0 for packed-bed of spheres). D_m is the effective molecular diffusivity of acid (m²/s). D_ex is the effective dispersion coefficient in the acid phase in the x direction (m²/s). D_ey is the effective dispersion coefficient in the acid phase in the y direction (m²/s).

2.5. Energy Model

In order to study wormhole propagation in carbonate matrix acidizing affected by temperature, the models considering reaction heat, which are used to describe a thermal non-equilibrium between the solid and acid phases, are introduced. We assumed that the exothermic heat of reaction between acid and rock was transported from the solid phase to the acid phase, and the reaction rate based on the temperature is mainly determined by the solid temperature, T_s, since the reaction happened at the surface of the rock.

(13)

$\begin{array}{l}\frac{\partial }{\partial t}\left(\left(1-\varphi \right){\rho }_{\text{s}}{C}_{\text{Ps}}{T}_{\text{s}}\right)\\ =\nabla \cdot \left(\left(1-\varphi \right)k_{\text{es}}\cdot \nabla T_{\text{s}}\right)+\Delta H_{\text{r}}\left(T_{\text{s}},1\text{atm}\right)a_{\text{v}}R\left(C_{\text{s}}\right)-h_{\text{c}}{a}_{\text{v}}\left(T_{\text{s}}-T_{\text{f}}\right)\end{array}$ (14)

where T_f is the temperature of the acid phase (K); T_s is the temperature of the solid phase (K); C_Pf is the heat capacity of acid (J/kg・K); C_Ps is the heat capacity of solid (J/kg・K); k_ef is the thermal conductivity of acid (W/m・K); k_es is the thermal conductivity of solid (W/m・K); h_c is the acid to solid heat transfer coefficient (W/m² K); ΔH_r (T_s, 1 atm) is the acid-rock molar reaction enthalpy (J/mol). As for the reaction heat depending on the temperature and pressure in porous media [29] , the standard molar reaction enthalpies of HCl acid and limestone at porous surface is given by:

(15)

where ΔH_r (298.15 K) is the molar reaction enthalpy of reaction at 298.15 K (J/mol); C_p,m is the constant-pressure molar heat capacity (J/K mol);

The standard molar reaction enthalpy (obtained from the Handbook of Physical Chemistry [30] ) of each parameter in Equation (15) is listed in Table 1.

Table 1. Parameters of thermodynamic in Equation (15).

The simplified expression of ΔH_r (T_s, 1 atm) can be achieved through substituting values of Table 1 into Equation (15).

(16)

2.6. Reaction Kinetics

The chemical reaction kinetics of HCl acid with various concentrations and temperatures on MISSAN limestone from Iraq has been studied in this paper with rotating disk apparatus at a constant pressure and shear rate. And the functions are used to calculate the reaction rate with the changing temperature and the consumption of hydrochloric acid, which are obtained by regressing experimental data.

Since the reaction rate decreases with the consumption of HCl acid, the function describing the relationship between the concentration of HCl acid and the reaction rate is obtained by regressing laboratory test data shown in Figure 1.

(17)

Figure 1. Relationship of concentration and reaction rate@ 100˚C; 5 MPa; 500 r/min.

The temperature is another factor which affects the reaction rate significantly. In the acid-rock reaction process, the reaction and heat transfer between acid and solid make the changes of reaction rate very complicated. By regressing lab tests on the temperature and reaction rate shown in Figure 2 and integrating

Figure 2. Relationship of temperature and reaction rate@ 20% wt HCl acid; 5 MPa; 500 r/min.

Equation (17), an approximation function describing the reaction rate vs. acid concentration and temperature is obtained as

$k_{\text{s}}=0.021{\text{e}}^{-18616/RT}{C}_{\text{f}}^{1.0341}$ (18)

2.7. Initial and Boundary Conditions

Initial condition

(19)

Pressure boundary condition

It is assumed that the total injection capacity is constant at the inlet of the carbonate core during acidizing. However, the injection capacity for each grid block changes with time because of the permeability improvement. Li (2004) listed out the way to solve the constant flux boundary [31] .

(20)

Boundary condition of HCl concentration

(21)

Boundary condition of acid temperature

(22)

Boundary condition of solid temperature

(23)

where Q_inj is the injection capacity (m³/min). Q_i is the injection capacity at each grid block at the inlet of core (m³/min). P_e is the constant pressure at the outlet of core (MPa). is the mass concentration of the acid at the inlet, (kmol/m³). is the initial solid temperature (K).

3. Computational Methods

The pressure, temperature and HCl concentration for each core segment at each time step are solved numerically by using boundary and initial conditions in Equations (20)-(23). The partial differential Equations ((3), (4), (13) and (14)) are performed with spatial and time discretization through implicit finite difference method to obtain the corresponding ODEs. And we programed an algorithm to save the coefficients of ODEs into a matrix by using MATLAB. The results of pressure, temperature and HCl concentration in each time step are calculated using Ax = B, where x is the vector to be calculated. The detailed solving procedures are listed as: 1) initialize the mesh of core and generate primary distribution of porosity through the method developed by Liu et al. (2012) [15] ; 2) calculate the pressure profile and velocity distribution in the core at time, n + 1, by Equations (1)-(3); 3) calculate the molar reaction enthalpy of reaction, ΔH_r (T_s, 1 atm), in each core segment at time, n, according to Equation (16); 4) substitute ΔH_r (T_s, 1 atm) into Equation (14), and calculate the temperature profile of acid phase and solid phase by solving the integration of Equation (13) and Equation (14) using coupled method; 5) calculate the reaction rate constant, k_s, the mass transfer coefficient, k_c, and molecular diffusivity, D_e, for each core segment; 6) calculate the HCl concentration profile at time, n + 1, through Equation (4); 7) solve Equation (5) and Equation (9) in sequence, the porosity, permeability, pore radius and interfacial area can be obtained; 8) judge: if the acid breaks the core (when the inlet pressure drops to 1% of its initial value, breakthrough phenomenon is considered [14] ), then end the calculation; else go to step 2).

4. Analysis of Simulation Results

The temperature is one of the important factors which influence the wormhole pattern in matrix acidizing in carbonate reservoirs, thereby affecting the optimum injection capacity and the acid volume to break the damage zone of formation [32] [33] [34] [35] . Notably, there are no reports that have simulated the wormhole propagation process under different injection temperatures of acid in carbonate reservoirs. To understand the injection temperatures on wormholing, we focus on studying the acid volume of breakthrough and optimum injection capacity under different temperatures of injected acid using the simulation method with the 2-D vision. The analysis of simulation results is based on the conditions listed in Table 2.

4.1. Analysis of Temperature Profiles under Non-Isothermal Conditions

Figure 3 illustrates the wormhole patterns under different injection capacities with 293.15 K temperature of injected acid, and 393.15 K initial temperature of

Table 2. The corresponding values of parameters in simulation.

Figure 3. Wormhole patterns under non-isothermal conditions using different injection capacity@ (a) Q = 0.2 ml/min; (b) Q = 1.5 ml/min; (c) Q = 10 ml/min; T_inj = 293.15 K; $T_{\text{s}}^0=393.15\text{K}$ .

solid phase. Comparing to the wormhole patterns under the isothermal conditions proposed by Maheshwari et al. (2013) [12] , the typical wormhole patterns can also be observed under the non-isothermal conditions as shown in Figure 3. The conical wormhole, dominant wormhole and uniform dissolution are generated under the corresponding injection capacity 0.2 ml/min, 1.5 ml/min and 10 ml/min in Figures 3(a)-(c).

Figure 4 and Figure 5 illustrate the corresponding temperature profiles of acid phase and solid phase when the acid breaks the core shown in Figure 3. With the colder acid entering the core, the solid phase is cooled down through the thermal conduction between acid phase and solid phase. Under the low injection capacity, the temperature of solid phase decreases more significantly than that of the moderate and high injection capacity by comparing Figures 5(a)-(c). It might be because the breakthrough volume of injected acid is much more under conical wormhole condition or face dissolution condition than that under dominant wormhole condition or uniform dissolution condition. By comparing the temperature in the wormholes for acid phase and solid phase, the temperature along the wormholes near the inlet of core nearly equals to the injection temperature, 293.15 K, and the temperature along the wormholes in the acid phase is slightly lower than the solid temperature because of the reaction on the solid surface. It is interesting to note that the maximum temperature is found at the tips of wormholes since the majority of the reaction occurs here, which agrees with the modeling results of Kalia et al. (2009) [24] .

Figure 4. Temperature profiles of acid phase when the acid breaks the core@ (a) Q = 0.2 ml/min; (b) Q = 1.5 ml/min; (c) Q = 10 ml/min; T_inj = 293.15 K; $T_{\text{s}}^0=393.15\text{K}$ .

Figure 5. Temperature profiles of solid phase when the acid breaks the core@ (a) Q = 0.2 ml/min; (b) Q = 1.5 ml/min; (c) Q = 10 ml/min; T_inj = 293.15 K; $T_{\text{s}}^0=393.15\text{K}$ .

4.2. Effects of Acid-Rock Reaction Heat on Effluent Temperature

A series of studies have shown that the reaction heat of acid-rock should be considered when simulating acidizing progress [23] [24] [26] . All of them assumed acid-molar reaction heat as a constant value while it will change with the temperature and pressure as illustrated in Equation (18). The difference of reaction heat treating methods will affect the predicted temperature to some extent, therefore affecting reaction rate. In order to compare the available model (with constant reaction heat) and present theoretical model (with dependent reaction heat), the temperature in the effluent of core changing with time under different injection rates and acid concentrations were analyzed. The values of the constant parameters in calculation are given by: T_inj = 298 K, $T_{\text{s}}^0=365\text{K}$ , ΔH_r = 4.86 kJ/mol [23] .

Figure 6 illustrates the effluent temperature calculated by constant and dependent reaction heat models under different injection rates: 0.2 m³/min, 1 m³/min and 10 m³/min. The curves indicate that the effluent temperature declines from 365 K to about 298 K under both reaction heat treating methods in different injection rates. This is because the colder acid supplying at the inlet of core persistently replaces the heated acid in the core. However, the temperature decreases faster with the increase of injection rate. It is contributed to two points. In the first place, the velocity of cold supplying acid is faster in the higher injection rate, leading to a faster replacing rate of heated acid. The second reason is that the high injection rate increases the heat transfer rate between acid and rock. The trends of temperature are the same under both reaction heat treating models. Nevertheless, the temperatures treated by dependent reaction heat model are much higher than that of the constant reaction heat model by comparing the line-graph and dash-graph in Figure 6. Moreover, raising the injection rate from 0.2 m³/min to 1 m³/min and 10 m³/min, the maximum differences of temperature for both simulation results experience an increase from 1.32˚C to 8.12˚C and 13.57˚C as shown in stripe-graph in Figure 6.

Figure 6. Effluent temperature of different reaction heat models under varying injection rate@ 0.2 m³/min, 1 m³/min and 10 m³/min.

Figure 7 shows the variation of effluent temperature declines solely caused by persistently injecting cold acid at the inlet of core under different acid mass concentrations, 5wt%, 15wt% and 30wt%. When the constant reaction heat model is chosen, the effluent temperature under strong acid concentration is a little higher than that of the weak one in the whole acidizing process as shown in the line-graph of Figure 7. It is caused by the phenomenon that the reaction heat will increase with the raising acid concentration. And Guo et al. (2014) also achieved the same conclusion in the study of acid fracturing [27] . When it comes to the dash-graph drawn by the simulation results from the dependent reaction heat model, the trends and conclusions are almost the same as the constant reaction heat model did before. However, the difference of effluent temperature between strong acid concentration and weak acid concentration is more remarkable than that of the line-graph. Furthermore, with the increase of the acid concentration from 5wt%, to 15wt% and 30wt%, the maximum differences of temperature for both treating models catch an increase from 1.12˚C to 6.17˚C and 12.16˚C as shown in the stripe-graph of Figure 7.

Figure 7. Effluent temperature of different reaction heat models under varying acid mass concentration@ 5wt%, 15wt% and 30wt%.

According to the discussions above, an obvious difference was observed by the simulation results of effluent temperature through both treating models, especially in high injection rate and strong acid concentration. This is mainly caused by the phenomenon that the high injection rate will stimulate the molar heat transfer rate, and the strong acid concentration will make the reaction heat increase. The difference of effluent temperature indirectly proves that the temperature profiles in core are almost different, therefore affecting reaction rate, wormhole propagating velocity, wormhole patterns in the matrix acidizing of carbonate reservoir. Moreover, the dependent reaction heat model takes more comprehensive factors into account and more conforms to reality, therefore, it is of great importance to consider to apply the dependent reaction heat model in the future work referring to the matrix acidizing simulation.

4.3. Effects of Injection Temperature on Pore Volume of Breakthrough

Under the isothermal conditions, we assume that the injected acid temperature equals to the solid temperature, and they are uniform throughout the core medium. Therefore, the concentration of HCl acid is the only driving mechanism for wormhole propagation under the certain initial temperature of acid/solid. Figure 8 plots the acid volume required to break the core for various acid injection capacities under 293.15 K, 338.15 K and 393.15 K. As the acid efficiency curves show, no matter what the isothermal conditions are, the PV_BT decreases as the injection capacity increases when it is smaller than the optimum injection capacity. However, it experiences a relatively slow growth with the increase of injection capacity when it exceeds the optimum injection capacity. On the other hand, the optimum injection capacity increases with the increase of initial temperature of acid/solid, and the corresponding optimum injection capacities are 0.7 ml/min, 1.2 ml/min and 1.7 ml/min for the temperatures 293.15 K, 338.15 K and 393.15 K, respectively. This is because that reaction rate decreases with the decreasing temperature from Equation (18), and the injection capacity should be lowered to make sure the convection and reaction become comparable to generate the dominant wormhole.

Figure 8. Volume of acid required to break the core for various acid injection capacity under isothermal conditions@ 293.15 K, 338.15 K and 393.15 K.

The wormhole propagation process is more complicated when the non-isothermal conditions are considered. The reaction rate in each core segment at each time step is influenced not only by the acid concentration, but also by the temperature significantly, which will affect the volume of acid breakthrough and optimum injection capacity. Figure 9 plots the acid efficiency curves with different injection temperatures under non-isothermal conditions. Different from the isothermal conditions shown in Figure 8, there is no such tendency that the optimum injection capacity increases with increasing injection temperature, although the optimum injection capacity exists in each group of simulation results. And they obtain a nearly similar value for the optimum injection capacity, 1.5 ml/min. However, the injection temperature affects PV_BT under any injection capacity, especially under low injection capacity. As we can see from the curves, the maximum difference of PV_BT is 46.67 PV which occurs at 0.2 ml/min, and the minimum value is 0.31 PV at 1.5 ml/min.

Figure 10 shows the pressure drop as a function of time for non-isothermal cases with different injection temperatures using the optimum injection capacity obtained from Figure 9. All of the pressure drops decrease almost linearly with time until the acid breaks the core, which means the dominant wormholes are generated in those cases [36] . On the other hand, the pressure drop decreases more slowly with increasing injection temperature under optimum injection

Figure 9. Volumes of acid required to break the core for various acid injection capacity under non-isothermal conditions@ T_inj = 293.15 K, T_inj = 338.15 K and T_inj = 383.15 K.

Figure 10. Pressure drop for non-isothermal conditions using optimum injection capacity.

Figure 11. Volumes of acid breakthrough for non-isothermal conditions using optimum injection capacity.

capacity by comparing the curves in Figure 10. However, it is observed that the trend of PV_BT does not agree with that for the pressure drop velocity. The PV_BT experiences a decrease with the increasing injection temperature, and it is followed by an increase with the increasing injection temperature as shown in Figure 11, so that an optimum injection temperature also exists.

5. Conclusions

In the present work, the dependent reaction heat model was developed and introduced into the integrated governing models of matrix acidizing in linear core. And the effluent temperatures were compared on both reaction heat treating methods. Finally, the acid wormholing influenced by temperature under isothermal and non-isothermal conditions was studied by this simulation method. The conclusions coming from the calculation and analysis could be summarized as: 1) Since there is an important difference of simulation results between the dependent reaction heat model and constant reaction model, and the former takes more comprehensive factors into account and more conforms to reality, therefore, it is suggested to use the dependent reaction heat model when modeling the acidizing process; 2) Under isothermal conditions, the optimum injection capacity decreases with the decreasing temperature; 3) Under non-isothermal conditions, the wormhole patterns, i.e. uniform dissolution, ramified wormhole, dominant wormhole, conical wormhole and face dissolution, can also be generated by changing the injection capacity, and the optimum injection capacity is slightly affected by the injection temperature; 4) Similar to the optimum injection capacity, there also exists an optimum injection temperature for developing dominant wormhole under non-isothermal conditions.

In the present work, the influence of supercritical CO₂ has been ignored, which will affect the simulation results to some extent. When the supercritical CO₂ is considered, it will affect the reaction heat definitely, and thereby affect the wormholing patterns. By the way, the two-phase flow should be considered in the modeling instead of one-phase flow in the present work when the supercritical CO₂ is considered.

The wormhole propagation process in acidizing carbonate is very complex, especially when using the unconventional acid system such as in-situ self-diverting acid and in-situ cross-linked acid etc. It requires further studies to capture the details of wormholing of unconventional acid system by considering the affection of temperature. The work presented in this paper can be improved by including the effects of supercritical fluid CO₂; the fracture-vug-matrix medium and the non-Newtonian acid. In the future studies, some of the extension work will be pursued.

Acknowledgements

The work is supported by grants from the National Natural Science Foundation of China, No. 51474182 and International Science and Technology Cooperation Research Program of Sichuan Province (China) under Grant, No. 2017HH0014.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Daccord, G. and Lenormand, R. (1987) Fractal Patterns from Chemical Dissolution. Nature, 325, 41-43. https://doi.org/10.1038/325041a0
[2]	Daccord, G. (1987) Chemical Dissolution of a Porous Medium by a Reactive Fluid. Physical Review Letters, 58, 479. https://doi.org/10.1103/PhysRevLett.58.479
[3]	Hoefner, M.L. and Fogler, H.S. (1988) Pore Evolution and Channel Formation during Flow and Reaction in Porous Media. AIChE Journal, 34, 45-54. https://doi.org/10.1002/aic.690340107
[4]	Golfier, F., Zarcone, C., Bazin, B., et al. (2002) On the Ability of a Darcy-Scale Model to Capture Wormhole Formation during the Dissolution of a Porous Medium. Journal of Fluid Mechanics, 457, 213-254. https://doi.org/10.1017/S0022112002007735
[5]	Huang, T.P., Hill, A.D. and Schechter, R.S. (1997) Reaction Rate and Fluid Loss: The Keys to Wormhole Initiation and Propagation in Carbonate Acidizing. International Symposium on Oilfield Chemistry, 18-21 February, Houston, SPE-37312-MS.
[6]	Gdanski, R. (1999) A Fundamentally New Model of Acid Wormholing in Carbonates. SPE European Formation Damage Conference, The Hague, 31 May-1 June 1999, 107-116.
[7]	Buijse, M.A. (1997) Understanding Wormholing Mechanisms Can Improve Acid Treatments in Carbonate Formations. SPE European Formation Damage Conference. Society of Petroleum Engineers, The Hague, 2-3 June 1997, SPE-38166-MS.
[8]	Fredd, C.N. and Fogler, H.S. (1999) Optimum Conditions for Wormhole Formation in Carbonate Porous Media: Influence of Transport and Reaction. SPE Journal, 4, 196-205. https://doi.org/10.2118/56995-PA
[9]	Liu, X., Ormond, A., Bartko, K., et al. (1997) A Geochemical Reaction-Transport Simulator for Matrix Acidizing Analysis and Design. Journal of Petroleum Science and Engineering, 17, 181-196. https://doi.org/10.1016/S0920-4105(96)00064-2
[10]	Panga, M.K.R., Ziauddin, M. and Balakotaiah, V. (2005) Two-Scale Continuum Model for Simulation of Wormholes in Carbonate Acidization. AIChE Journal, 51, 3231-3248. https://doi.org/10.1002/aic.10574
[11]	Kalia, N. and Balakotaiah, V. (2007) Modeling and Analysis of Wormhole Formation in Reactive Dissolution of Carbonate Rocks. Chemical Engineering Science, 62, 919-928. https://doi.org/10.1016/j.ces.2006.10.021
[12]	Maheshwari, P., Ratnakar, R.R., Kalia, N., et al. (2013) 3-D Simulation and Analysis of Reactive Dissolution and Wormhole Formation in Carbonate Rocks. Chemical Engineering Science, 90, 258-274. https://doi.org/10.1016/j.ces.2012.12.032
[13]	Smith, M.M., Sholokhova, Y., Hao, Y., et al. (2013) CO 2-Induced Dissolution of Low Permeability Carbonates. Part I: Characterization and Experiments. Advances in Water Resources, 62, 370-387. https://doi.org/10.1016/j.advwatres.2013.09.008
[14]	Kalia, N. and Balakotaiah, V. (2009) Effect of Medium Heterogeneities on Reactive Dissolution of Carbonates. Chemical Engineering Science, 64, 376-390. https://doi.org/10.1016/j.ces.2008.10.026
[15]	Liu, M., Zhang, S. and Jianye, M. (2012) Fractal Nature of Acid-Etched Wormholes and the Influence of Acid Type on Wormholes. Petroleum Exploration and Development, 39, 630-635. https://doi.org/10.1016/S1876-3804(12)60086-X
[16]	Izgec, O. (2010) Numerical and Experimental Investigation of Acid Wormholing during Acidization of Vuggy Carbonate Rocks. Journal of Petroleum Science & Engineering, 74, 51-66. https://doi.org/10.1016/j.petrol.2010.08.006
[17]	Cohen, C.E., Ding, D., Quintard, M., et al. (2008) From Pore Scale to Wellbore Scale: Impact of Geometry on Wormhole Growth in Carbonate Acidization. Chemical Engineering Science, 63, 3088-3099. https://doi.org/10.1016/j.ces.2008.03.021
[18]	Bulgakova, G.T., Kharisov, R.Y., Pestrikov, A.V., et al. (2013) Experimental Study of a Viscoelastic Surfactant-Based in Situ Self-Diverting Acid System: Results and Interpretation. Energy & Fuels, 28, 1674-1685. https://doi.org/10.1021/ef4019542
[19]	Liu, M., Zhang, S., Mou, J., et al. (2013) Diverting Mechanism of Viscoelastic Surfactant-Based Self-Diverting Acid and Its Simulation. Journal of Petroleum Science and Engineering, 105, 91-99. https://doi.org/10.1016/j.petrol.2013.03.001
[20]	Mou, J., Cai, X., Liu, M., et al. (2013) The Mechanism of Leak off Reduction of Clean Self-Diversion Acid in Acid Fracturing. IPTC 2013: International Petroleum Technology Conference, Beijing, 26-28 March 2013.
[21]	Liu, P., Xue, H., Zhao, L.Q., et al. (2015) Analysis and Simulation of Rheological Behavior and Diverting Mechanism of in Situ Self-Diverting Acid. Journal of Petroleum Science and Engineering, 132, 39-52. https://doi.org/10.1016/j.petrol.2015.04.042
[22]	Ratnakar, R.R., Kalia, N. and Balakotaiah, V. (2013) Modeling, Analysis and Simulation of Wormhole Formation in Carbonate Rocks with in Situ Cross-Linked Acids. Chemical Engineering Science, 90, 179-199. https://doi.org/10.1016/j.ces.2012.12.019
[23]	Liu, M., Zhang, S., Mou, J., et al. (2013) Simulation of Temperature Field in Near-Wellbore Region in Carbonate Acidizing. Ciesc Journal, 64, 3542-3548.
[24]	Kalia, N. and Glasbergen, G. (2009) Wormhole Formation in Carbonates under Varying Temperature Conditions. 8th European Formation Damage Conference, Scheveningen, 27-29 May 2009, 19. https://doi.org/10.2118/121803-MS
[25]	Bousri, A., Bouhadef, K., Beji, H., et al. (2012) Heat and Mass Transfer in Reactive Porous Media with Local Nonequilibrium Conditions. Journal of Porous Media, 15, 329-341. https://doi.org/10.1615/JPorMedia.v15.i4.30
[26]	Kalia, N. and Glasbergen, G. (2010) Fluid Temperature as a Design Parameter in Carbonate Matrix Acidizing. SPE Production and Operations Conference and Exhibition, Texas, 8-10 June 2010.
[27]	Guo, J., Liu, H., Zhu, Y., et al. (2014) Effects of Acid-Rock Reaction Heat on Fluid Temperature Profile in Fracture during Acid Fracturing in Carbonate Reservoirs. Journal of Petroleum Science and Engineering, 122, 31-37. https://doi.org/10.1016/j.petrol.2014.08.016
[28]	Alkattan, M., Oelkers, E.H., Dandurand, J.L., et al. (1998) An Experimental Study of Calcite and Limestone Dissolution Rates as a Function of pH From 1 to 3 and Temperature from 25 to 80 °C. Chemical Geology, 151, 199-214. https://doi.org/10.1016/S0009-2541(98)00080-1
[29]	Medeiros Jr., F. and Trevisan, O.V. (2006) Thermal Analysis in Matrix Acidization. Journal of Petroleum Science and Engineering, 51, 85-96. https://doi.org/10.1016/j.petrol.2005.11.011
[30]	Perry, R.H. and Green, D. (1997) Chemical Engineers Handbook. 7th Edition, McGraw Hill Professional, London.
[31]	Li, C. (2004) Fine Scale Sandstone Acidizing Coreflood Simulation. The University of Texas, Austin, 37-38.
[32]	Wang, Y., Hill, A.D. and Schechter, R.S. (1993) The Optimum Injection Rate for Matrix Acidizing of Carbonate Formations. SPE Annual Technical Conference and Exhibition, Houston, 3-6 October 1993, 13. https://doi.org/10.2118/26578-MS
[33]	Frick, T.P., Kurmayr, M. and Economides, M.J. (1994) An Improved Modeling of Fractal Patterns in Matrix Acidizing and Their Impact on Well Performance. SPE Production & Facilities, 9, 61-68. https://doi.org/10.2118/23789-PA
[34]	Fredd, C.N. and Fogler, H.S. (1999) Optimum Conditions for Wormhole Formation in Carbonate Porous Media: Influence of Transport and Reaction. SPE Journal, 4, SPE-56995-PA.
[35]	Bazin, B. (2001) From Matrix Acidizing to Acid Fracturing: A Laboratory Evaluation of Acid/Rock Interactions. SPE Production & Facilities, 16, 22-29. https://doi.org/10.2118/66566-PA
[36]	Alghamdi, A.H.A. (2010) Experimental and Theoretical Study of Surfactant-Based Acid Diverting Materials. Texas A&M University, College Station.