Utilization of nature-inspired algorithms for gas condensate reservoir optimization

Hydrocarbon production from a gas condensate reservoir has a major concern to the petroleum industry around the world. Due to pressure decline in a reservoir, liquid drops out of solution and cause a significant reduction of well productivity. As a result, huge amount of valuable oil stays underground. The proper production management for a gas condensate reservoir minimizes or eliminates the blocking problems. Optimization of condensate recovery can be achieved through appropriate selection of production and injection rates. However, determination of the optimal rates is a complex problem, involving many factors including geological uncertainty, rock petrophysical characteristics, fluid properties, economic costs, and technical ability. In this research, nature-inspired optimization algorithms are employed to evaluate the potential of production boost. Numerical testing and comparative study revealed that investigated natural-inspired algorithms outperform widely reported in the literature optimization methods and provide a higher quality solutions. The primary concern of this research is to develop an optimal gas condensate production strategy in terms of economic efficiency by utilization of nature-inspired optimization algorithm and providing their ability to generate better quality solution than engineering approach as well as gradient optimization.


Introduction
Oil and gas reservoirs can be divided into five categories based on fluid types (McCain 1990): black oil-consists of heavy, non-volatile hydrocarbons; volatile oil-contains fewer heavy molecules and more intermediate components (C 2 H 6 -C 6 H 14 ); condensate-condensate gas is very similar to volatile oils in terms of produced oil; however, the reservoir temperature of the condensate gas reservoir is greater than the critical temperature of the fluid; and volatile oil is liquefied at original reservoir pressure and temperature, and a condensate gas is in gaseous state; wet gas-natural gas that is rich in significant heavy hydrocarbons such as propane, butane and other liquid hydrocarbons; dry gasnatural gas that appears in the absence of condensate and liquid hydrocarbons or gas from which condensable hydro-Communicated by V. Loia.
B Damian Janiga janiga@agh.edu.pl 1 Faculty of Drilling, Oil and Gas, AGH University of Science and Technology,al. Mickiewicza,Poland carbons were removed. Hydrocarbons produced from gas condensate reservoir have a major concern to the petroleum industry around the world. Gas condensate reservoirs are significantly different from conventional reservoirs in terms of thermodynamic properties and flow characteristics. Reservoir performance prediction and optimization require careful analysis and project management (Fasesan et al. 2003;Bozorgzadeh et al. 2006). During production process, reservoir pressure falls below dew point and a liquid phase rich in heavy ends drops out of solution in near wellbore zone. As a result, volume of condensate becomes immobile and does not take part in the flow stage (App et al. 2007). The phase behavior of gas condensate reservoir type is strongly dependent on phase envelope and thermodynamic conditions of reservoirs fluid. A typical P-T (pressure-temperature) phase diagram ( Fig. 1) has the critical temperature lower than the reservoir temperature; however, the cricondentherm temperature is larger than the reservoir temperature. The microscopic condensate droplets tend to be trapped in pores causing additional pressure drop in near wellbore zone (Fig. 2). In the zones far from well, the liquid is immobile due to low saturation and capillary forces. Moving through the well center liquid saturation is higher than critical value (S or ) and mobility increases. Condensate blockage problems are registered in the numerous hydrocarbon fields. Arun field operated by Mobil lost over 50% of its initial productivity due to condensations issues (Afidick et al. 1994; Barnum et al. 1995). Shell reported 67% productivity loss (Smits et al. 2001) because of the same aspect.
Various remedial actions to alleviate the pressure drop caused by gas condensate blockage through years were proposed. The most common activity is gas injection to maintain reservoir pressure above dew point line. Dry gas injections were implemented in Sleipnerost Ty (Huerta Quinones et al. 2010) and Margham condensate gas field in Dubai (Ernster et al. 1988). In the field practice, dry gas (consisting of mostly methane) has been often used as injection gas; nevertheless, some researchers propose more widely available and cheaper substitute in the form of CO 2 , nitrogen, air, C 4 H 10 + C 5 H 12 injection, rich gas injection (Siregar et al. 1992;Sanger et al. 1998;Sadooni and Zonnouri 2015;Xiong et al. 2007;Lewis and Couples 1993;Asar et al. 1988). However, due to mixing of nitrogen, CO 2 and others with reservoir fluid, significant liquid dropout occurs, even though the pressure is maintained above the dew point value of the gas. Abel et al. (1970) introduced two schemes of gas injection: full and partial pressure maintenance. The most widely used technique to improve gas condensate reservoir performance is gas cycling (Kleinsteiber et al. 1983). The objective of the cycling in gas condensate reservoir is to maintain pressure close to the dew point and recover more wet gas with a maximum condensate yield. The problem is that during the life of the reservoir, breakthrough appears and the condensate yield decreases in some wells located close to the gas injectors.
Efficiency prediction of gas condensate field development strategy is usually carried out by enumeration of certain amount of possibilities (trial and error method). Over the past few years, several investigators proposed approaches to evaluate gas condensate reservoirs and optimize condensate recovery. Kamari et al. (2017) researched the influence of reservoir pore size on condensate PVT properties and concluded that increase in pore sizes and reservoir temperature increases dew point pressure. Sakhaei et al. (2017) developed the study of enhancement condensate production through wettability alteration. Su et al. (2017) performed experimental investigation and modeling simulation. Results revealed that the gas injection technique is more effective than water flooding method. Zhou et al. (2016) evaluated the critical gas rate necessary to prevent liquid loading in deep condensate gas wells. Kalugin et al. (2015) proposed to use the steepest descent method with fractional step to determine the optimal production rates which can increase volume of recovered gas condensate. For clarity, only depletion mode was considered. Ghorbani et al. (2017) suggested to use the firefly optimization algorithm for prediction of gas flow rates from gas condensate reservoirs through wellhead chokes. Kaydani et al. (2016) established the multi-gene genetic programming model for dew point determination with sensitivity analysis. Mahdiyar and Jamiolahmady (2014) developed the new methodology for the optimization of fracture design in gas condensate fields.
Because of reservoir heterogeneity, the configuration of exposure well and different factors, it is dubious that solution found during scenario evaluation process is really optimal. For that reason, several nature-inspired optimization algorithms, which are metaheuristic models, offer high-quality solution for the given optimization problem, which are applied to develop optimal gas condensate production strategy. Nature-inspired algorithms mimics nature mechanisms and principles which increased the evolutionary adaptability of individuals to the environment, enabling survival and reproduction. For years, humans have utilized the guidance of nature in finding the most appropriate solution for complex engineering problems (Yazdani and Jolai 2016). However, due to principle of "No free lunch theorems" there is no optimization technique for solving all optimization problems (Wolpert and Macready 1997). Different algorithms possess capabilities for solving different types of optimizing issues, since it is difficult to predict the best algorithm for every optimization problem in terms of objection function value and computational cost.
Motivated by the above consideration, in this paper, particle swarm optimization (PSO), genetic algorithm (GA), grasshopper optimization algorithm (GOA), moth-flame optimization (MFO) colony and Grey Wolf optimization (GWO) were in detail analyzed to find the best solution which influence and determine the right gas condensate field management. PSO and GA were chosen because of large amount of real case optimization studies, not only for petroleum engineering field. The GOA, MFO and GWO are relatively new algorithms and were chosen due to their benchmark results comparing with different methods (Heidari and Pahlavani 2017;Sulaiman et al. 2015;Mei et al. 2017;Medjahed et al. 2016;Kamboj 2016;de Moura Oliveira et al. 2016). Obtained results were compared with expert knowledge approach for the gas condensate recovery process and with adjoint gradient optimization. From the experimental results, the investigated natural-inspired algorithms outperform knowledge-the base gradient optimization approach.
The remainder of this paper is organized as follows: Section 2 briefly presents the basic concept of gas condensate field production optimization, work flowchart, objective function definition and description of the utilized algorithms; in Sect. 3, experimental results and performance comparison are given; in Sect. 4, conclusions and future work are described.

