Optimization of Horizontal Well Direction and Length Considering Geomechanics Properties and Drainage Area Using Genetic Algorithm in A Gas Field

To maximize a horizontal well production, we need to determine the optimum direction and horizontal well length that maximizes the gas field recovery for a certain constant flow rate called by plateau rate. This problem is conventionally solved by using a reservoir simulation model and trial and error procedure that consumes considerably a lot of time and efforts. This study uses a random search method, i.e., Genetic Algorithm (GA), to solve this optimization problem and it very much eases to find the best well location with less time and efforts consumed. Along the general technique in directing a horizontal well towards the least principal stress of rocks, this study considers the geomechanics effects that influence the gas production performance. And also, the drainage area of horizontal well will be considered in this study to obtain the optimum horizontal well direction and length. In order to do this, a new proposed objective function for the GA has been constructed based on basic reservoir properties (i.e., porosity, permeability and gas saturation) and geomechanics properties (i.e., Young’s modulus and Poisson’s ratio). The results of the proposed method are validated using a reservoir model and economics evaluation. It may be concluded that the applying GA, with the appropriate objective function, can give accurate and faster results compared with the trial and error method using reservoir simulator, technically and economically, and also the proposed method is able to reduce the amount of works considerably time.


Introduction
Until now, the determination of the optimum horizontal well direction and length has usually been conducted by using a conventional trial-and-error method, when placing coordinates at a geological or reservoir model and then running with a simulator to obtain the gas recovery; this process requires experience, considerable time, and a relatively high cost. Especially when addressing oil or gas fields that have relatively large size, we are faced with the possibility of a remarkably large number of solutions (locations of wells), which causes the conventional methods to become ineffective and inefficient.
Generally, the exhausted trial and error procedures have to be run to obtain the optimal production scenario in forecasting future reservoir performance. It spends quite a lot of time and cost of applying some case sensitivity of the development scenario with trial and error procedure.
Basically, longer horizontal well section will give higher gas recovery. However, after certain length, increasing the horizontal well length causes the plateau time period and recovery factor of the reservoir decreased. Hence, there are some considered parameters in applying the horizontal well length, such as: (1) Reservoir boundary. Reservoir production performance will be lower in the area near the reservoir boundary. This phenomenon happened because of the drainage radius achieved by the well near the reservoir boundary is narrower. Thus, increasing the horizontal well length too close to the reservoir boundary will not give higher gas recovery. (2) Low reservoir properties quality. If the horizontal well section has reached the bad reservoir properties zone, it will reduce the gas production. Furthermore, if the bad zones penetrated by the horizontal section are the high pressure zone, then the field will experience the significant decreasing in production.
Logically, one would like to create fractures which are perpendicular to the horizontal wells in the vertical plane. The purpose is to enhance vertical permeability and to increase the productivity of horizontal well. The direction of the fracture orienting from the horizontal well is the same as that for vertical wells, i.e., parallel to the plane of minimum principle stress (Joshi, 1990). Therefore, if a horizontal well is drilled along the low principle stress direction, then the stimulated fractures will be perpendicular to the horizontal well. However, if the horizontal well is also drilled along the maximum principle stress direction, then the stimulated fractures will be parallel to the horizontal well. Therefore, if the horizontal well is to be stimulated, it is important to consider the local directions (Joshi, 1990). Horizontal well drainage area is one of the important variables needed to be determined in the production forecast of horizontal wells. By knowing the drainage area, we can determine the wells spacing required so that there will be no production interference between the wells. The drainage area can be calculated by using either of the methodologies proposed by Joshi in 1990 and1991. These two methods consider only the geometry of the drained surface, and no other factors such as reductions in productivity that triggers pressure loss through the productive section (Saavedra et.al, 2001).
As a directed random search technique, Genetic Algorithm (GA) might be offered as an approach to select best production scenario without having a deal with the optimization parameter directly and to replace the conventional trial and error procedure efficiently. A genetic algorithm model is developed in order to execute effectively reservoir simulation results to find an optimal development scenario in a significantly reduced time. GA technique can simplify the optimization problem because it is not directly working with much information about the optimization problem. In this study, GA as a directed random search technique is applied to select an optimum development scenario for X gas field. The X gas field is developed from one horizontal well with a crucial parameter of horizontal well direction and length that influence much on the plateau time production, as well as on the recovery factor. GA has a role to replace the conventional reservoir simulation technique of trial and error procedure in finding the best horizontal well direction and length for maximizing plateau time production in X gas field, technically and economically.
In previous study results , we can see that the best direction of horizontal well is 315 0 with the optimum length 700 meter which yields 7.173 years of plateau time. As shown in the simulation results in previous study , sensitivity in horizontal well direction and length gives different production performance from the reservoir. Basically, longer horizontal well section will give higher gas recovery. However, after certain length, increasing the horizontal well length causes the plateau time period and recovery factor of the reservoir decreased. Hence, there are some consideration parameters in applying the horizontal well length: (1) Reservoir boundary: Reservoir production performance will be lower in the area near the reservoir boundary. This phenomenon happened because of the drainage radius achieved by the well near the reservoir boundary is narrower. Thus, increasing the horizontal well length too close to the reservoir boundary will not give higher gas recovery.
(2) Low quality of the reservoir properties: If the horizontal well section has reached the bad reservoir properties zone, it will reduce the gas production. Furthermore, if the bad zones penetrated by the horizontal section are the high pressure zone then the field will experience the significant decreasing in production.

