Study of Hydraulic fracturing in water dominated Geothermal field using experimental design and numerical simulation

As a water-dominated geothermal field continues in operation, production rate and pressure will decrease. To enhance production rate, hydraulic fracturing stimulation is normally conducted. This fracturing stimulation is not always successful due to reasons related to the uncertainty of numerous design and reservoir parameters. This study is conducted to answer the uncertainty of hydraulic fracturing parameters problems specifically in a water-dominated geothermal reservoir by means of extending the experimental design for geothermal field applications. In this study, experimental design will be applied to study the uncertainty parameters of hydraulic fracturing stimulation in a water dominated reservoir using CMG STARS simulator followed by statistical analysis using MiniTab software and performing Monte Carlo probabilistic simulation. The result provides necessary steps to conduct hydraulic fracturing study in CMG STARS and proxy polynomial equation that well describes mass rate response output. This study also concluded that fracture half-length appears to have the highest effect on mass rate even it is not too different from the other two parameters. Last, the most probable mass rate output is estimated to be 124 kg/s which is 7.1 times higher than the base case mass rate. The first section in your paper


Introduction
Geothermal energy as one of the alternative energy today continue to develop and grow over time. Using earth heat source to generate electricity or for other direct uses, it is said that geothermal energy is more sustainable and more environmentally friendly compared to fossil fuels, e.g. oil and gas, even though the principles are mainly similar, differ in the fluid being extracted. When a geothermal well is drilled obtaining a low production rate, the well could be immediately abandoned. However, when it has high production rate, produced steam can be used either for direct uses or indirect uses such as electricity utilization. When a geothermal field continues its production through years, the reservoir performance such as production rate and pressure will eventually decrease over time. If this situation is not properly handled, the production rate could fall below the economic limit, which could then lead to the abandoned of the well. The low production rate that occur directly after drilling and after years of production can be caused by several things such as formation damage, completion effects or lack connectivity to main fluid conduits [1].
The stimulation in geothermal wells was started in 1970 in New Mexico. This was adopted from the application in petroleum industry, which has been developed since 1947 in the United States of America. Since then, hydraulic fracturing has been attempted in geothermal formations as means to stimulate both production and Injection wells [2]. Hydraulic fracturing can be described as a stimulation method that is used to increase the performance of a well by increasing the productivity index or the injectivity index by increasing the permeability around the wellbore, eliminate formation damage and connecting fractures that already exist in the reservoir.
Experimental design is described as a scientific approach that allows the researcher to comprehend the process and to determine how an input affects the response or the desire output by changing inputs. This method is commonly used for resource assessment in the oil and gas industries [3], yet there are still very few applications in the geothermal industry and stimulation methods. Although there are authors starting to use this method for resource assessment in geothermal fields such as Quinao and Zarrouk [4], Riswan [5] and Quinao and Zarrouk [6].  In this study, the experimental design will be applied to assess the uncertainty parameters of hydraulic fracturing stimulation in a water-dominated reservoir where in this case, Wairakei geothermal field. Using CMG STARS software, the mass rate of production is selected as the objective function. The main goal of this research is to analyze the significances of the input parameter to the mass rate production response after one year of production period. The next objective is to develop a proxy polynomial equation describing the mass rate production for the period of one year after the hydraulic fracturing treatment. This proxy is as a function of the uncertain parameters or controlled parameters. The last objective is to forecast the most probable mass rate using Monte Carlo simulation. This study will also offer the schematic of how to conduct a hydraulic fracturing in CMG STARS. The workflow of this study can be seen in figure 1.

