Optimization and Analysis of CO 2 Huff-n-Puff Process in Shale Oil Reservoirs Using Response Surface Methodology (RSM)

The objective of this study is to determine the inﬂuence of operational parameters of CO 2 Huﬀ-n-Puﬀ EOR process in the Eagle Ford shale oil reservoirs using response surface methodology (RSM). Single-factor analysis was ﬁrst conducted for establishing the Box-Behnken model in RSM. We selected the primary depletion time, gas injection time, cycle number, production time per cycle, and injection rate as the primary input variables using RSM. The cumulative oil production and net present value are optimized as the output factors. After that, Design of Expert 12 software was used to design the experimental table for the above setting factors. Corresponding to the results taken from the optimization, the most signiﬁcant factor is injection rate, followed by injection time, cycle number, and primary depletion time, and production time per cycle is the least signiﬁcant. Additionally, the optimum responses were found as primary depletion time of 2.37 years, injection time of 3.4 months, cycle number of 3, production time per cycle of 2.2 years, and injection rate of 5000MSCF/D. Moreover, correlation coeﬃcient ( R 2 ) of the regression model is 0.9977, and adjusted model R 2 is 0.8898, which further indicates that the model has good reliability. Results show that RSM is a useful technique for optimization, and it also provides insights into optimizing and designing the CO 2 Huﬀ-n-Puﬀ process in the shale oil reservoirs.


Introduction
In recent years, more attention has been paid to the specific potential of unconventional reservoirs. Tight oil and shale oil contribute to the oil production all over the world. e global energy landscape was changed by the successful development of shale oil reservoirs in USA [1,2]. Shale oil reserves are rich with huge development potential in China [3]. However, due to the low porosity and permeability of tight reservoirs, conventional waterflooding is not ideal. CO 2 injection has been proved to be a good choice to improve the oil recovery in the unconventional reservoirs [4,5].
Researchers carry out a series of studies to determine the applicability of CO 2 -EOR process in the unconventional reservoirs [6][7][8]. Song and Yang [9] have conducted several experiments to evaluate the performance of CO 2 injection in enhanced oil recovery from tight oil reservoirs. Results indicate that CO 2 flooding results in a significant increment in the oil recovery. Ding et al. [10] performed core flooding tests to examine the performance of CO 2 Huff-n-Puff process, concluding that CO 2 Huff-n-Puff injection is capable of effectively recovering tight matrix oil, but the performance is impacted by matrix size and permeability. Alfarge et al. [11] and Song et al. [12] found that molecular diffusion is the main mechanism controlling the development of CO 2 -EOR for shale oil reservoirs by combining numerical simulation methods with field production data. Hejazi et al. [13] optimized the cycle number for the CO 2 -EOR used in the Bakken.
A number of numerical simulation studies have assessed the impact of CO 2 concentrations in tight oil reservoirs on development [14][15][16][17]. Gamadi et al. [18] studied the effect of circulating gas injection on shale oil reservoirs to improve oil recovery and found that CO 2 injection was better than N 2 injection. Yu et al. [19,20] conducted numerical simulation of CO 2 Huff-n-Puff EOR process in Bakken tight oil reservoir and found that CO 2 -EOR process significantly improved the oil recovery. ey also pointed out that CO 2 diffusion coefficient is a key parameter affecting the ultimate oil recovery. In the studies by Ma et al. [21], injection rate and injection volume had a great influence on the performance of enhanced oil recovery.
Zhang et al. [22] investigated the factors affecting the CO 2 Huff-n-Puff process in tight oil reservoirs; however, only injection rate and CO 2 diffusion coefficients are considered. Wang et al. [23] adopted the orthogonal experimental design to analyze the influence of parameters in the CO 2 Huff-n-Puff process; however, injection time, soaking time, and cycle number were not included in their models. Moreover, they did not consider the interaction between different parameters. Tang et al. [24] conducted a numerical simulation study on CO 2 injection in tight oil reservoirs, and single-factor analysis was applied to optimize the CO 2 Huffn-Puff schemes. However, it can only screen the listed schemes, and the optimal design cannot be obtained in their model.
As aforementioned, current studies focus on the singlefactor sensitivity analysis, and the mutual interactions between different parameters are always neglected. Additionally, the economic cost is also significant in the design of CO 2 Huff-Puff process [25,26]. erefore, a comprehensive model for optimizing the operational parameters is essential.
At present, the main optimization algorithms include RSM, ANN, and Taguchi orthogonal analysis. In the previous study, Taguchi orthogonal analysis is often used to evaluate the influence factors of CO 2 Huff-n-Puff process, which can determine the general influential factor. However, it cannot optimize the parameters of the actual projects, which means that engineers still need some "empirical" to set the value of parameters. Compared to the traditional Taguchi orthogonal analysis, the RSM can get the multiple regression equation relation according to the result of the design of experiment analysis with multiple inputs and outputs [27]. Within the scope of the experimental design of planning, RSM can find the best combination of the dependent variable in the experiments through the analysis of regression equation. In this way, RSM could determine the optimal combination plan and fill the limitation of Taguchi orthogonal analysis. ANN has a similar optimization ability to RSM, and, in most application cases, the fitting equation and optimization scheme are slightly better than the results of RSM. However, the information processing system in ANN that simulates the biological human nervous system needs more experimental samples before it can be successfully modeled than RSM. Hence, RSM is thought to be a good choice in the optimization analysis [28][29][30][31].
Among the experimental models in design of experiment, RSM is outstanding in many fields such as the chemistry, food science, and biology science [32][33][34][35]. In recent years, RSM has gradually gained favor in the engineering field [36][37][38][39], including in the oil-gas field development [40][41][42]. RSM can not only build a continuous variable surface model with limited experimental data points, but also evaluate the interaction between the influential factors, which fills in the lack of research on the interaction of factors in the current study. erefore, RSM is a useful tool for design and optimization in the CO 2 -EOR operation.
In this work, we first conducted a series of single-factor analysis of operational parameters in the CO 2 Huff-n-Puff process, including injection time, injection rate, soaking time, production time per cycle, and cycle number. e cumulative oil production and economic benefits are selected as the optimization outputs. Box-Behnken design was then utilized to analyze the interaction between operational parameters. Finally, the optimal CO 2 Huff-n-Puff process was obtained, and the economic benefits were also evaluated. By introducing RSM into the optimization of CO 2 Huff-n-Puff process, the impacts of different parameters can be effectively distinguished, and the optimal scheme can be determined, which provides useful guidance for the design and operation of CO 2 -EOR process.