Method and Materials
X gas field that will be used in this study is located about 7.5 km east of Java in the part of South Madura Basin with area ± 1.6 km 2 . South Madura Basin is located around East Java Sea structural area. The fluid sample was taken during Drill Stem Test (DST) that performed in the gas zone of X-1 exploration well, which it shows that the methane composition was 95.83%. Thus, the X field reservoir fluid can be categorized as dry gas.
The reservoir thickness of X field is 443 ft. The Gas Water Contact (GWC) is located at depth -1362 mSS based on log analysis and RFT plot. The GWC is used to estimate the gas reserves in X gas field. The volumetric calculation shows the Initial Gas in Place (IGIP) in X gas field based on 2P reserves model is 128.3 BSCF.
A full field 3D numerical simulation model of X gas field (as shown in Figure 1) has been constructed in order to analyze and investigate reservoir fluids distribution and to optimize reservoir fluids production. The seismic interpretation result, exploration and development well log data were used in static reservoir model construction.
Other data, such as fluid properties, rock properties (permeability, relative permeability and capillary pressure), and tubing performance were added to make the dynamic reservoir model.
The dimension of X gas field fine grid is 82i x 69j x 40k, with the size of each grid is 50 x 50 meter. Because the reservoir geomechanics data of X gas field which are required in reservoir simulator (Young's modulus, Poisson's ratio, Biot's constant) are not available, then we will generate those required data by using normal distribution and data from some references. In GA, the random searching technique can be derived from the evaluation of the objective function, or it also might be known as fitness function. The objective function will be generated by using the basic reservoir and geomechanics properties, which are porosity, permeability, gas saturation, Young's modulus and Poisson's ratio.
Figure 2 explain about how GA approach applied in this study. However, reservoir simulator is still required to validate whether the GA can give some good solutions on this production optimization problem. a. This study only observes one development scenario, which is one horizontal well which is drilled only in one layer of reservoir model b. Production constraint applied in this production optimization is 25 MMSCF/day of gas rate and 250 psia of Tubing Head Pressure (THP) c. Because the geomechanics data are not available, then we will generate it using Gaussian distribution (normal distribution) d. The vertical well radius r ev for horizontal well drainage area calculation is 50 meter.

Drainage Area of Horizontal Well
Horizontal well drainage area is one of the important variables needed to be determined in the production forecast of horizontal wells. By knowing the drainage area, we can determine the wells spacing required so that there will be no production interference between the wells. The drainage area can be calculated by using either of the methodologies proposed by Joshi in 1990 and1991. These two methods consider only the geometry of the drained surface, and no other factors such as reductions in productivity that triggers pressure loss through the productive section (Saavedra et.al, 2001). Calculation methods are illustrated in Figure 3(a) would represent the Method 1. It is assumed that the area is consisted of a rectangle with base L and height 2r ev and two semi-circles mas.ccsenet.org Modern Applied Science Vol. 11, No. 9;2017 with a radius of r ev (Saavedra et.al, 2001). The drainage area of horizontal well A hw is calculated with the following equation: where L is the length of the horizontal well and r ev is the drainage radius of a vertical well in the same area. Figure 3 (b) represents the Method 2. Method 2 uses the assumption that the drainage area is an ellipse, calculated with the following equation: Where a is half the major axis of the ellipse and b is half the minor axis of the ellipse.

Geomechanicss
Described model of Darcy flow imply that the flow assumed in laminar condition, since geomechanics affect reservoir permeability that will be smaller during production lifecycle, and if it not applied in fractured reservoir so the Darcy's law still valid .
The geomechanics model considered to use an isothermal, poroelastic medium with small strains, while three basic relationships for a poroelastic theory (Shu, 2003) are: a) Stress equilibrium c) Stress can be defined in strain terms (poroelastic theory), i.e.: where Young's modulus ( ) and Poisson's ratio ( ) as variable in shear modulus = ( , ) = ( ) , and . is Biot's poroelastic constant which represents the interactions between rock and fluid. Biot's constant value is ranging from 0 to 1. = 0 means that no interaction between rock and fluid. And also, the pore pressure does not affect the rock stress. From those stress equations, we could conclude that stress is related with Poisson's ratio and Young's modulus. These two geomechanics parameters will be used later for the input of reservoir simulator. Table 1 shows the average value of geomechanic properties that will be used in this GA study.

