Numerical Simulation of Low Salinity Water Flooding on Core Samples for an Oil Reservoir in the Nam Con Son Basin, Vietnam

Low-salinity water flooding (LSWF) is environment-friendly and operates similarly to conventional waterflooding without the need for synthetic chemical materials. The application of LSWF makes sense in Vietnam as HC production has steadily declined since 2002, and the majority of main oil fields have become near mature and mature fields. In the next years, Enhanced Oil Recovery (EOR) should be a top priority for Petro Vietnam to boost its oil production, for which the key issue is how to select a suitable EOR technology. In this study, LSWF of the Lower Miocene sand using low salinity water from Lower Oligocene sand was investigated. Previously at the Ruby field in the Cuu Long Basin, an LSWF feasibility study was carried out based on a conventional core flooding experiment, which is time-consuming and costly. This study targets the Chim Sao field in the Nam Con Son Basin, for which a cheaper and faster assessing method is required. As a result, a numerical code written in Matlab was developed and successfully validated with the core flooding experiment results obtained at the Ruby field. The LSWF simulation was conducted using the multiple ion-exchange mechanisms (MIE), and the results obtained showed an increase in the oil recovery factor by 2.19% for the Lower Miocene Sand. Another important outcome of this study is the innovative proposal and successful simulation to use the abundant low salinity water from the underlying Lower Oligocene sand as a natural LSW source to inject into the Lower Miocene oil reservoir that can be a decisive factor to help apply LSWF in practice on a wide scale not only for Chim Sao but also other similar oil fields in southern offshore Vietnam.


Introduction
For a long period since 1987, the main oil production of Vietnam has come from the famous White Tiger fractured granite reservoir in the Cuu Long Basin, which peaked in 2006 and then started to decrease gradually. On the other hand, contributions from the clastic Miocene and Oligocene reservoirs in the Cuu Long Basin and the Nam Con Son Basin become increasingly significant (see Figure 1). While the petroleum production from the Miocene sands is second to that from the fractured granite basement reservoirs since 2013, it also has started to reduce. The recovery factor of the Miocene reservoirs ranges from 15 to 35%. As a matter of fact, most major oil fields in Vietnam are in a declining phase, which requires an intensive study on Improved Oil Recovery (IOR) or Enhanced Oil Recovery (EOR) [1] Thang et al. (2021) in [1] reported that many major oil fields in Vietnam are in the declining phase and can be classified as mature or near-mature fields, which need the application of IOR/EOR to sustain their production life. Thang and Giao (2017) [2] made an overview of EOR applications in Vietnam. Many studies related to the application of EOR using polymer, water-alternative-gas injection, thermally stable surfactant and nano-surfactant, namely [3][4][5][6], have been carried out. An EOR master plan was proposed based on a recent comprehensive study by Quy et al. 2020 [7], in which LSWF is considered as one of the viable methods, as shown in Figure 2 below: Figure 2. Low-salinity water flooding (LSWF) as one applicable method in the enhanced oil recovery (EOR) master plan of Petro Vietnam after Quy et al. 2020 [7].
The main objective of this study is to develop a numerical code in Matlab to simulate the 1D LSWF flooding process on core samples to investigate its feasibility in a mature oil field, i.e., the Chim Sao field, in the Nam Con Son Basin (see Figure 3a). The simulation in this study focused, however, only on the multiple ion exchange (MIE) mechanism so that the simulation results can be compared and validated by the laboratory core flooding experiment results obtained at the Ruby oil field in the Cuu Long Basin [8]. The scopes of this study include the followings: (i) performing a comprehensive overview of LSWF mechanisms and applications; (ii) screening and selection of a suitable study location in the Nam Con Son Basin; (iii) development of a one-dimensional numerical simulation of LSWF and validating it with the core flooding results from the Ruby field in the Cuu Long Basin, where an experimental study of LSWF was recently conducted using the conventional laboratory core flooding experiment; (iv) Performing a one-dimensional LSWF simulation using the newly developed code for the Lower Miocene sand using the natural low salinity water source from the below Oligocene sand for the Chim Sao oil field, the Nam Con Son Basin.