Response Surface Methodology.
In 1951, RSM was formally proposed by British scholars Box and Wilson. With mathematical statistics, design of experiment modeling and statistics, the multiple input variables can be analyzed in RSM.
e objective of this method is to determine and optimize the effects on the dependent variable and degrees of effects of input factors. By RSM, we can use the least number of operations to evaluate the uncertain factors on the quantitative evaluation of the influence of target size, which can save time and money. is method optimizes the responses by deriving the relationship between the input parameters and outputs.
In order to obtain the optimal parameter value of the experimental scheme and optimize the whole project implementation, the following equation can be drawn: where y is the output, x i (i � 1, 2, 3, . . . , z) are the values of input parameters, f(x 1 , x 2 , . . . , x z ) is the regression equation obtained by regression analysis of the impact factors, and ε is the error of the regression equation. e RSM needs to meet the requirement that the test area is close to the optimal area; that is, the designed experimental site should include the optimal experimental conditions. If without the single-factor analysis, improper selection of experimental points will be caused, and good optimization results cannot be obtained by using RSM. erefore, the single-factor analysis is needed to determine the reasonable factors and their ranges of value in the RSM, which is benefit for screening significant factors and reducing the number of experiments.
Box-Behnken design is a spherical design with rotation, which greatly reduces the number of experiments compared with the Central Composite Design experiment. For an experiment with 6-factor and 3-levels, it only needs 54 experiments, and the interaction between parameters can still be considered [43]. With these characteristics, Box-Behnken design is often acted as the preferred experimental design in RSM.
In this study, we applied Box-Behnken design to determine the influential factors of the CO 2 Huff-n-Puff project in the Eagle Ford shale oil reservoirs. According to ANOVA, the significance and applicability of the regression model can be tested and the degrees of the influence on the dependent variable can be determined. e RSM flowchart is shown in Figure 1.