Genetic Algorithm
Optimization problem solving using conventional search methods are often not reliable to determine the optimum solution in non-linear multimodal functions. Even, some of the conventional search methods apply trial and error procedure by using a lot of amount of data and efforts to obtain the optimum solution. This method is exhausting to be done. In several cases, an efficient random search method might be required to avoid trial and error procedures and to decrease time and efforts consumed to obtain the optimum solution.
A Genetic Algorithm (GA) is a directed random searching technique inspired by natural evolution mechanism in optimization problem, invented by Holland (Holland, 1975). There are four conditions that influence the evolution process (Holland, 1975):  Reproduction ability of the organism  Reproduction process  Individual diversity in a population  Difference in survival ability GA can find the global optimum solution in complex multi-dimensional search spaces. Survival of the optimum fittest among a population of individuals, selection criteria, and a reproduction strategies are inspired from the natural evolution theory and used as the genetic operators in the artificial environment.
A simple GA flowchart methodology is described in Figure 4. Selection mechanism in GA uses the genetic operators employing to manipulate individuals in a population over several generations to improve their fitness in optimizing the solution gradually. Individuals in a population can be formed into chromosomes link and usually represented as a strings of binary numbers. In this case, conversion of real numbers into binary numbers and re-conversion techniques is very important in selecting the optimal solution.
Genetic algorithms do not use much knowledge about the problem to be optimized. GA work with codes to represent the parameters influencing the optimization problem that they do not deal directly with those parameters. Thus, some matters to be solved in using GA application include: 1. how to represent the problem parameters, 2. how to create the initial population of possible solutions, 3. how to select or devise a suitable set of genetic operators, and 4. how to design the interface between the problem environment and the GA model to know the quality of the found solutions whether they are the optimal solution for the problem.  In this study, the control parameters for the GA are listed as shown in Table 2. Choosing a large size of population results in a bigger probability of converging to global optimum, however, it consumes more computation time, and thus, we chose 100 as the population size. To give a quick convergence to the global optimum solution, we chose a high value of the crossover probability that is 0.95. On the other hand, we chose a very low value of the mutation probability that is 0.05 to avoid instability and prominent diversity of the population. In this study, we used a maximum number of iteration as a stopping criterion for the GA calculation.
Here we used 1000 iterations. We have observed that this number of iterations is sufficient to achieve a global optimum solution (Ariadji et.al, 2012).