Methodology
In this section, the optimization workflow process is presented in detail. Moreover, synthetic numerical model and involved optimization algorithms are also described.

Optimization problem definition
Hydrocarbon flow through porous media can be written as a set of non-linear functions ( f ()) representing conservative relations (1) (Ramirez 1987): where x is reservoir state vector and u is control variable. Flow equations are solved in volume-time domain (Ω × t) of problem where start and end of simulation are t 0 = 0 and t f , respectively. Optimization gas condensate reservoir depletion problem can be written in continuous form as (2): where y M is solution quality after performing M iterations of optimization loop and J () is an objective function depended on the state and control vector in time t. Detailed objective function takes the form of profit coming from condensate sales reduced by cost of gas injection (3): Cumulative condensate production and gas injection rates from all wells are described as Q P and Q I , where q I and q P are individual rates for N p production and N I injection wells, respectively (4): Equation (1) was solved numerically by ECLIPSE software; therefore, continuous volume-time domain was divided into control steps. Following this, simulation time was split into N t discrete time steps t. Elements of control vector are: individual gas production (w) and injection rates (q I ), which for every discrete time step cover simulation time domain (5): Condensate withdraw rate (w(t)) is a function of gas production (q P (t)), reservoir state (x(t)), petrophysical characteristics and fluid properties as well as fluid-rock interaction (h(x(t))) (6): Individual gas production and injection rates were bound by their lower and upper limit: 0 Mscf/day ≤ w ≤ 5000 Mscf/day and 0 Mscf/day ≤ q I ≤ 2000 Mscf/day, respectively.
Given problem of gas condensate flow in porous media was evaluated by numerical simulation of the synthetic field. Reservoir parameters are in detail described in the below subsection.