LSWF Mechanism
The mechanisms of LSWF have been studied in the laboratory scale more than two dedicates and can be classified as: (1) fines migration [9,10]; (2) mineral dissolution [11]; (3) increased pH effect and reduced interfacial tension [12]; (4) emulsification [8,13]; (5) multiple ion exchange [13][14][15]; (6) double-layer effects [16,17]; and (7) wettability alteration [18][19][20]. The following part will focus more on the multiple ion exchange (MIE) as one of the main mechanisms of LSWF employed in this study (see Figure 4), which is highly dependent on compositions of the formation water and injected brine. For LSWF to be effective, it is necessary for the formation (connate) water to contain sodium, calcium and magnesium ( [13][14][15]21]). In an MIE process, the ions like K, Na will replace the divalent ions (e.g., Ca, Mg) that are adsorbed on the rock surfaces, leading to the liberation of oil in the form of calcium carboxylate complex (see Figure 4). In addition, the ion-exchange process explains why LSWF does not work when a core is acidized and fired when the cation exchange capacity (CEC) of the clay minerals is destroyed and why LSWF has little effect on refined oil because no polar compounds are present to strongly interact with the clay minerals ( [16]).
An important observation reported from core flooding experiments, and field implementations ( [22][23][24]) is that the effluent pH tends to increase during LSWF. This phenomenon can be explained by the dissolution of carbonate minerals, such as calcite and dolomite. Ionic exchange during LSWF leads to the adsorption of divalent ions (Ca 2+ , Mg 2+ ) and promotes mineral dissolution, which provides further calcium and magnesium for ion exchange. This explanation is consistent with the results of Evje and Hiorth in 2011 [25] and Hiorth et al. 2010 [26].

LSWF Governing Equations
Integrated simulation of geochemical reactions and fluid flow with active ions in the porous media, especially LSWF for EOR, is essential for understanding the behavior of each chemical component in the fluid. PhreeQC, a USGS's geochemical reaction package [28], can help simulate various geochemical processes, including equilibrium between water and minerals, ion exchangers, surface complexes, solid solutions, and gases. Consequently, PhreeQC can be used to couple with the other multiphase flow source codes, one of which is the Matlab reservoir simulation toolbox (MRST) developed by SINTEF to prototype new models and computational methods for reservoir simulation [29].
In this study, first, the black-oil model, as described by Bao et al. (2017) [30], is a special multicomponent, multiphase model with no diffusion among the fluid components. The name "black-oil" refers to the assumption that various hydrocarbon species can be lumped together to form two components at surface conditions and, namely, a heavy hydrocarbon component called "oil" and a light component called "gas." At the reservoir conditions, these two components can be partially or completely dissolved in each other, depending on the pressure, forming a liquid oleic phase and a gaseous phase. In addition, there is an aqueous phase, which herein is assumed to consist of only water. Similar to the polymer model given by Bao et al. 2017 [30], the continuity equations for the LSWF model are given in the Equation (1) as follows: where: φ is rock porosity; s α , ρ α and v α are saturation, density and velocity of phase α, respectively; c is the concentration of a certain chemical component given in mass per volume of water; φ re f and ρ r are reference porosity and density of the rock.
The phase fluxes velocity, v α , are estimated from Darcy's law given as follows: where: K is absolute permeability of the reservoir rock; λ α = k rα µ α is the mobility of phase α; k rα is relative permeability and µ α is the phase viscosity.
The fully implicit simulation with low and high order schemes have been proposed and solved Equations (1) and (2) by Bao et al. 2017 [30], and Mykkeltvedt et al. 2019 [31] and the finite volume method was used to solve numerically.