Economic Evaluation.
e main purpose of CO 2 Huff-n-Puff process is to obtain economic benefits from oil and gas extraction and incremental oil recovery. erefore, operational costs should also be considered when designing the relevant parameters of stimulation process. In this study, cumulative oil production and economic benefits are set as the objective to determine a feasible and reasonable CO 2 -EOR project in the Eagle Ford shale oil reservoirs. e cumulative oil production at the end of 21 years is predicted and the dynamic evaluation method is adopted to evaluate the net present value (NPV).
e NPV method refers to the net present value of the net cash flow discounted to the initial construction period. e NPV is greater than 0 when the scheme generates additional benefits, which is expressed as where n is the total prediction calculation period of the investment plan; (C a − C 0 ) t is net cash flow in t th year; e is the discount rate; (1 + e) t is the discount coefficient in the t th year. e dynamic payback period, which is defined as the time required when the sum of investment and income is 0, is calculated. e investment payback time after the project put into production can be determined as where j is internal rate of return.
In the economic evaluation, the following assumptions are made: (1) Excluding the interest of the construction period, depreciation is divided into ten years. (2) Only crude oil is produced from wells. (3) e life period, construction period, and production period of oil wells are determined values. (4) e benchmark discount rate is set as 9%.
e sales revenue is the total revenue calculated by combining the annual oil production at the current crude oil price and the crude oil commodity rate (96%). Based on the oil price in 21 years, the crude oil price of $50 per barrel is now used for this calculation.
As reported, capital expenditure in oil and gas field development includes one-time investment cost of horizontal well and fracturing construction, ground supporting cost of injection equipment, gathering and transportation equipment, power and off-site system, and operating cost in production. It mainly includes comprehensive costs such as employee salary, test maintenance cost, engineering material cost, and oil and gas processing cost. e taxes are mainly calculated as federal income taxes. Expenditure items are shown in Table 1.

Reservoir Profile.
In this study, a typical block in the Eagle Ford is selected for the simulation study. e CMG-GEM module is used to establish the corresponding reservoir model with the size of 7785 ft × 1300 ft × 100 ft. ere are 173 grids in x direction, 65 grids in y direction, and 1 block in z direction. One horizontal well is set in the middle, which can act as the injection well or production well during CO 2 Huff-n-Puff process, as shown in Figure 2. e basic data of the reservoir model are listed in Table 2. It should be noted that the reservoir heterogeneity is not included in this model.
As detailed fluid composition data for the Eagle Ford well is not available, the WINPROP module is used to establish a tight oil fluid model for the reservoir, using six pseudocomponents including CO 2 , N 2 , C 1 , C 2 -C 5 , C 6 -C 10 , and C 11+ . eir critical properties and interaction parameters are shown in Tables 3 and 4, respectively. e formation volume factor of 1.65rb/STB, bubble point pressure of 3446 psi, and gas-oil ratio of 1900 scf/STB were determined using the Peng-Robinson equation of state and flash calculations. Relative permeability curves are taken from Yu et al., which are obtained after history matching, as shown in