The Optimization Problem
The objective function at grid (x,y) is a simple formula which is generated from the basic reservoir and geomechanics properties. In this study, we will consider the drainage area of horizontal well. Thus, we need to define and discretize which grids are included as horizontal well drainage area. Then, we will use those grids in constructing the objective function. The purpose of the discretization is to simplify the horizontal well drainage area and to ease the making of GA code. The illustration of dicretization of horizontal well drainage area is shown in Figure 5 and Figure 6. The gray grids represent the horizontal well drainage area.  ( , ) = gas saturation at grid (x,y), fraction; ( , ) = permeability at grid (x,y), mD; ( , ) = Young's modulus at grid (x,y), bar; ( , ) = Poisson's ratio at grid (x,y), fraction.
The numerator part of f consists of the multiplication of basic reservoir properties, i.e., porosity, permeability and gas saturation. As we know that those reservoir properties have great influence in production performance. Higher value of porosity and gas saturation at (x,y) means that (x,y) will have more gas saturated in porous rock instead of the properties with lower value. And also, the fluid will be easier to flow through the rock that highly permeable. Thus, higher numerator part of f means better production performance. The denumerator part of f consists of the multiplication of geomechanics properties, i.e., Young's modulus and Poisson's ratio. From the mas.ccsenet.org Vol. 11, No. 9;2017 sensitivity analysis that we have done before, we could see that the lower of those geomechanics properties will give better production performance. Then, we put it as the denumerator of f.
Because of the horizontal well will be drilled through more than a grid in reservoir model, then the objective function F that we constructed is a summation of f(x i ,y i ) where (x i ,y i ) is the horizontal well drainage area . Note that f and F have no specific physical interpretation. However, we believe that it will be sufficient to represent a model for optimizing the production rate. Note that all of reservoir and geomechanics properties at (x,y) have been normalized for keeping the fairness of the influence of every properties.
Hence, we can obtain the optimum direction by using the (x o ,y o ) and (x,y) coordinates, as mathematically explained below. Suppose that vector u is on a vector that is constructed by the horizontal well, and vector v is a north direction vector. Thus, the optimum direction of the horizontal well α in the counter clockwise direction could be obtained from where u = (x 0 -x y 0 -y) and v = (0 1). We can also determine the length of horizontal well L (in meter) by using Pyhtagoras theorem, which is = ( − ) + ( − ) .
(15) In this study, we have known that the initial location (kick-off point) of horizontal well is located in grid = 42 and = 34. The constraint that we used in this study is that the horizontal well cannot be drilled over the reservoir boundaries (i.e., 26 ≤ ≤ 54 and 18 ≤ ≤ 53 ). And also, we would optimize the optimum direction α for a certain length L of horizontal well (i.e., for = 500 meter). Thus, the optimization problem of this study is

Economic Evaluation
To determine the most optimum length, we should do the economic analysis. Figure 7 shows the Production Sharing Contract (PSC) that is used for the economic calculation. : 20% For gas market, the economic analysis will be only done at the plateau period (based on the contract of plateau rate with the buyers). From the assumptions above, we could tell that the difference investment for each case of horizontal well will be determined by drilling cost per foot. The vertical lengths of wells before Kick Off point for all cases are similar, i.e., 1589 meter. Drilling cost depends primarily on well location and well depth (Bourgoyne, 1991). The location of the well will govern the cost of preparing the wellsite, moving the rig to the location and the daily operating cost of the drilling operation (Bourgoyne, 1991).
Drilling costs tend to increase exponentially with depth (Bourgoyne, 1991). In this study, we examined two types of drilling cost. First, we used the constant drilling cost per foot (CPF), i.e., 6000 US$/ft. Second, we used an = constants which depend on the well location. Because the empirical drilling cost data are not available, we will use the assumptions to estimate the value of a and b to give realistic horizontal well drilling cost estimations in recent time, i.e., a = US$ 5,000,000 and b = 3.45 x 10 -4 ft -1 . The figure below shows the comparison between two types of drilling cost. We can see that drilling cost of constant CPF gives linear profile of drilling cost while the exponential drilling cost gives exponential line profile of drilling cost. The exponential drilling cost tells us that if we drill longer horizontal well, the cost will increase exponentially due to technical considerations, such as cutting transport (longer horizontal well causing harder to transport the cutting), pipe sticking, mechanical problems, etc. which is more realistic. However, the drilling cost with constant CPF is frequently used sometimes because it is more simple and easier.

