A simulation study on hydrogel performance for enhanced oil recovery using phase-field method

Hydrogels are increasingly applied in oil recovery processes. This leads to more controlled flow of fluids in porous media. In this process, hydrogel is injected to the reservoir to block the high permeability areas. The trapped oil in low permeability regions, is then swept by water flooding. pH‐sensitive hydrogel microspheres were synthesized in another work of the authors, which effectively increased the oil recovery factor in experimental studies. In this communication, phase-field approach was used to simulate this process and to obtain the tuning parameters of the model including thickness of the contact surface (є), phase transform parameter (M0), and excess free energy (∧). Diffusion of hydrogels was studied by Cahn–Hilliard conservative approach and the breakage, deformation, and plugging mechanisms were analyzed, based on pressure drop variations in micromodel. Moreover, Effective parameters on oil recovery factor were analyzed. Results indicated a good agreement between experimental and modeling studies of oil recovery factor in water and hydrogel flooding with absolute errors of 2.29% and 4.06%, respectively. The recovery factor was calculated using a statistical method which was in good agreement with the modeling results. The tuned parameters of the model were reported as, є = 111.7 µm, M0 = 5 × 10−13 m3/s, ∧=-0.0003\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\wedge = - 0.0003$$\end{document} J/m3.

Pre-cross-linked polymers. This category of polymers includes microgels, performed particle gels (PPG), colloidal dispersion gels (CDG), bright water, and pH sensitive gels. Polymers are prepared at surface to be injected in the reservoir. Dispersed aggregates are injected to the reservoir and swell in the matrix in area of high permeability. Unlike in-situ cross linked polymers, they do not have the limitations of changing the gelant composition, its dilution, and reaction control problems. On the other hand, they might be filtered within the pores of the reservoir or mechanical retention might lead to high pressure drop in the system. This might lead to poor injectivity, swelling problems and plugging undesirable regions of the reservoir. Various studies have been implemented for application of pre-cross-linked polymers in the near well-bore region and deep in the reservoir rock 27 .
Gel and polymer flooding simulation studies are divided into macroscopic and microscopic categories. In macroscopic approach, there are various studies based on classical percolation theory. Yuan et al. 28 proposed a 2D multi-component model for polymer cross-linking in the reservoir. Chen et al. 29 performed a simulation study on oil displacement mechanism. Based on this model, Xu et al. 30 proposed a 3D model with eight components for gel flooding. The results for saturation were improved using high-order explicit techniques. A 3D streamline method was established by Feng et al. 31 to consider the capillary and gravity forces, and applied it for gel flooding simulation studies. Wu et al. 32 and Cui et al. 33 improved the proposed gel flooding models by coupling the chemical reaction and mass transfer between phases. An optimization study was performed on injection parameters for gel flooding, using starting pressure gradient method. This study was then developed for fractured reservoir by Zhao et al. 34 using fracture-matrix coupled technique in a 3D two-phase flow. An optimization study was also performed for profile control in this reservoir. Herzig and Payatakes 35 established a classical percolation theory. In their study, percolation coefficient was a function of percolation velocity and concentration. Particle deformation and expansion were considered by Zhang et al. 36 in which a profile control model was proposed for fluid-solid coupling. PPG displacement was modeled by Wang et al. 37 , considering the relationship between pore throats and particles. Restarting pressure and deformation were modeled in this study. A power law model was proposed to calculate the restarting rate. Moreover, a distribution function was recommended to estimate the probability of pore plugging. Surfactant migration model was coupled with PPG flooding to describe the interaction between surfactant and PPG. One of the limitations of classical percolation theory, is that the breakthrough time is always considered one pore volume. Moreover, the particle size, deformation, and strength are not studied in this theory. Therefore, size exclusion theory was introduced to develop more accurate results. Flow modeling in pore scale is very useful to understand the transport of the fluids in the reservoir. Interfacial force could be modeled as a continuum force using diffuse interface models. In this technique, numerically resolvable layers are applied for distributing the interface discontinuities on them 38 . These methods are numerically attractive due to advantage of easier solving the Navier-Stockes equation on fixed grid in comparison to exact equation, which requires adaptive interface fitting grids. Diffuse interface models are classified into three main groups including tracking/ distributed force model, continuum surface force technique [39][40][41] , and phase-field based models [42][43][44] . Santos 45 and Bedrikovetsky 46 established a population equilibrium model in which particles penetrate to large pore throats or might be captured by the reservoir rocks. A modified colloidal model was proposed by You et al. 47 in which inlet concentration was increased by the particles in accessible pores. This model led to more appropriate results for fitting the experimental data. Liu et al. 48 upgraded the single-phase flow to two-phase flow model in which PPG migration was modeled based on the pore distribution.
Another approach for studying the displacement mechanisms in microscale, include molecular dynamics. In a recent study, snapshots of molecular dynamic simulation was applied for analyzing the whole displacement process. Laplace equation was applied to configure the oil-water interface 49 . This technique gives an insight to residual oil displacement in reservoir and provides a vision in EOR studies 50 www.nature.com/scientificreports/ of migration process in turbulent flow for irregular shaped particles. In this study, flow was analyzed using LBM and Smagorinsky turbulent technique. Particles interaction was characterized using DEM and particle-fluid interaction was studied using intrusive boundary conditions. Ohtsuki et al. 54 applied a hybrid simulation model for studying the particle concentration effect on permeability. The fluid-solid interaction was studied using a coupled system. LBM was incorporated to model, for studying the fluid flow in mesoscopic level. Mechanical and geometric characteristics were analyzed using DEM at microscopic level. Interaction between solid and liquid was modeled based on boundary method. A new solver was introduced by Xiong et al. 55 for simulation of particle-fluid interactions in large scale. Zhou et al. 56,57 introduced a new approach for simulation of deformable PPG migration using LBM, IBM, and DEM. In this study, a functionality was developed between particle-throat diameter, elastic modulus, and critical pressure gradient. Phase-field method (PFM) is an interface capturing technique, that is widely applied due to its simplicity and accuracy for complex geometries and moving interfaces 58,59 .
In this paper, an attempt was made to use the phase-field model for simulation of hydrogel injection. This model makes it possible to track particles during hydrogel injection in micromodel. Phase-Field model together with Conn-Hilliard and Navier-Stokes flow equations makes it possible to simultaneously study on various effective parameters in flooding. These parameters include, interfacial tension, wettability, density and viscosity of various fluids, hydrogel swelling, and etc.,. Deformation of hydrogel particles and their breakage is well modeled during the movement in micro-model. Results indicated that large particles are broken at inlet of the micro-model and deformed particles with smaller diameter are moved more deeply. Also in this study, we tried to show the real conditions of the reservoir, in which the high-permeability and low-permeability areas are not considered as separated sand packs, but are continuously connected to each other and fluid transfer is possible between them. This model also allows to study the effects of temperature and gravity. Moreover, in this study swelling ratio is represented as follows: where W s represents the weight of the swollen hydrogel at a specified time and W d denotes the weight of the dried hydrogel. The swelling ratio was obtained using experimental studies. In simulation studies, the free energy of the system was defined in such a way that increasing the volume of hydrogels due to swelling was applied in the model.