Finite Volume Formulation of Ion Transport Equation
To simulate an LSW core flooding experiment, a one-dimensional LSW flow model coupled with chemical reactions will be employed in this study. Figure 5 shows a finite volume (FV) discretization of the core with specified boundary conditions at the core end. The sequence of core flooding to be simulated include the followings steps: (i) First, (flow simulation), low salinity water is injected at a constant pressure higher than the reservoir pressure at the left end (x = 0), while at the right end (x = L) the pressure is kept constant as the reservoir pressure; (ii) Second, (simulation of occurring chemical reactions), the chemical reactions of transported ions of the formation (connate) water with the ions in the injected LS water will be simulated, and as the results, the salinity of the mixed fluid will be calculated, based on which the relative permeability of hydrocarbon will be determined as a function of salinity that is discussed in the literature by Shojaei et al. 2015 [32]. By using the 1D model and neglecting compressibility, we can rearrange the flow Equation (1) with chemical component transport in Buckley-Leverett form as follows: in which f w is the fractional flow, which is defined as: where p c is the capillary pressure, ∆ρ wo the density difference between water and oil, g the gravitational acceleration, ϑ the angle inclination of the reservoir. For the case of core flooding, the capillary pressure is negligible, and the core sample is placed horizontally, Equation (4) can be reduced to: The Buckley-Leverett equation is discretized in one-dimension space of length L with a mesh of N + 1 grid points x i for i = 0, 1, 2, . . . , N is shown in Figure 5. All grid cells have ∆x = l/N width and are of the cell-centered type.
Taking integration of Equation (3a) over a single grid cell, one gets: where s w,i is average water saturation at each grid cell i.
Using the first accurate Euler forward method for the left-hand side (LHS) of Equation (7), then for the spatial discretization, the first-order upwind is applied to get the fully explicit form of the water saturation as follows: In the matrix form, Equation (8) is written as follows: where s n w is N × 1 vector of the water saturation of all grid cells at a time step n; K v is N × N matrix that is defined as follows: and f n w and f n bc are N × 1 vectors given as: A second-order time discretization was introduced using the modified Euler method [33], which can capture both information at the beginning of a time step (equal to the Euler forward method) and the information at the end of a time step. Additionally, this approximation also mitigates the small oscillations happened when the relative permeability curve switch at the shock front. Equation (9) can be rewritten according to the modified Euler method as follows [33]: where f n+1 w can be determined by first calculating the predictor step of water saturation, s n+1 w as given in Equation (9).
In fact, the above-mentioned modified Euler method is a special case of Runge-Kutta method for solving ODE (ordinary differential equations) [34]. For the single chemical component transport equation (Equation (3b)) by taking the derivative of a functioning product, Equation (3b) becomes: It is noted that, physically, the water saturation of each cell (S w ) is higher than that of the connate water saturation (S wc ), i.e., s w ≥ s wc > 0, so that Equation (15) is valid.
Theoretically, Equation (15) can be discretized in a similar way the Buckley-Leverett advection equation was done as mentioned above. A problem arising with first-order spatial discretization is that it often exhibits large numerical dispersion. Especially for multi tracers augmented reservoir simulation, this has been discussed in the works by Alsofi and Blunt (2010) [35], Jerauld et al. (2006 and2008) [36,37]. In fact, one of the ways of limiting numerical dispersion is to increase the number of grid cells, which is often inapplicable as it leads to a very large amount of grid cells in augmented reservoir simulations. In order to overcome the dispersion in a first-order upwind system, the higher numerical schemes can be used while keeping the realistic amount of grid cells as given in the following equation [38]: However, higher-order numerical schemes suffer from oscillation near discrete steps. That is, if one tries to simulate a slug of low salinity concentration with the numerical scheme from Equation (14) will be heavily distorted and lead to faulty conclusion and non-physical results [39].
The removal of oscillations from higher-order schemes can be done by using the so-called flux limiters [40]. The procedure consists of two steps. In the first step, the 2nd order upwind scheme is rewritten as a correction to the monotone first-order upwind difference and the correction term is written in the form of successive gradients as in the followings: In the second step, both non-monotonous terms are multiplied by the flux limiter functions Ψ(r i ) and Ψ(r i−1 ) where There is no suitable flux limiter function for all cases. A full review of all mathematical requirements and their implementation is summarized by Hirsch, 2007 [40].
Equation (17) can now be fully solved numerically. However, it is still used for the first-order accuracy in the time domain. As second-order accuracy is desired, the same modified Euler method that was applied to solve Equation (3a) as mentioned above will be applied to Equation (17), and the final numerical prediction of a single time step, as shown in the following: the conversion criteria known as Courant-Friedrichs-Lewy (CFL) condition as described by Chen, 2008 [41], is given as follows: Consequently, for the case of the Buckley-Leverett equation (Equation (3a)), the derivative of the fractional flow with respect to the water saturation is maximum in order to satisfy the CFL condition, which means: For the cases of chemical components (ions) transport (equation 3b), the stability condition will be: Figure 6 shows the flowchart in this study to simulate the process of two-phase fluid transport with chemical reactions in a porous media that was used for simulating a core flooding experiment with constant porosity and permeability. In this workflow, the simulation starts at time t once the initialization is made. After this, the fractional flow ( f w ) and water saturation (s w ) are implicitly solved by Equation (12), and the concentration equation of each chemical component will be further solved by Equation (18). When the concentrations of chemical components at a grid cell are calculated, they will be transferred to the IphreeQC [42], which is designed to increase the flexibility in interfacing with the widely used PhreeQC. Consequently, the newly obtained concentrations of these components/ions after chemical reactions will be compared with a certain salinity threshold to determine whether the injected water is of high or low salinity and to make the changes of relative permeability of water and oil as represented by Corey's coefficient [43]. In the studies by Shojaei et al. 2015 [32]; Jerauld et al. 2008 [37] and Tripathi and Mohanty, 2008 [4], Corey's coefficients (n o ) are a function of salinity, X c , but above a certain threshold of salinity (high salinity X HS c ) and below a certain level of salinity (low salinity, X LS c ) it has no salinity dependency as described below:

Validation of the LWSF Simulation Results with Those of a Core Flooding Experiment
Based on the FV formulation of equation 18 and its coupling with chemical reactions, as shown by the workflow in Figure 6, a fully functional code was developed in Matlab with IphreeQC functions being called via Com-server for the chemical reaction's simulation to implement the simulation of LSW core flooding. To validate the new code, we used X c LS = 1000 ppm, X c HS = 7000 ppm for low and high salinity, respectively, following Tripathi and Mohanty's (2008) model [44]. The validation of the new code was done by comparing the simulation results with those of a core flooding experiment conducted for the Ruby field in the Cuu Long Basin. The core sample and fluid properties and ion concentrations, including oil, high salinity (HS) and low salinity (LS) water, are presented in Tables 1 and 2. The core flooding experiment first started with HS water injection; after that, at a pore volume of 20 PV, an LS water will be injected. The setup of an LWS core flooding is shown in Figure 7 [8], in which the reservoir pressure (2570 psi) and temperature (93 • C) in an oven are kept constant. First, the cores were saturated with FW and cured for 7 days; then, the n-decane was injected until irreducible water saturation was achieved. In the next step, stock tank oil (STO) was injected until the n-decane was perfectly replaced, and the cores were cured by the STO for 40 days. The injection rate was 0.1 cc/min, equivalent to 1 ft/day of flow velocity in the reservoir. The stable outlet pressure of the core was controlled by using a back pressure regulator (BPR), and the differential pressure was measured from changes in the inlet and outlet pressure gauges. The recovered oil and effluent were collected at fixed intervals. After a flooding test, recovery was confirmed using a Dean-Stark apparatus. In addition, cation concentrations (Ca 2+ , Mg 2+ , Na + and K + ) and pH in the effluent collected were analyzed using the inductively coupled plasma atomic emission spectrometry (ICP-ES) method. Figure 8 shows a chart of relative permeability for high and low salinity water. The Corey oil and water coefficients are presented in Table 1. Figure 9 compares the recovery factors by LSWF simulation and experiment, which match quite well, in particular after the low salinity water is injected. Both simulation and experiment confirmed that an increase of about 2.5% in oil production due to the low salinity water injection method. Similarly, the simulation results of ion transport for Na + , K + , Ca 2+ , Mg 2+ , Cl − and pH, using Equation (21) with and without chemical reactions for Cations Na + , K + , Ca 2+ , Mg 2+ , anion Cl − and pH are shown in Figure 10a-f, respectively. There are two stages of water injection, including formation water injection (high salinity) and low salinity water injection, whose fluid properties and ion concentrations are shown in Tables 1 and 2. In this core flooding experiment and simulation, the divalent ions (calcium and magnesium) are expected to exchange with sodium and potassium that creates an effect of the MIE mechanism for LSWF. To examine the effect of ion exchange, the simulations were run with and without an ion-exchange reaction. The simulated and experimental concentrations of Na, K, Ca and Mg are shown in Figure 10a-d, respectively. Figure 10d helps to show that the effluent Mg++ concentration in the case with ion exchange (blue line) matches better with the experimental one in comparison with the case without ion exchange (orange line). The simulated and experimental concentration of passive ion (chloride) also matches quite well, as seen in Figure 10e. In addition, Figure 10f shows an increasing trend of pH, which is commonly mentioned in the literature.