Numerical model of gas condensate reservoir
A numerical 13-well model was constructed using gas condensate fluid and reservoir properties based on SPE 10 benchmark model. For this purpose, Schlumberger ECLIPSE E300 compositional reservoir simulator (version 2016.2) was utilized. This simulator was successfully applied for the compositional problems by Maleki et al. (2012) and Nasriani et al. (2015). Fully implicit solution is selected for solving the mathematical models and governing equations for objective problem.
Heterogeneous reservoir simulation model ( Fig. 3) contains six layers and has a three-dimensional system and consists of 30 × 30 × 6 grid cell arrangement, corresponding to the x, y and z direction. Simulation model covers area of 9842.52 ft × 9842.52 ft × 92.42 ft. Initial reservoir pressure (324 bar) is above dew point (dew pressure P d = 308 bar) at temperature of 355 K. PVT fluid properties were studied using three-parameter Peng-Robinson equation of state. Production wells exposure layers 3-6, while injection wells are open in 1-2 layers. Wells distribution is presented in Fig. 4. Simulation time was set for 15 years with a half year control step. Minimum bottom hole pressure for production wells was 68.9 bar. Production wells were controlled by gas production rate (in surface conditions). Pressure maintenance was performed by gas injection with maximum bottom hole pressure equal 413 bar. Reservoir temperature was expected to be constant during gas injection (no thermal effects).

Optimization methods
To evaluate best gas condensate production strategy, comparison between base case scenario (expert knowledge), gradient optimization as well as nature-inspired algorithm was carried out. The following section is introduced in order to let readers to better understand solution methodology.

Expert knowledge approach for gas condensate recovery process
Gas condensate recovery process should be investigated for every reservoir independently. The petrophysical and fluids properties are unique and affect development strategy. In case of direct solution copied from similar field or closely located, some engineering problems may occur. As indicated by Coats et al. (1985), gas cycling with full pressure maintenance is most promising technique to enhance condensate production; therefore, base simulation was performed with the following assumptions: -individual injection rate for 9 injection wells was set to 1000 Mscf/day, -individual withdraw rate for 4 production wells was set to 2250 Mscf/day.
Injection and production rates were kept at the constant level during recovery time. The conventional simulation approach is referred as ENG on the result graph.

Gradient optimization method
The problem is to maximize the objective function over time period T (from t 0 to t f ); thus, the main constraints are described by the system equations (1). The optimization variable (u) is referred as controlling input in terms of production and injection gas rates. Equations set (1) can also be expressed as (Ramirez 1987;Brouwer and Jansen 2004): By introducing Lagrange multiplier λ, the problem can be transformed into an equivalent unconstrained problem (Ramirez 1987;Brouwer and Jansen 2004): Derivative of objection function was used to determine the search direction for the improved vector control in every iteration: where d m = −d J /du and ω is step size multiplier. Given search direction, d, the standard line search algorithm is to determine some step size multiplier, ω, which maximizes the objective function value Passino (2005): As a initial solution for gradient optimization u m=0 engineering approach with adequate assumption was used. Solution gradient was obtained by internal algorithm of the utilized reservoir numerical simulator.