Materials and methods
Experimental setup of hydrogel injection process is explained in the work of Jamali et al. 60 and illustrated in Fig. 1. Oil samples were obtained from Siri reservoir with specific gravity of 0.85, °API gravity of 32.27, and viscosity of 16 cp. pH-sensitive hydrogel microspheres containing silica nanoparticles were synthesized. The hydrogel solution namely, poly(acrylamide-co-methylenebisacrylamide-co-acrylic acid) was prepared with concentration of 1000 ppm and homogenized with an ultrasonic probe with power of 30 kW for 15 min. Micromodel tests were designed for injection of hydrogel with optimum formulation based on the work of Jamali et al. 60 in which the mixer round was 443.5 rpm, AAcid/AAm wt.% was 62%, and Si nano-particle content was 0.25%.

Micromodel experiments.
A two-layer homogeneous micromodel was designed to study on increasing the oil recovery factor using the synthetized hydrogels. To do tests, water was injected to a saturated micromodel with oil. Then water was injected to micromodel. The injected water was channelized through high permeable zone. Oil recovery factor in this step was 51.46%. Hydrogel solution was then injected to micromodel. Since hydrogel particles were swelled, and high permeable region was blocked, the flow was diverted to the low permeability area, and the pressure difference between input and output was increased. For the last step, water was injected to recover the trapped oil in low permeability region. Experimental values of the recovery factor from the work of Jamali et al. 60 are represented in Table 1.