Geology of the Study Site
The Cuu Long Basin is separated from the Nam Con Son Basin to the SE by the Con Son Swell (Figure 11). The geologic history of the Cuu Long Basin can be divided into four phases as follows: (i) The pre-rift phase (Late Jurassic-Paleocene) with widespread emplacement of granitic plutons and NE-trending swarm dykes. These structural fabrics constitute the basis for the granitic basement play in the Cuu Long Basin; (ii) Rifting phase (Eocene ?-early Late Oligocene): in this period, the basin was opened in NW-SE direction, which reactivated earlier NE-SW trends as listric normal faults. These normal faults form grabens and half-grabens filled initially with coarse-grained sediments. As rifting continued and waned at the late Early Oligocene-early Late Oligocene, the basin expands into a widespread lacustrine environment; (iii) Compression phase (early Late Oligoceneearly Early Miocene): Combined with pre-rift and syn-rift trends, the compression phases created pervasive fracture network in the granitic basement, which is crucial to the success of the basement play. Furthermore, during this time, lacustrine deposition is widespread and enabled the formation of the D shale, which is both an excellent source and seal for the basement play. Deep lake environment prevails throughout most of the D sequence and part of the C sequence; (iv) Thermal subsidence (late Early Miocene-Recent): After the compression phase, the basin entered a period of tectonic quiescence. The first major marine influence is marked by the Bach Ho shale ("Rotalid shale"), which can be correlated over the whole basin and serves as an excellent regional seal for the basin. Thermal subsidence and marine influence continued to the present day [45,46].  Figure 11 shows a long interbasin geological cross-section from the Cuu Long Basin, through the Con Son Swell, to the Nam Con Son Basin. The pre-Tertiary basement may be 6 to 8 km deep. The Ruby Field is located in the most northern part of Cuu Long Basin. The Ruby structure is an anticlinal drape that overlies one of these intra basin horst features. Hydrocarbon in the Ruby Field was discovered in the Early Miocene, Oligocene as well as the basement. The Nam Con Son Basin is separated from the Cuu Long Basin by a broad SW-NE trending basement ridge called Con Son Swell. It is about 550 km long and 200 km wide across its widest southern part. From the granitic outcrop on Con Son Island, it is believed that some part of Con Son Swell is built of Late Cretaceous granitoids. On the other hand, the Chim Sao Field is located towards the southwestern end of the Nam Con Son Basin. The stratigraphic section penetrated in the Chim Sao Field ranges in age from the recent strata to the Oligocene formation. Elsewhere in the Nam Con Son Basin, a basement of Late Jurassic to Early Cretaceous weathered and fractured granites and granodiorites have been penetrated. Hydrocarbon in Chim Sao Field was found in Early and Middle Miocene sands.