Nature-inspired algorithms
Five nature-inspired algorithms (GOA, MFO, GWO, PSO and GA) were applied to solve the optimal gas condensate production management and algorithms were compared each other in terms of solution quality and convergence speed. For every algorithms solution, population consisted of N = 20 individuals. All algorithms were terminated after 30 generation. Due to used numerical simulation for the objective function evaluation, each objection function called demand complete run; therefore, convergence speed is the crucial parameter. Efficiency comparison of the same initial population of solutions was implemented to keep reliability of nature-inspired algorithms. In the worldwide literature database, there are only few papers about implementation nature-inspired algorithms to solve optimization problem connected with oil and gas industry. Best to the authors knowledge, this is the first application of MFO and GOA in reservoir engineering. Widely reported PSO and GA can be found in the related works (Eberhart and Kennedy 1995;Song et al. 2016).
Grasshopper optimization algorithm (GOA) was developed by Saremi et al. (2017). Algorithm mimics social behavior of insects, where individuals join the largest swarm (Simpson et al. 1999). Life cycle of grasshopper consists of the three phases: eggs, nymph and adulthood. The unique aspect of every form can be found. The main characteristics of grasshoppers in the larval phase are slow movement and small steps. In contrast, long range and abrupt movement are essential features of the swarm during adulthood. In mathematical position of ith, GOA particles form can be written as (11) Saremi et al. (2017): where ub d and lb d are upper and lower bounds in the Dth-dimensional search space, s() function indicates social forces, c is a decreasing coefficient to tighten the comfort zone, the repulsion zone, and the attraction zone andT is the best solution found. Parameter c is a tuning operator describing:(I) exploration balances and exploitation balances of the entire swarm around the target, (II) the attraction, the comfort and repulsion zones between grasshoppers. Moth-Flame Optimization Algorithm was established by Mirjalili (2015), and inspiration was drawn from moth navigation skill during night, which is called transverse orientation. In this method of orientation, moths maintain a fixed angle on a distant source of light with respect to the Moon. It is very effective mechanism for traveling long distances in a straight path (Gaston et al. 2013). Despite transverse orientation effectiveness, it can be observed that moths fly spirally around the light because of artificial light tracking. Mathematical form assumes positions of each moth (particle) in search space, flame position (source of light) and value of objective function for moths and flame. It has to be noted that both moths and flame are solutions. The difference between them is the way that investigator treats and updates them in each iteration. Moths are the actual searching agents that move around the search space, whereas flame is the best position of moths to be reached. In order to mathematically model moth-flame optimization, position of every moth has to be updated with respect to flame using the following spiral equation (12) (Mirjalili 2015): where D i indicates distance between ith moth for jth flame, b is a spiral constant, and t is a random number in range [−1, 1]. Mirjalili inspired by gray wolves behavior developed Gray Wolf Optimization (GWO) algorithm Mirjalili et al. (2014). The GWO is an evolutionary computing-based technique which simulates the hierarchy of gray wolves pack. Detailed description of haunting mechanism was described by Muro et al. (2011). In mathematical convention counting, surrounding, and assailing the prey (updated solution) can be written as (13) Mirjalili et al. (2014): where X i is a vector position of wolves, X p is a position of prey, A is an attracting parameter, and C is a random vector introducing exploration element. In mathematical model, the social hierarchy of wolves when designing GWO the best solution is called alpha (α). Consequently, the second and third the best solution are named beta (β) and delta (δ), respectively. The rest of the solutions are assumed to be omega (ω).

Optimization framework
To find the optimal gas condensate reservoir developing plan, every possible solution u from N -elements population is evaluated by reservoir simulation. Control variables were transferred to simulator by functional keywords. The results from simulation (total condensate production and total required cumulative lean gas injection) were used to calculation of objection function value (3). Optimization flowchart for population-based algorithms is presented in Fig. 5. Optimization framework was implemented in MATLAB software in conjunction with Schlumberger Eclipse Simulator which is industry reference reservoir simulation tool. To keep reliability of nature-inspired algorithm, efficiency comparison of the same initial population was implemented. For every of proposed nature-inspired algorithms, elitism was introduced, what allowed to pass from the best solution to the next generation.