Simulation of hydrogel flooding
Simulation studies was performed in two steps. Water flooding that was followed by hydrogel injection. Comsol Multiphysics software was applied for simulation of gel flooding process using phase-field approach. Phase-field model. Phase-field method (PFM) is a widely applied technique for complex geometries and moving interfaces, due to its simplicity and high accuracy 58,59 . Since phase-field method take use of advection techniques, they are easily implemented in three dimensions by finite element method and unstructured grids. Dissipative energies and surface tension could be easily implemented in this method. Constitutive relations including relative permeability and capillary pressure would be beneficial for flow simulation 61 . The interfaces should have four to eight cells to calculate surface tension and interface energy. There are two categories for phase-field model including conservative Cahn-Hilliard model in which diffusion is driven by the chemical potential 62 , and non-conservative Allen-Cahn model 63,64 . In phase-field model, advection term is implemented without distortion in three dimension in unstructured grids.
Two-phase flow equations. In this method, the two-fluid contact surface is considered as a diffusive thin layer represented in Fig. 2. In fact, the diffusive contact surface model provides a way for continuous modeling of contact surface forces, whereby the discontinuities at contact surface are propagated on a thin, numerically soluble layer. In this method, transfer equations at contact surface are replaced by continuous transfer-diffusion  www.nature.com/scientificreports/ equations, where the emission rate is determined by the chemical potential gradient (ϕ) and interfacial tension is calculated based on the interfacial disturbance energy. Using these equations, contact surface displacement and deformation are calculated on a fixed network. Using the surface diffusion model for Navier-Stokes multiphase flow equations is much easier than solving them with precise equations. Cahn-Hilliard method is applied for studying conversion, formation and dissolution at contact surface 65 . In this method, the function of the free energy is minimized, which is represented in Eq. (1).
where ϕ represents the chemical potential and u is the fluid velocity. The mass concentration and chemical potential are related by the following Eq. 38 : Chemical potential is equivalent to the surface tension in curvature for phase-field model. m represents the mass concentration of the component. G denotes the chemical potential and is represented as follows: where M c is the specific mobility and describes the transmission stability due to diffusion. є describes the thickness of the contact surface between the two fluids. λ is the mixing energy density In phase field equations, viscosity and density are calculated as follows: where ρ 1 and ρ 2 are densities of fluids 1 and 2. μ 1 and μ 2 are the viscosity of fluids 1 and 2 and v is the volume fraction of each phase. Since the fluids are considered incompressible, mass and volume conservation laws are applied 62 . Average curvature (1/m) is calculated using the following equation: Cahn-Hilliard equation has a curvature-dependent solution that makes it possible to use it for simulation of nucleation, evaporation, and swelling processes. Areas of high curvature at contact surface are generally high potential areas, and soluble materials move through these areas to low potential surrounding zones 38 . While the fluid transfer is involved, the chemical potential would not necessarily be uniform. Precision analysis of phase-field model is complicated, since convergency is provided by three factors including mesh size, interface surface thickness, and mobility. Mobility influences the thickness and turbulence of the boundary layers for the chemical potential. The convergence rate is determined by a dual limit, which combines the approximation of phase-field model to the real sharp surface contact physics as well as approximation of the numerical methods of phase-field model to its real solution.
Three-phase flow equations. In this model, the contact surface is tracked between the phases by solving the continuity and Navier-Stokes equation. The contact surface of the fluid is analyzed using chemical potential and phase-field models. The free energy is minimized to define the liquid-liquid contact surface motion. Navier-Stokes is represented as follows for the three-phase flow: where P is the pressure, µ represents viscosity, F shows the external force. I is the identity matrix. F st is the surface tension calculated as follows: η i denotes the chemical potential of each phase. Boyer technique 66 was applied to track the three phase flow contact surface.
(1) ∂ϕ ∂t + u · ∇ϕ = ∇ · (M∇G). where σ ij is the surface tension between phases, i and j. ∧ denotes the excess free energy of the system. i is defined as follows: Model implementation. The properties of the two parts of micromodel are represented in Table 2. In this study, the modeling approach was finite element in Comsol Multiphysics software. For walls as impermeable surfaces with no-slip boundary condition, we have: where ṫ is the tangent unit vector to the surface, u f is the tangential fluid velocity and u w is the tangential wall velocity. Since the walls are stationary, the fluid velocity along the impermeable walls is zero ( u fṫ = o).
n is the normal vector to the surface and the fluid velocities at interface,re equal ( u f 1 = u f 2 ). Boundary layer mesh generation was applied in this study. This type of mesh is suitable where the fluid flow is coupled with energy and mass transfer equations. For solving the Navier-Stokes equations, the modified Newton's method was implemented. The solution starts with the initial guess for u 0 .