Result and Discussions
3.1. Single-Factor Analysis. In this section, we conducted a series of sensitivity analysis to determine the input parameters and their reasonable ranges in the CO 2 Huff-n-Puff process. Bottom-hole pressure, CO 2 diffusion coefficient, primary depletion time, CO 2 injection time, CO 2 injection rate, soaking time, production time per cycle, and cycle number are initially selected as the influencing parameters, while cumulative oil production is set as the response value. We only change one factor at one time while keep other parameters as constant. Detailed results are illustrated in the following sections.
Additionally, a base case of CO 2 Huff-n-Puff process is set with bottom-hole pressure of 1500 psi, diffusion coefficient of 0.01 cm 2 /s, and primary depletion time of 2 years. In each cycle, injection time of 4 months, CO 2 injection rate of 3000 MSCF/D, soaking time of 1 month, and production time per cycle of 2 years are set. e well will experience 5 cycles during the total production time. Figure 4 shows the difference between depletion and CO 2 stimulation under different bottom-hole pressure of 1000 psi, 1500 psi, 2000 psi, and 2500 psi, respectively. Results indicate that the production gradually decreases in the case of depletion drive, and the deviation between the CO 2 Huff-n-Puff process and depletion derive is more significant with the increase of Geofluids 3

Bottomhole Pressure.
bottomhole pressure. As shown in Figure 4, the production of CO 2 Huff-n-Puff process reaches the maximum when the bottom-hole pressure is 1500 psi, after that, it gradually decreases with the continuous increase of bottom-hole pressure. Since bottomhole pressure has such a great influence on the production, the bottomhole pressure is set as 1500 psi in the following simulation, which could obtain the highest cumulative oil production. Figure 5, we set the gas injection time as 2 months, 4 months, and 6 months, respectively. It can be seen that cumulative oil production increases with the increase of the injection time. However, the increment with injection time increasing from 2 to 4 months is larger than that from 4 to 6 months. erefore, the injection time is determined as 4 months for the following design, and the minimum and maximum levels are set as 2 months and 6 months, respectively. Figure 6 shows the improvement of oil production with the CO 2 Huff-n-Puff process with CO 2 diffusion coefficients of 0.0001 cm 2 /s, 0.001 cm 2 /s, and 0.01 cm 2 /s, respectively. As shown, the CO 2 diffusion coefficient has noticeable impact on the production; additionally, for this specific block, only when the CO 2 diffusion coefficient is 0.01 cm 2 /s, CO 2 Huff-n-Puff process is more effective than depletion production. e CO 2 diffusion coefficient is finally set as a constant value of 0.01 cm 2 /s, not considered in the following design.

Primary Depletion Time.
Primary depletion time is another important factor affecting the performance of CO 2 Huff-n-Puff process. As shown in Figure 7, the depletion production time is set as 1 year, 2 years, and 3 years, respectively. We compared the production of depletion drive and the CO 2 Huff-n-Puff process under different depletion production time of 1 year, 2 years, and 3 years. Results illustrate that the depletion production time has a little influence on the increment of the cumulative oil production. However, the length of depletion production time will affect the economic costs. Considering the effect of increase cumulative oil production and economic benefits, the primary         Figure 8 compared the well performance with the gas injection rates of 1000 MSCF/D, 3000 MSCF/D, and 5000 MSCF/D, respectively. As shown, increase in the gas injection rate results in the higher cumulative oil production. When the injection rate increases from 3000 MSCF/D to 5000 MSCF/D, the increment with injection rate increasing from 1000 to 3000 MSCF/D is much larger than that from 3000 to 5000 MSCF/D. at is, the injection rate has a reasonable range. In the following design, the injection rate is finally set as 3000 MSCF/D. e minimum level is 1000 MSCF/D and the maximum level is 5000 MSCF/D. Figure 9, the difference between depletive production and CO 2 Huff-n-Puff process under soaking time of 1 month, 2 months, and 3 months indicates that soaking time has a minimal effect on the improvement of cumulative oil production. In the previous studies, soaking time has an effect within a certain range, and beyond the upper limit of the range, the effect becomes limited and can be ignored [28,30]. Finally, the soaking time was set as a fixed value of 1 month, which was not analyzed as an engineering factor in the following.