Performance of nature-inspired algorithms
Aforementioned algorithms were examined in terms of convergence speed and exploration ability. Convergence (σ ) in mth generation for populated-based algorithm can be expressed as (14) Fig. 5 Optimization workflow for population-based optimization algorithm, due to elitism best solution can be pick from last population where y * is a global optimum solution, y m max is the best solution found in population in mth generation and y m=0 max is the best solution in initial population. Population diversity measures the distribution of individuals. The diversity changing rate is a way to monitor the degree of convergence/divergence. Diversity of solution (Ψ ) in mth generation can be written as (15) (Ni and Deng 2014) : whereȳ m is average solution for mth generation.

Results and discussion
This section provides experiment results with discussions. Results which are put forward in this section are straightforward optimization indicator in terms of objective function value and parameters which have significant impact on solution. Effect of well rate control was carried out by comparison of solution matrix. Exploration and exploitation abilities of algorithms were investigated in conjunction with statistical analysis of obtained results.

Optimization results
As a result of engineering approach simulation 172.55 MM$ of objection function value was obtained. By including gradient optimization, proper production and injection policy can increase solution quality to 178.14 MM$.  Fig. 6. Comparison of the obtained results indicates that all of proposed nature-inspired optimization algorithms outperform engineering approach of gas condensate recovery problem. Summary results of objective function values are presented in Fig. 7. Genetic algorithm delivered better results by 4.15%, GOA 2.90%, GWO 7.17%, MFO 23.88% and PSO 26.78% in respect to engineering knowledge, thus providing the necessity of optimization work. Objective function value increased from 172.55 MM$ to 178.14 MM$ (+ 3.24%) for base case study using gradient optimization. The results revealed that PSO generates the best project profitability. The results are in line with previous investigation where the authors employed PSO and four different nature-inspired optimization algorithms to find optimal polymer injection strategy to enhance oil recovery. We indicated that of all the five implemented nature-inspired algorithms are comparable in magnitude and among them, PSO allowed to find the highest net present value (NPV) (Janiga et al. 2017).

Effect of wells rate control on optimization process
Investigation of control parameters (injection and production rates) was carried out in three generation points and is presented in Fig. 8 (after second and third). Diversity of objective function value (− 1.71 MM$ in second generation between MFO-PSO and − 6.61 MM$ in third generation) has reference in control matrix and is illustrated in Fig. 9, where in case of injection, 49.62% of controls changed in positive direction which increase rates are between MFO and PSO. Similar trend is observed for production control matrix where 48.33% rates has positive sign change. Thus, relative improvement in objection function value is related with increase average injection and production rates. Change of Fig. 12 Dependencies between injection and production rate including function fitness value control vector between second and third for PSO and MFO is presented in Fig. 10, revealing which part of solution matrix u has significant impact on improving solution. Increase of objective function (+ 6.63%-PSO, + 3.74-MFO) between following generations is associated with injection and production control. Generally, enhanced injection and production rates give better results. Taking into consideration both MFO and PSO algorithms, which are preferable for the proposed optimization problem, they present relatively different optimization process. As was assumed, decision variables are production and injection rates which are closely related to cumulative production and injection. Analysis of optimization process for MFO and PSO algorithm is presented in Fig. 11 proving differences. PSO algorithm at initial stage, where exploration was a major aspect, frequently changed possible solution. After that more schematic select of variables was observed-turnaround from exploration to exploitation. At the final point of optimization run, over 200 and 300 trials of around 70,000 MMscf of total production and 68,000 MMscf of total gas injection were For MFO algorithm, balance between exploration and exploitation is harder to be noticed. Histogram of possible solution of cumulative production and injection is more uniform and unimodal. Analysis of relationship between control data (injection and production cumulative volume) in regard to function fitness value reveals linear dependency, especially for total gas injection-project profitability (R 2 = 0.9351 for MFO and R 2 = 0.9832 for PSO), and detailed analysis is presented in Fig. 12. Optimal control matrix (u) for PSO and gradient algorithm is presented in Fig. 13 and consists of production and injection rates for all wells for every discrete time step. After high-rate production period, every production well rate decreases and injection oscillations rates are periodic.