Result and discussion
In this study, gel flooding process was simulated using Comsol Multiphysics, based on the phase-field method to study the mechanism of hydrogel performance for increasing the oil recovery factor. The developed model was validated with experimental results. Sensitivity analysis was performed on various parameters regarding micromodel and hydrogel, to propose the functionality of the recovery factor with these variables for water flooding and gel flooding processes.
Validation of experimental results for water flooding. The modeling result is represented in Fig. 3 where oil and water are specified in blue and red, respectively. The results are in good agreement with experimen-  www.nature.com/scientificreports/ tal data from the work of Jamali et al. 60 . Comparison between recovery factor values for modeling and experimental studies are presented in Table 3. As it is obvious, in layer of low permeability, capillary force is increased and the injected fluid flows to high permeability area. Water injection is continued until the breakthrough time.
Sensitivity analysis on oil recovery factor. A statistical model was developed for prediction of oil recovery factor in water flooding process. In this section, dual effects of a number of parameters were considered.
Statistical model for prediction of recovery factor in water flooding. Sensitivity analysis was performed to determine the effective parameters on oil recovery factor including injection rate, wettability, and viscosities of oil and water. Using these parameters, a statistical model was proposed. The results for ANOVA analysis of variance are reported in Table 4. As it is obvious, P value is less than 0.0001, which denotes the high accuracy of this model. Coefficients of the model are reported in Table 5.
In this equation, q represents the flow rate (cc/min), Ɵ is contact angle (degree), µ o denotes oil viscosity (cp), µ w is water viscosity (cp), and T represents temperature (°C). As it is obvious in Fig. 4, the modeling results are in good agreement with the proposed model using design expert, since most of the data are nested near the unit slop line curve.
Capillary number. Capillary number (Ca) defined as the ratio of viscous force to surface tension, is one of the dimensionless parameters that is used in flow interpretation. For capillary numbers more than 1, it is increased with increasing the velocity and is negligible for low Reynolds numbers. Capillary number is in the order of 10 −6 and is about 1 in the near wellbore region and is defined as follows: where γ is the surface tension. Ca is increased with increasing the oil recovery factor. Figure 5 illustrates the capillary number diagram in terms of normalized saturation for the modeling results in the present study.
Capillary and viscous forces are in contradiction with each other during oil recovery processes. The capillary force traps the oil in the cavities, while viscous force acts to displace oil by the injected fluid. High sweeping efficiency is an important parameter for successful EOR process which is a function of capillary number. In fact, to improve the recovery factor, capillary number is increased and the mobility ratio is reduced. According to Fig. 5, for low capillary numbers, the diagram is linear with zero slope. In critical point which is a function of    www.nature.com/scientificreports/ wettability, permeability and pore structure, a sharp reduction is observed. Accurate representation of capillary number functionality with normalized saturation is an essential tool for simulation of oil sweeping and minimizing the residual oil saturation. Capillary number diagram in terms of normalized saturation is mostly related to reservoir heterogeneity, but it is recommended to use its model for investigation on sweeping process.