Results and Discussion
As the coreflooding experiment for LSWF is quite expensive and time-consuming in general, and for the Chim Sao field in particular, we have investigated the feasibility of LSWF using a numerical simulation approach in this study. First, the availability of low salinity water in this oil field is checked based on the reservoir petrophysical properties and production information. Figure 12 shows the ion concentrations of formation water of the main reservoir in the Chim Sao field, indicating that the Lower Oligocene formation water can be a good candidate for LSW source as its total salinity is essentially below 2000 ppm. In the first step of LSWF modeling for the Chim Sao field, the new code was applied to simulate a core flooding experiment with the input data of reservoir and fluid properties being shown in Table 3. It is noted that the real ion concentrations of LSW from Lower Oligocene sand (the possible natural LSW source) and the formation water of the injected reservoir (Lower Miocene sand) were used in our LSWF simulation for the Chim Sao field. Figure 13 shows the results of simulated ion concentrations during the HS and LS injection stages, which show a very similar trend to those in Figure 10 for the Ruby field, where the injected water was an artificial one. Figure 14 shows the results of the recovery factor calculated for a two-stage LSWF simulation in the Chim Sao field, first with high salinity (HS) water flooding up to 4.4 PVI for maintaining reservoir pressure, then followed up by an LS water flooding up to 8 PVI. The simulation results obtained indicated that an EOR taking the natural LSW from Lower Oligocene to inject into the Lower Miocene oil sand could result in an increase of the recovery factor by 2.19%.

Conclusions and Recommendations
In this study, a comprehensive review of LSWF mechanisms was done. In line with a master EOR plan of PetroVietnam [7], the low salinity water flooding is considered as one of the possible EOR options to be applied for the mature and near-mature oil fields in Vietnam. To meet the need to quickly assess the LSWF feasibility of these oil fields, a new computation code using the finite element volume (FV) method was developed in Matlab to simulate a core flooding experiment, and it was successfully applied for the Chim Sao oil field in the Nam Con Son Basin. The major conclusions drawn from this study are as follows: 1.
Based on the results of a laboratory core flooding experiment conducted for the Ruby oil field in the Cuu Long Basin concerning a low salinity water flooding (LSWF) of the Lower Miocene sand, the multiple ion exchange (MIE) was identified by [8] as one of the most suitable mechanisms for this oil sand. The simulation conducted in this study also has confirmed the same for the Lower Miocene sand in the Nam Con Son Basin; 2.
A new computation code based on the finite volume element (FV) method was successfully developed in this study to simulate the ion transport in the aqueous in coupling with chemical reaction. The numerical simulations were run, and the results were compared with and are well-validated by the experimental results of the laboratory core flooding experiment for the Ruby field in the Cuu Long Basin in term of oil recovery factor, multi-ion concentrations (Na + , K + , Mg 2+ , Ca 2+ , Cl − and H + ); 3.
This study's major outcome was the successful simulation of LSWF using the natural low salinity water from the Lower Oligocene formation to inject into the overlying Lower Miocene sand, a main oil productive unit in the Chim Sao oilfield, the Nam Con Son Basin. The simulation results obtained show that such an LSWF could increase the oil recovery factor by 2.19%; 4.
It is recommended the application of the new code be applied for assessing the other mature oil fields in terms of applicability of LSWF to enhance their oil recovery, and some LSWF pilot full-scale field tests should be planned in the near future in the Nam Con Son Basin, Vietnam. In addition, further development of the code should consider the general form of flow equation (Equations (1a) and (1b)) in coupling with PhreeQC for chemical reactions the simulation to take care of the capillary pressure effect.