Model Description
An existing model of Wairakei water-dominated geothermal field, New Zealand developed by Lukmana [7] is used as the reference of this study. Wairakei field starts production at the 1950s with an almost linear decline rate of 3% per year until the early 1970s followed by a 1% decline rate until the late 1980s [8]. Up to the 2000s, production rate from the main production area at Wairakei geothermal field is decreasing, except from the Te Mihi area due to injection, infill wells and reduction of the operating wellhead pressure [9]. The production profile of the Wairakei Field is depicted shown at figure 2. Figure 2. Wairakei average monthly production rate for the Eastern Borefield (EB), Western Borefield (WB) and Te Mihi area including injection rate [9]. A CMG STARS model was built adopting the conceptual model of Wairakei geothermal field in Figure 3 as a reference. From the conceptual model, Wairakei field is located on a flat terrain categorized as a water-dominated geothermal system. The reservoir temperature is estimated to be ranging from 200˚C until 260˚C using hot rocks as it heat source. It can also be seen that the reservoir condition generated a thin layer of steam cap at shallow depth.
The synthetic model has 11 horizontal layers consisting of 968 grid blocks with different elements covering equal to of the field. The 2 upper layers represent the surface, the 3 following layers represents the cap rock, followed by 5 layers of the reservoir and 1 layer at the bottom as heat source basement. The rocktype configuration for each layer are based on Elfajrie and Syihab [10]. Figure 4 shows the layers and grid blocks of the model. The natural state condition is than calibrated and compared for the temperature and pressure vs depth between the actual field data, the TOUGH2 model by Elfajrie and Syihab [10] and the CMG STARS model by Lukmana [7] is represented on figure 5.

Hydraulic Fracturing Parameters
Three parameters were investigated in this study which covers the fracture geometry and the fracture permeability. These are fracture half-length, fracture width and fracture permeability. Fracture halflength describes the length of 1 wing of fracture, the fracture width describes how wide the fracture is opened and the fracture permeability describes the permeability of the fracture itself. The fracture height was not investigated due to the inability of the model and CMG STARS to adjust the fracture height. Therefore, it is assumed that the fracture height is the same along the fracture

Experimental Design
A full factorial design for the three parameters (n=3) at two levels (max [+1] and min [-1]) was used to identify the significant parameters affecting fluid production due to hydraulic fracturing treatment after 1 and 30 years of production. It is worthily noted that the range of parameter values will limit the results within the values identified. An engineering senses has guided these values. The values of max and min variables used in this study were based on [12] and [13] presented in table 1. Eight numerical models were then built using Mini-Tab according to the full factorial design presented in table 2.

Model Development
As mentioned before, hydraulic fracturing option in CMG STARS is not active for some reasons differing from CMG GEM and CMG IMEX. Hence, in order to do the further study, a pre-study is needed and conducted to determine the right steps to have a hydraulic fractured grid in CMG STARS which is still a representative model of the numerical model performance. In short, hydraulic fracturing option in CMG reform the existing grids by splitting the grid into the desired size of fracture and change the value of permeability where the hydraulic fracturing is conducted to form a new grid representing the fracture inside the reservoir. The validation is conducted by creating two numerical models from the same base model but with a different process. The first model is formed by changing the grids and permeability of the base model in the CMG STARS simulator while the second model is formed by changing the grids and permeability in the same exact place and value like the first model using CMG GEM. For the second model, grid and spatial data from the base model is imported toward the CMG GEM where the change in grid and permeability is conducted and then reimported back into CMG STARS to be compared with the previous model. Both models are then compared by the pressure versus depth condition, temperature versus depth condition and the fluid flow profile for 30 years of production from a well under identical coordinates and condition. If the results are not exactly the same, model b will then be re-evaluated and modified to gain the exact same output. These modifications will then be used as the procedure for the further study to gain a representative result.
After the result of pre-research is obtained, using the CMG STARS model [7] as the base model, hydraulic fracturing is conducted to the model according to the scenarios values based on table 1 and table 2. As it has been explained before that hydraulic fracturing option cannot be conducted in the CMG STARS, CMG-GEM is used by importing the grid and spatial data from the base model and re-importing the fractured grid and special data to CMG STARS base model. Followed by re-entering manual and automatic modification to the input data according to the pre-research. The hydraulic fracture itself was performed to the J direction at the centre of the reservoir at each layer of the reservoir which are layer 6,7,8,9 and 10. The 3D model of the fractured reservoir can be seen in figure 6. Simulation of 30 years production is then conducted in the exact similar condition for the 8 scenarios.