Hydrogel flooding simulation.
As represented in previous sections, after water breakthrough, the prepared hydrogels are injected to micromodel. By injection of hydrogels, the high permeability region is blocked and water is transferred to non-swept areas of low permeability. Figure 6 illustrates the modeling result using Comsol Multiphysics software which are in good agreement with experimental study. In simulation studies, density of hydrogels was equal to that of displacing fluid and its viscosity was 80 cp. The injection time was 1260 s in this step. The injection temperature was 25 °C, the outlet pressure was 10 5 pa. Injection rate was 0.0005 ml/min and wettability was 40°.
During injection of hydrogels, they might break or transform. Table 6 represents the RF% for experimental and modeling studies. Thickness of the contact surface in simulation studies was obtained 111.7 µm. The phase transform parameter in Cahn-Hillard equation (M o ) was equal to 5 × 10 −13 m 3 /s, and the excess free energy ( ∧ )  Figure 7 illustrates the variations of oil recovery factor for the whole process.
Mechanisms for hydrogel movement in porous media. The mechanisms of hydrogel movement in porous media include breaking, deformation, and blocking, which are illustrated in Fig. 8a-c, respectively. All of these mechanisms are effective on the performance of the injected hydrogels, which are reviewed in this section.   www.nature.com/scientificreports/ Hydrogel breakage is mostly occurred within the throats near injection point. The broken hydrogels are then penetrate more deeply in porous media. Hydrogels are broken by increasing the pressure to more than a threshold value, otherwise hydrogels would block the throats. Pressure variations at two sides of micromodel is represented in Fig. 9. Two regions are defined in this figure. In the first zone, hydrogel particle blocks the throat, temporarily. The pressure is increased at inlet and reduced at outlet. This leads to increasing the pressure difference in micromodel. As a result, hydrogel particle is deformed and broken (Fig. 8a-part b). The two generated particles are get to the next pore (Fig. 8a-part c) and one of the particles enter to the next throat (Fig. 8a-part d). As a result, the outlet pressure is slightly increased (Zone II). In all the steps for breakage and deformation, gels would plug the throat, if the pressure is not increased enough. Hydrogels might also cross the throats by deformation, which actually happen at injection point and in medium depths of the reservoir. Figure 10 shows pressure variations at two sides of micromodel for deformation of hydrogels. Two regions are defined in this figure. In the first region, hydrogels are blocked the throat. This causes increasing the pressure difference in this region. This pressure difference leads to deformation and movement of hydrogel particles in the throat. The movement of hydrogels is continued toward the next pore and throats (Fig. 8b-part d). This is obvious with an incremental increase in outlet pressure (region II).
While the pressure variations is not enough for deformation or breakage of hydrogels, the blockage is the predominant mechanism. Blockage is mostly occurred at depth of the porous media by reducing the pressure gradient. Figure 11 shows pressure variations at two sides of micromodel for the mechanism of blockage. Hydrogel deformation is occurred at the first region by pressure variations (region I). For pressures less than the threshold value, blockage is occurred at the throat and pressure value at two sides of micromodel remains almost constant (region II). For the first region of Fig. 11, the pressure difference between two sides of micromodel is increased since hydrogels block the throat. The pressure increase leads to deformation of hydrogel, but the pressure is not high enough for breaking the hydrogels or passing them through the throats. Since the throat is blocked, the inlet and outlet pressure values, remain constant which is obvious in region II.
In this section, pressure difference between the two ends of micromodel is discussed for hydrogel flooding process. The first section of Fig. 12 represents the pressure variations in micromodel for water flooding at various www.nature.com/scientificreports/ time steps. The pressure drop is reduced during water flooding. The next part represents the pressure variations for hydrogel flooding. After water breakthrough, during hydrogel injection, the pressure difference between inlet and outlet, represents considerable increase. Blockage of some of the throats during hydrogel injection causes this phenomenon. Moreover, there is distinctive pressure fluctuations after hydrogel injection. This is due to permanent or temporary throat blockage during hydrogels injection.
Pressure profile diagram. To obtain pressure profile along the micromodel, the pressure variations along line A is investigated (Fig. 13). The pressure changes during water flooding is represented in Fig. 14 at various time intervals. As shown in the figure, when the initial parts of micromodel are swept, the pressure variation along the micromodel is decreased at various time steps and water is replaced by oil. At 1320 s, the pressure profile with almost monotonous slope is obtained. This means that breakthrough has occurred.
Sensitivity analysis on oil recovery factor for hydrogel injection. In this section, effective parameters on oil recovery factor in hydrogels injection process are analyzed. These parameters include mean diameter of hydrogel particles, hydrogel injection rate, and hydrogel viscosity.    www.nature.com/scientificreports/ Hydrogels mean radius. Hydrogel size is an important parameter in determining the penetration depth and hydrogel performance in increasing the oil recovery factor. The final radius of the injected hydrogel is a function of its initial size and the swelling ratio. Depending on the method of hydrogel synthesis and its chemical compounds, the blocking ability is different. Bai et al. 17 stated that particles four times the throat diameter, are able to pass the throat. Larger particles should be broken before passing the throat. However, the operating conditions and hydrogel chemistry are of important parameters. Depth of penetration is a function of particle size, since it is effective on injectability of the fluid. In the present study, effect of hydrogel radius on oil recovery factor was studied. The mean radius of hydrogel particles in experimental study was 69.9 µm and the average size of throats in high and low permeability areas were 350 and 240 µm, respectively. Sensitivity analysis are performed at constant injection rate of 1e-3 cc/min and viscosity of 80cp. The results are represented in Fig. 15. Increasing the size of hydrogels led to increasing the percentage of blocked throats in high permeability region, which results in alteration the direction of the injected fluid toward the low permeability region. Therefore, the oil recovery factor would increase in this section and the amount of produced water is controlled. Of course, the size of hydrogels should not be so large that could not enter the pores. By increasing the hydrogels radii, their ability to penetrate to micromodel depth is decreased, and blockage is occurred in areas close to the injection point. This would reduce the sweeping efficiency in highpermeability region. Variations of the recovery factor in the whole micromodel, is also represented in this figure.
Hydrogel injection rate. The functionality of recovery factor with injection rate is analyzed from different standpoints. The first one is the effect of injection rate on functionality of hydrogel particles, and the other is the effect of injection rate on displacement process similar to that of two-phase flow. Result of the two parameters is represented in Fig. 16. In this study, the mean radius of hydrogels was 69.9 µm and the solution viscosity was 80 cp. As represented in this figure, the recovery factor of the high permeability region is increased with increasing the injection rate up to an optimum value (1 × 10 −3 cc/min). The recovery factor in low permeability region is at its minimum value at this point. The total oil recovery factor is increased with increasing the injection rate.
Hydrogel viscosity. Hydrogels are considered as a high viscosity liquid for simulation studies. To simulate solid particles, Mirzaei 67 estimated their viscosity to be about 100 times that of a liquid. In Saghafi et al. 68 study, the viscosity of hydrogels for simulation studies was 80 cp. The flow direction could be changed due to high viscosity of the injected hydrogels. In the present study, viscosity of the injected hydrogels was 80 cp. This value leads to an acceptable penetration depth for hydrogel particles to change the flow direction and increasing the oil recovery factor. Effect of hydrogel particles viscosity was studied on oil recovery factor in high and low permeability regions and the whole micromodel. In this study, the injection rate and hydrogel radius were 1e-3 cc/min and 69.9 µm, respectively. Effect of this parameter on oil recovery factor, is illustrated in Fig. 17. With increasing the viscosity of hydrogels, the flow in high permeability region is stopped near the injection point and leads to deviation of the flow to low permeability region. This leads to reducing the oil recovery factor in high permeability region, since hydrogels could not penetrate to the acceptable depth in micromodel. The flow resistance is increased in high permeability region and the flow is directed to low permeability part which is represented in www.nature.com/scientificreports/ Fig. 18 (path 1). Moreover, hydrogels are accumulated at inlet of high permeability region which steers the flow to low permeability area (path 2). With increasing the viscosity of hydrogel particles, the penetration and placement of the particles in high permeability region is reduced. This leads to reducing the oil recovery factor, since hydrogels are unable to direct the flow to low permeability area.
Statistical model for hydrogel injection. In this section, a statistical model is represented to determine oil recovery factor in hydrogel injection process. This model was proposed based on the effective parameters including, fluid injection rate, hydrogel size and viscosity of hydrogel solution. The results of ANOVA analysis for the high permeability region is represented in Table 7. F-values and P-values represent that the model is significant. The proposed model is represented in Eq. (17). www.nature.com/scientificreports/ where Q represents the injection rate of hydrogels (cc/min), r denotes the hydrogel radius (µm), and μ shows the viscosity of hydrogels (cp). The fitting statistics are reported in Table 8. The same procedure was applied for the low permeability region and oil recovery factor was proposed in terms of injection rate, viscosity, and hydrogel radius by the following equation. The fitting statistics are also reported in Table 9. The results of ANOVA analysis for the low permeability region is represented in Table 10.   www.nature.com/scientificreports/ Figure 19a and b, represent the results of modeling and prediction of oil recovery factor based on the proposed statistical models for high and low permeability regions, respectively. As represented in these figures, most of the data are nested near the unit slope line curve which states the high accuracy of the proposed statistical model.