Trial and Error
Before we use Genetic Algorithm (GA) to determine the optimum horizontal well length and direction, we have also done trial and error procedure in X gas reservoir model with geomechanic effects. The purpose of this trial and error is to give some comparisons for GA method whether this GA method could give a better and faster solution.
We will use trial and error procedure to gain only the best direction of horizontal well. The horizontal well section will be in fixed length, i.e. 500 meter. Figure 9 shows the horizontal wells that will be checked by trial and error. The results summary of the trial and error is shown in Figure 10 and Table 3. Based on the trial and error, we obtained that the optimum direction for 500 meter of horizontal well section length is 45 0 which yields 7.45 years of plateau time.   Table 4 shows the result of GA simulation. Figure 11 shows the new proposed horizontal well that is resulted by GA method with L = 500 meter. The two types of drainage area give similar result of the optimum direction angle, i.e., 69 0 . Figure 11. The new proposed well that is resulted by GA method

Genetic Algorithm Results
After we obtained the optimum direction, then we will validate the new proposed well using commercial reservoir simulator. Figure 12 shows the result from the reservoir simulator for L = 500 meter and α = 69 0 . From Table 4, we could see that the GA method yields a good result, although the plateau time is shorter compared with trial and error results (the difference is only 0.02 years). However, the GA could give the optimum result faster without using a lot of reservoir simulator runs. Moreover, the GA is a very reliable and efficient tool to determine the optimum direction and length of horizontal well.   Table 5 shows the results from GA simulation for all horizontal well length L. It shows that two types of the drainage area give similar results, except for L = 100 meter. Figure 13 shows the visualization of horizontal wells resulted from GA simulation. We can see that all of horizontal wells have slightly different direction angles. All of horizontal wells have direction angles ranging from 0 0 -90 0 . It means that that area has the best reservoir and geomechanic properties. After we obtain those results from GA simulation, then we will validate it using commercial reservoir simulator to see how long the plateau time that will be resulted. Figure 14 show the results from reservoir simulator. We can see that all of the horizontal wells give quite good results. The plateau time is ranged between 5 -7.5 years.  Vol. 11, No. 9;2017 Table 6 and Figure 15, we could see that the longer horizontal wells give the longer plateau time. However, when the horizontal well length L is longer than 200 meter, the increment of plateau time is quite small and not significant (i.e., less than 0.5 years). Thus, we can conclude that the horizontal well with α = 56 0 and L = 200 m resulting 6.59 years of plateau time is the optimum horizontal well scenario in X gas field. However, we still need to evaluate all of these results economically to determine which is the best horizontal well direction and length.

Economic Evaluation
All the tables and figure below shows the result of the economic analysis calculation.   From the above results, we could see that all of scenarios from both of drilling cost estimation methods give interesting results for the investor. For example, the Net Present Value (NPV) for all scenarios are ranged from 19 million US$ to 24 million US$. In addition, the internal rate of return (IRR) of all observed horizontal wells length give good results, i.e., above 50% (above MARR = 20%). It shows that all of these projects are very interesting to be implemented.
From Figure 16, constant CPF drilling cost gives more optimistic NPV than exponential drilling cost. It happens because based on Figure 8, the exponential drilling cost has higher cost value than constant CPF drilling cost. Then, with same gross revenue, the exponential drilling cost will give lower NPV results.
However, based on Figure 16, both of drilling cost estimation methods lead to the same conclusion. If we drill longer horizontal well then the increment of NPV will decrease, and moreover, the NPV will decline after a certain horizontal well length. Then, we can decide which the most optimum horizontal well length is. We could see that L = 300 meter gives the largest Net Present Value, i.e., 24,331.78 MUS$ for constant CPF drilling cost and 22,826.12 MUS$ for exponential drilling cost. Thus, from the above economic analysis, we can conclude that the most optimum horizontal well direction and length to be developed in X gas field is the horizontal well with α = 59 0 with L = 300 meter.

Conclusion
The proposed objective function in GA method without loosing accuracy could significantly reduce the amount of reservoir simulation work in finding the best horizontal well direction and length. Using the limited reservoir properties data and few reservoir simulation runs, GA could obtain acceptable and accurate results. The proposed GA procedure can be an effective and efficient alternative tool to determine the optimum horizontal well direction and length in field development scenario, technically and economically.

Copyrights
Copyright for this article is retained by the author(s), with first publication rights granted to the journal.
This is an open-access article distributed under the terms and conditions of the Creative Commons Attribution license (http://creativecommons.org/licenses/by/4.0/).