Production Time per Cycle.
Appropriate production time during each cycle is essential for pressure maintenance. As shown in Figure 10, we set the production time per cycle as 1 year, 2 years, and 3 years, respectively.
Comparing the production of depletive drive and CO 2 Huff-n-Puff process under different production time per cycle, it can be found that production time per cycle has a little effect on the increase in the cumulative oil production. However, the length of the production time per cycle will affect economic costs, which needs to be considered as engineering factor. When this factor is 2 years, the dependent variable reaches the maximum. As a result, the production time per cycle is set at 2 years, with a minimum level of 1 year and a maximum level of 3 years. Figure 11, cycle number is set as 3, 5, and 7, respectively. When the number of cycles increases from 3 to 5, the cumulative oil production increases, but it begins to decrease with the cycle number from 5 to 7. erefore, the cycle number has an appropriate range to maximize its effectiveness. For this specific block, the production reaches the highest at cycle number of 5. In the following design, the cycle number is set to 5, the minimum level is 3, and the maximum level is 7.

Box-Behnken Design.
With the single-factor analysis, primary depletion time, injection time, cycle number, production time per cycle, and injection rate are selected as the input variables in the Box-Behnken design. e input variables and their levels are shown in Table 5. "−1"and "1" represents the minimum and maximum of the variables. Design of Expert 12 was used to select Box-Behnken design module and input the prescreened single factor and its level range. A response surface analysis table was designed automatically by Design of Expert 12, and 46 experiments were obtained. After setting the Box-Behnken design models, 46 groups of reservoirs simulation models are totally run by CMG-GEM module, and the cumulative oil production and accumulated cash flow of each scheme are selected as the two output parameters. Since economic evaluation has its own formula, it is not included in the regression analysis, but it still needs to be taken into consideration in Box-Behnken design as a response value.
Regression analysis was conducted on the cumulative oil production to screen out the applicable model. Table 6 illustrates the fitness for different regression methods. As listed, the quadratic model is used for the regression analysis in this study, which provides higher correlation coefficient (R 2 ). After that, variance analysis (ANOVA) was applied to offer the numerical information for the parameters and their combination terms. Table 7 shows ANOVA results and the most important parameter in ANOVA results is p-value.
e p-value is described as the probability of accidental occurrence of dependent variable results, and if the p-value of a factor is less than 0.05, it indicates that this factor has an important effect on the developed model [45]. F-value is the ratio of two variances and it is always combined with p-value to determine the importance of each factor. Factors that have higher F-value result in more significant influence.
As shown in Table 7, the p-value of the model is less than 0.0001, indicating that the cumulative oil production has a very significant regression relationship with its respective variables, which can better reflect the real linear relationship. In terms of quadratic coefficients, p-value of production time per cycle is higher than 0.05, meaning that it is less important for the cumulative oil production. e value of p for other factors are less than 0.0001, indicating that they have the greatest impacts on the cumulative oil production. Figure 12 also presents the significance of each influential factor. e Y value of 0 is taken as the reference; the larger the deviation of each curve is, the greater the influence of the corresponding factor on the response value will be. As shown, the injection rate is the most significant, which has the largest deviation, while the production time per cycle with the smallest deviation is the least significant. Hence, the order of influence of various factors on the cumulative oil production of CO 2 Huff-n-Puff process is injection rate-> injection time > cycle number > turn time > production time per cycle. Additionally, the p-values of interaction parameters AC, AD, BC, CD, CE, DE, A 2 , B 2 , C 2 , D 2 , and E 2 are also less than 0.05, which demonstrates that they have extremely significant influence on the cumulative oil production difference. Figure 13(a) shows the normal probability plots of residuals for the cumulative oil production, which is a diagnostic tool to assess the validity of the model. It also can be    evaluated by Figure 13(b), which shows a scatter plot of real and model predicted values. e coefficient of determination R 2 is utilized to check the proposed model. R 2 value closing to 1 means that the compatibility between experiments and the results obtained from the proposed model is high.
us, a reasonable regression equation for the influence of the input variables on the response value is generated, as shown in (4). In this equation, y is the value of cumulative oil production.