Conclusion
The present survey, focused on the performance mechanism of hydrogel for enhanced oil recovery. A simulation study was conducted for hydrogel flooding using Comsol Multiphysics software, based on the phase-filed approach. The simulation results for oil recovery factor were validated using the experimental data from another work of the authors. Sensitivity analysis was performed on oil recovery factor, based on the effective parameters in this process. Various mechanisms of hydrogel performance in the porous media were discussed and oil recovery factor was quantitatively analyzed, based on ANOVA analysis. The main results of this study are summarized as follows:  www.nature.com/scientificreports/ • Experimental and modeling results for water flooding and hydrogel injection processes were in good agreement with each other with absolute errors of 2.29% and 4.06%, respectively. • In gel flooding process, the tuning parameters of the model was reported as follows: Thickness of the contact surface (є) was 111.7 µm. Phase transform parameter (M 0 ) was obtained 5 × 10 −13 m 3 /s and excess free energy ( ∧ ) was reported − 0.0003 J/m 3 .
• Increasing the hydrogel radius and viscosity has the same behavior on oil recovery factor. Hydrogel ability for in-depth penetration is reduced, and the blockage is occurred near the injection point. This would reduce the sweeping efficiency in high-permeability region and alteration the flow direction toward the low permeability region. • The recovery factor of the high permeability region is increased with increasing the injection rate, up to an optimum value (1 × 10 −3 cc/min). The recovery factor in low permeability region is at its minimum value at this point. • Results of sensitivity analyses indicated the direct functionality of sweeping efficiency in low permeability region with injection rate, viscosity, and hydrogel radius. This expresses the acceptable performance of injected hydrogels for increasing the total recovery factor.