Result and Analysis
As a complimentary goal, the procedure to build hydraulic fractured grids in CMG STARS is obtained and it is presented in figure 7. As it is expected, several properties are changed when the grid and spatial data is re-imported into CMG STARS. The other issue is that the CMG STARS cannot calculate the well index because of the grid is too small so it must be manually configured.  Figure 7. Recommended procedure to conduct hydraulic fracturing study using CMG STARS.
After 30 years production simulation using constant bottom hole pressure for the base model and each scenario of experimental design is conducted, the result is as shown in figure 8. As predicted, the mass rate during 30 years of production of all the scenarios is greater than the based model. Maximized at all maximum parameters which are scenario 8 showing a significant difference compared to the other scenarios. The mass rate after 1 year of production for all scenarios is presented in table 3 including the fold of increase.  Figure 8. Eight scenarios mass rate production profile over 30 years of production graph.
From the result in table 3, a proxy polynomial model is generated using MiniTab to represent the mass rate after 1 year of production for the Wairakei geothermal field as shown in equation (1) (1) Where A shows fracture half-length ( ), B shows fracture width ( ) and C shows fracture permeability ( ). Since the experimental design chosen to be used at the start of the analysis is the full factorial design, the interaction of variables is also analyzed resulting two-way interactions and threeway interactions.
Based on the statistical analysis using Minitab, from the three parameters there are no significant parameters. This means that the effects of three parameters on the mass rate after hydraulic fracturing simulation are nearly alike. In figure 9, the Pareto chart ranks the linear effects of the parameters on the mass rate by decreasing order. It shows that A, which is fracture half-length have the highest effect on mass rate even it is not too different from the other two parameters. Although A has the highest effect, it cannot be categorized as a significant parameter because it doesn't pass the red line. The normal plot in figure 10 shows that all the parameters have positive effects on the mass rate of Wairakei geothermal field even it is categorized as a not significant parameter.  Furthermore, the proxy polynomial model was used to carry out Monte Carlo probabilistic analysis resulting the graph in figure 11. Using 10000 iterations, the parameters were randomly sampled from each parameter range of the max (+1) and the min (-1). From the graph it can be concluded that hydraulic fracturing stimulation in Wairakei geothermal field can increase the mass rate production into P10 of 47 kg/s which is 2.7 times higher than the base case, mean of 124 kg/s which is 7.1 times higher than the base case and P90 of 260 kg/s which is 14.8 times higher than the base case. It must be noted that the spread of the probabilistic mass rate depends on the range of parameters mentioned before. The positive skewness shows that the mass rate will probably have a value under the average data.  Figure 11. Monte Carlo simulation result for mass rate one year after hydraulic fracturing stimulation.

Conclusion
From the results and the analysis of the research performed, it can be concluded that: • hydraulic fracturing can be conducted and analyzed using CMG STARS simulator software by the procedure that has been obtained from this research in the result and analysis section. • From the three parameters that were investigated in this study, fracture half-length appears to be the most sensitive parameter affecting the output of mass rate production, even though it is just slightly more sensitive compared to the others. From the three parameters, there are none that is classified as a significant parameter according to the Pareto chart and the normal plot. • The proxy polynomial equation describing mass rate production 1 year after hydraulic stimulation on Wairakei geothermal field is performed as a function of the three parameters is obtained as shown in equation 1. • From probabilistic Monte Carlo simulation, it is obtained that the most probable mass rate obtained is 124 kg/s which is 7.1 times higher than the base case flowrate. A better proxy polynomial equation must be obtained to get a better representative result as it is find that the range of the probable result is still large and the result of the Monte Carlo simulation did not produce a normal distribution.