Evaluation of Influencing Factors and Prediction of the Optimal Scheme.
e graph of the RSM method is a threedimensional spatial graph and a contour map on a twodimensional plane composed of factors corresponding to the specific response value (cumulative oil production). is typical graph can intuitively reflect the influence of each factor on the response value. e shape of the contour map reflects the strength of the interaction. A steeper response surface curve indicates a stronger interaction. Figure 14 shows the response surface diagram of several factors with significant interaction judged by the F-value (as higher as possible) in Table 7. Cumulative oil production is set as the response value. As shown in Figure 14, it can be seen that the surface has a partial convex arch, which corresponds to the innermost circle region of the lower two-dimensional contour graph. If the shape of the innermost line deviates to the circle line, it means the obvious interaction. We can see the inner circle of Figure 14(b) is more elliptical meaning that the interaction of factor C and D is more obvious, which is also confirmed by their maximum F-value in Table 7.    Furthermore, the denser the contour line in a certain direction, the greater the influence of the factor on the response value. As shown in Figure 14(c), the number of contour lines in the direction of factor E is 5, which is more than that in the direction C, indicating that factor E has a greater influence on the response value in the interaction between factor C and E. e collective effect of injection time and cycle number on the cumulative oil production are presented in Figure 14(a). e cumulative oil production increases at the beginning and then decreases with the rising injection time and cycle number. e maximum cumulative oil production is found as 618 Mbbl. Figure 14(b) provides the interaction of cycle number and production time per cycle. Increasing the cycle number and production time per cycle will lead to a reduction in the cumulative oil production. Moreover, the cumulative oil production rises as expected with the incremental injection rate, and it also can be seen that surface color change is rich, which means these factors has a good interaction, as seen in Figures 14(c) and 14(d).
Based on the data obtained from the single factor experiment, the combination of factors and levels used in the design of the Box-Behnken model was determined. After a comprehensive analysis of the two response values of cumulative oil production and cumulative cash flow, the following optimization conditions can be obtained: primary depletion time is 2.37 years, injection time is 3.4 months, cycle number is 3, production time per cycle is 2.2 years, and injection rate is 5000 MSCF/D. For this optimal design, the forecast of cumulative oil production is 614.663 Mbbl, and the cumulative cash flow model forecast is 4476.86 thousand dollars.
We also run the simulation with these optimized parameters, and the simulation results of the cumulative oil production is 600.124 Mbbl, with an error of 3%; the cumulative cash flow is 4234.626 thousand dollars, with an error of 5%. Hence, the prediction model is reliable and reasonable in the design of CO 2 Huff-n-Puff operation.

Conclusion
In this study, RSM is applied to analyze the influential factors and optimize the operating parameters of CO 2 Huff-n-Puff process in the Eagle Ford shale oil reservoirs. Box-Behnken design is used to determine the input variables and ANOVA is utilized to analyze the compatibility of the proposed model and the importance of operating factors. According to the experimental results, the following conclusions can be drawn for this typical block: (1) With the single-factor analysis, the effect of soaking time on cumulative oil production is negligible. (2) Corresponding to the results taken from the optimization, the order of influential factors in this specific CO 2 Huff-n-Puff process is injection rate-> injection time > cycle number > primary depletion time > production time per cycle. (3) Based on the RSM results, the optimal operating parameters for the CO 2 Huff-n-Puff process are determined as primary depletion time is 2.37 years, injection time is 3.4 months, cycle number is 3, product time per cycle is 2.2 years, and injection rate is 5000 MSCF/D. (4) e little deviation of R 2 between the predicted and adjusted model illustrates the accuracy of the proposed model. Adeq Precision (93.1149) which is larger than 4 further indicates the good reliability. (5) is study provides a comprehensive tool to efficiently optimize the CO 2 Huff-n-Puff process in the unconventional reservoirs, and it also can successfully determine the best operating parameters to achieve the optimal solution.

Data Availability
All the data supporting the results of this study have been presented in the paper.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.

Authors' Contributions
Yinqing Wang conceptualized the study, developed the methodology, and wrote the original draft. Weiwei Xie suprevised the study. Jinghong Hu supervised the study and developed the methodology. Yuan Zhang performed projection administration and reviewed and edited the article.