Condensate production analysis
Gas condensate reservoir production optimization allowed to recover from 4.54 MMstb (GOA) to 5.61 MMstb (PSO). Additional hydrocarbon production can raise reservoir management profitability. The right managing policy in terms of production and injection rates eliminated condensate dropout in near well zone. Based on PSO results, condensate saturation in bottom part of reservoir was reduced from (average) 21.22 to 11.20%. Condensate saturation is presented in Fig. 14. Additionally, more uniform drainage zone can be observed, which resulted in increment of reservoir recovery factor.

Algorithms efficiency
Going through analysis of algorithm efficiency (Fig. 15), GWO offers higher convergence ratio at initial stage and exceeds other algorithms. GOA gives monotonic convergence increase. Search ability of solution space can be expressed as population diversity (Ψ ). The best exploration ability appears for MFO algorithm. For GOA algorithm, some premature convergence is observed, where population diversity decreases to zero at final optimization stage and this fact is referred as algorithm stuck in local minimum. Other algorithms are in the right balance between exploration and exploitation ability. The algorithms created the unique solutions (far to each other in search space) during optimization

Routine and advanced statistical analyses
Based on carried iterations, for every fifth generation minimum, maximum, average, and standard deviations of objective function are listed in Table 1. In order to provide better analysis of the results, the following hypothesis are investigated to test which algorithm (alg. 1) statistically outperform other algorithm (alg. 2) in terms of convergence and spacing metric according to the single tail t test with α = 0.05: If the p value of the statistical test satisfies condition p > α, then it indicates that the hypothesis H 0 is right and there is no significant difference; otherwise, e.g., algorithm 1 is better than other algorithms. Results of test are listed in Table 2. The Wilcoxon signed ranks test was used to determine similarities between two population (Derrac et al. 2011). Nonparametric procedure was employed to examine hypothesis, which involved designing two samples. This research is analogous to the paired t test in nonparametric statistical procedure; therefore, it is a pairwise test that aims to detect significant differences between two samples, i.e., the behavior difference between two algorithms. Evaluation of statistical significance difference between every population was carried out by tests following hypothesis with α = 0.05: If statistical test was less than or equal to the value of the distribution of the Wilcoxon for n degrees of freedom, the null hypothesis of equality of medians was rejected. That fact is means that some algorithm outperforms the other one, with the p value associated. The Wilcoxon test was used to evaluate all investigated algorithms with full range of generations (gen:i, where i = 1, . . . , N ). Test results are presented in Fig. 16, where similar populations are marked. Based on statistical test, GA and GWO have comparable population in over 37% of generation (omitting initial population) and it can be stated algorithms have not any exploration abilities. Detailed results are listed in Table 3.

Conclusions and future work
Imitation of nature for years has been the inspiration to solve immensity optimization problems. In recent times, natureinspired algorithms become the standard in researching complex engineering and non-trivial problems. Notwithstanding, individual characteristic of algorithms leads to disparate solutions, which can be hard for the prediction. Management of gas condensate reservoirs is complicated due to their phase and flow behavior; furthermore, geological and technological aspects affect possible solutions. Assumption of mentioned conditions can be aggregated in the form of optimization problem. In this research, the effort was put to combine GA, GOA, GWO, MFO and PSO algorithms with full-scale simulation model to determine optimal reservoir policy. The results of the conducted numerical experiment indicate that the investigated nature-inspired algorithms have various degrees of the robustness. PSO and MFO algorithm engender better results in case of solution quality and convergence speed. Application of nature-inspired optimization method (PSO algorithm) can improve condensate production up to 42.36% in conjunction with reduction of liquid genera- tion from 21.22% to 11.20%. Project profitability increased about 52.78% based on production and injection well rates and change is connected with low operational cost. Further work will investigate other features of utilized algorithms in different real oil and gas field case stud-  ies. Moreover, improved convergence ratio and optimization results will be analyzed. In this study reservoir temperature was assumed to be constant during gas injection process; therefore, in the future work, thermal effects will be taken into consideration.