Hydraulic Fracturing Treatment Optimization for Low Permeability Reservoirs Based on Unified Fracture Design

Hydraulic fracturing optimization is very important for low permeability reservoir stimulation and development. This paper couples the fracturing treatment optimization with fracture geometry optimization in order to maximize the dimensionless productivity index. The optimal fracture dimensions and optimal dimensionless fracture conductivity, given a certain mass or volume of proppant, can be determined by Unified Fracture Design (UFD) method. When solving the optimal propped fracture length and width, the volume and permeability of the propped fracture should be determined first. However, they vary according to the proppant concentration in the fracture and cannot be obtained in advance. This paper proposes an iterative method to obtain the volume and permeability of propped fractures according to a desired proppant concentration. By introducing the desired proppant concentration, this paper proposes a rapid semi-analytical fracture propagation model, which can optimize fracture treatment parameters such as pad fluid volume, injection rate, fluid rheological parameters, and proppant pumping schedule. This is achieved via an interval search method so as to satisfy the optimal fracture conductivity and dimensions. Case study validation is conducted to demonstrate that this method can obtain optimal solutions under various constraints in order to meet different treatment conditions.

This paper focuses on the treatment optimization for low permeability reservoirs. The conventional treatment optimization is to obtain optimal treatment parameters in order to maximize the production with some constraints such as the cost and treatment capability. Usually, the net present value (NPV) [24] is used as the optimization objective. Taking various treatment permeability reservoirs [59], horizontal wells [60,61], and coalbed methane wells [16]. Due to the high flow velocity of gas in fractures, the non-Darcy flow effect is combined with the UFD method to optimize the design of gas reservoirs [13,16]. The UFD method is also used for optimizing the design of fractured horizontal wells in heterogeneous tight gas reservoirs [19]. However, the above works only focus on the optimization of fracture number, fracture spacing and fracture geometry dimensions. Treatment optimizations are not mentioned [62][63][64][65]. In order to satisfy the optimal fracture dimensions, the treatment parameters should also be optimized further.
To this end, this paper proposes a method coupling the fracture geometry optimization and treatment parameter optimization to maximize the dimensionless productivity index, given a certain mass or volume of proppant. By introducing a desired proppant concentration, both the optimal fracture dimensions and treatment parameters can be solved rapidly. The methods are detailed in Section 2. A case study and some discussions are described in Section 3. Moreover, the proppant mass or volume can also be optimized further considering the fracturing cost or the treatment capability, which will be discussed in Section 3. Conclusions are drawn in Section 4.

Fracture Geometry Optimization Based on UFD
UFD method gives the analytical solution for the optimal dimensionless fracture conductivity in order to maximize the dimensionless productivity index, given a certain mass or volume of proppant. The optimal dimensionless fracture conductivity is the function of a proppant number. The details about UFD method are given in Appendix A.
When the optimal dimensionless fracture conductivity is determined, the optimal fracture half-length and width can be obtained accordingly: x opt = k f V f C f Dopt kh p 0.5 (1) w opt = C f Dopt kV f k f h p 0.5 (2) where: k f is the permeability of the propped fracture (md); k is the permeability of the reservoir (md); h p is the thickness of the reservoir, in this model, it equals to the fracture height (m); C f Dopt is the optimal dimensionless fracture conductivity; V f = V p 2 is the single wing volume of the propped fracture (m 3 ); V p is the total volume of the propped fracture (m 3 ).
As mentioned above, since the volume and permeability of the propped fracture are unknown, the fracture dimensions cannot be solved directly. Note that the volume and permeability of the propped fracture are both related to the proppant concentration, which reflects the uniformity of the proppant distribution. In theory, if the distribution of the proppant in the fracture is completely uniform, then the value of the proppant concentration can reach the bulk density of the proppant on the ground. However, due to proppant transport in the pipe and fracture, the distribution of the proppant is never completely uniform. Thus, the value of the proppant concentration should be less than the bulk density. Nevertheless, the proppant concentration can be maintained at a suitable level by optimizing the treatment parameters. So, we can preset a desired proppant concentration value and then calculate the volume and permeability of the propped fracture. The corresponding optimal fracture dimensions can be obtained accordingly, which can be reached through treatment parameter optimization. Thanks to the corresponding relationship between the optimal fracture dimensions and the proppant concentrations, once the optimal fracture dimensions are reached, the desired proppant concentration can be obtained. This paper proposes an iterative method to solve the optimal fracture dimensions given a preset proppant concentration. Suppose all the proppants stay within the range of the thickness of the reservoirs; then the relationship of the mass and concentration of the proppant in the fracture is as follows: where M p is the proppant mass (kg), which is constant if the fracturing scale is fixed. C s is the proppant concentration (kg/m 3 ). Then, the fracture width Equation (2) can be rewritten as follows: The permeability of the fracture is related to the proppant concentration per unit area and closure pressure. The proppant concentration per unit area is defined as: Then, function of the permeability of the fracture can be written as: where P c is the closure pressure (MPa), which can be determined by field tests. There is no explicit equation for the function F, which can be obtained from the permeability curve of different proppant concentrations per unit area and closure pressures, given one certain type of proppant. The permeability curves can be obtained through lab experiments. Given the proppant mass M p and the desired proppant concentration C s in the fracture, based on the optimal dimensionless fracture conductivity C f Dopt , the optimal fracture dimensions and permeability can be solved by the fracture half-length Equation (1), width Equation (4) and permeability Equation (6) using the iterative method.
Firstly assume an initial permeability k f , and substitute it into the fracture width equation to solve w opt . The proppant concentration per unit area C f can then be solved by w opt , and the new permeability k f can be obtained from the permeability curve. If the difference between the new value and the assumed one is within a certain small range, then the new permeability and corresponding fracture width can be obtained. Otherwise, replace the assumed permeability with the new value and repeat the above procedure until convergence is achieved. Once the fracture width and permeability are calculated, substitute them into Equation (1) to solve the optimal fracture half-length.

Treatment Optimization through a Fast Semi-Analytical Fracture Propagation Model
Once the optimal fracture dimensions are determined, the treatment parameters also need to be optimized. There are many parameters that influence the fracture dimensions, and the value of each should be within a certain range considering the treatment feasibility. Moreover, the treatment optimization result should satisfy the fracture length and width simultaneously. In order to solve this optimization problem, the treatment parameters need to be tuned and the fracture dimensions calculated repeatedly. Thus, the computation of the fracture propagation model should be rapid. Currently, the widely-used PKN [44][45][46] analytical model is very fast. However, its solution does not rigorously satisfy the flow continuity equation. Pseudo-3D [47][48][49][50] and fully-3D [37] models can solve the geometric dimensions of the propped fracture accurately. However, these methods are time-consuming and rarely suitable for treatment parameter optimization. This paper proposes a fast 2D semi-analytical fracture propagation model, which assumes that the fracture height is constant and the proppant is transported in the fracture at the same speed as the sand-laden fluid.
Energies 2018, 11, 1720 5 of 23 In this method, the total pumping time is divided into 1, 2, . . . k, . . . k end segments. The length of each time segment is ∆t and the pumping fluid volume during each segment is constant, ∆V. Then, after k time segments, the total pumping time is t = k∆t and the fluid elements in the current fracture can be recorded as 1, 2, . . . i, . . . k accordingly. The following values of each fluid element are calculated according to the total pumping time t: the length L i (t), the distance of its back interface to the fracture entry d i (t), the cumulative fluid loss volume V li (t), the remaining volume V i (t), the proppant mass M pi , and the proppant concentration C si (t). The proppant mass in each element is related to the pumping schedule and supposed to be invariable during fracture propagation. Additionally, we record the current total fracture half-length L(t) and width distribution W(x, t) at the end of each pumping time segment.
During one new pumping time segment, the length of each existing fluid element L i (t) varies continuously due to the change of the fracture width at its location and the fluid loss within it. As more fluid is injected, the existing fluid element moves forward. Hence, the distance of the back interface to the fracture entry d i (t) also changes. Due to fluid loss, the volume of the fluid element decreases continuously. For the pad fluid element, its volume decreases to 0 eventually. For the sand-laden fluid element, thanks to the existence of proppant, its final volume will be greater than 0. The maximum concentration C smax should be set to ensure proppant transport in the fracture. Otherwise, a sand plug will occur. C smax can be determined empirically. Figure 1a shows the fluid element distribution after k time segments. Figure 1b is the fluid element distribution after k + 1 time segments.  [44]: where: where ν is the rock Poisson's ratio, E is the rock elastic modulus (Pa), H is the fracture height (m), C is the fluid loss coefficient (m/s 1/2 ), µ a is the apparent viscosity of the power law fluid (Pa·s), K is the where: The volume of the fluid loss for each element is limited as follows: for sand − laden fluid (12) Because the fracture half-length L(t) cannot be solved explicitly, an iterative method is used. First, an initial fracture half-length is assumed, and then the fracture width distribution is solved by Equation (8). We can substitute it into Equation (10) to check whether the continuity equation is satisfied. If it is not, then we assume another fracture half-length and calculate again. In general, the fracture half-length increases with time. So, the assumed value can be set as the previous value plus a small increment until it finally satisfies the continuity equation. When both the fracture width and half-length have been determined, the fluid element space in the fracture should be reassigned. The remaining volume of each fluid element after pumping for k time segments is: According to the remaining volume of each fluid element, starting from the fracture tip, we calculate in turn the distance of the back interface of element 1, 2, . . . i, . . . k to the fracture entry: As the model satisfies the continuity equation, then the following condition will be satisfied naturally: Thus, the length of each fluid element becomes: The proppant concentration in the fluid element can be calculated as: Because fluid loss occurs during fracturing, sand-laden fluid that was injected early is subject to a long loss time; hence, the sand ratio of the injected fluid at early times should be relatively small. As the loss time of sand-laden fluid injected at a later stage is short, the sand ratio of the injected fluid can be relatively large. Thus, at the end of injection, the proppant in the whole fracture will be evenly Energies 2018, 11, 1720 7 of 23 distributed. For the convenience of operation, the staged proppant pumping schedule is usually used. It means that one constant sand ratio is maintained over a period of time. When the injection of these portions of sand-laden fluid is finished, the sand ratio is then changed. In this paper, the sand ratio is increased according to the following power function to ensure that the proppant is distributed as evenly as possible, so as to reach the desired proppant concentration: where S is the sand ratio (representing the ratio of proppant to fluid volume), a and b are coefficients (in this paper, b is referred to as the proppant pumping curve index), and t is the proppant pumping sequence. If all the proppants are to be injected through 10 times, then the value of t is an integer from 1-10. In this paper, the proppant pumping times and the maximum sand ratio will be determined empirically. Then, for a given value of coefficient b, the corresponding coefficient a can be determined. Thus, the optimization of the complicated proppant pumping schedule can be simplified as the optimization of the single coefficient b.
After the injection has finished, the pad fluid continues to be lost and the front of the fracture will close. Thus, the half-length of the propped fracture can be set as the length in which proppant is present. If the proppant concentration in the direction of the fracture width has not reached a certain value, then the fracture will close and the fracture width will decrease until the proppant concentration reaches the desired value. The final width is the propped fracture width. The average fracture width is used in the UFD method regardless of the fracture shape. Figure 2a shows the state of the fracture immediately after the end of injection. Figure 2b shows the state of the propped fracture.
Energies 2018, 11, x FOR PEER REVIEW 7 of 23 paper, the sand ratio is increased according to the following power function to ensure that the proppant is distributed as evenly as possible, so as to reach the desired proppant concentration: where is the sand ratio (representing the ratio of proppant to fluid volume), and are coefficients (in this paper, is referred to as the proppant pumping curve index), and is the proppant pumping sequence. If all the proppants are to be injected through 10 times, then the value of is an integer from 1-10. In this paper, the proppant pumping times and the maximum sand ratio will be determined empirically. Then, for a given value of coefficient , the corresponding coefficient can be determined. Thus, the optimization of the complicated proppant pumping schedule can be simplified as the optimization of the single coefficient .
After the injection has finished, the pad fluid continues to be lost and the front of the fracture will close. Thus, the half-length of the propped fracture can be set as the length in which proppant is present. If the proppant concentration in the direction of the fracture width has not reached a certain value, then the fracture will close and the fracture width will decrease until the proppant concentration reaches the desired value. The final width is the propped fracture width. The average fracture width is used in the UFD method regardless of the fracture shape. Figure 2a shows the state of the fracture immediately after the end of injection. Figure 2b shows the state of the propped fracture. This simulation needs to find the optimal treatment parameters in order to satisfy both the fracture half-length and width. So, the objective of the treatment optimization is to minimize the following error function: where is the calculated fracture half-length (m), is the optimal fracture half-length (m), is the calculated fracture width (m), and is the optimal fracture width (m). The basic flowchart of the proposed fracturing treatment optimization method based on the UFD is shown in Figure 3. This simulation needs to find the optimal treatment parameters in order to satisfy both the fracture half-length and width. So, the objective of the treatment optimization is to minimize the following error function: where x is the calculated fracture half-length (m), x opt is the optimal fracture half-length (m), w is the calculated fracture width (m), and w opt is the optimal fracture width (m). The basic flowchart of the proposed fracturing treatment optimization method based on the UFD is shown in Figure 3.

Basic Parameters
To validate the proposed method, it was applied to a case study reservoir: that of the Daniudi Gas Field located at the junction of the city of Yulin, Shaanxi Province and the city of Ordos, Inner Mongolia Autonomous Region, China. It contains lower Paleozoic marine carbonate rock as one of the main gas-bearing strata. Its basic parameters are listed in Table 1. The homogeneous rock is presented in this case study for illustration. Nevertheless, the heterogeneous rocks are usually divided into several homogeneous segments [19] and the proposed method can also be used for each segment respectively.

Optimization Method and Results
A proppant volume of 18 m 3 is used to illustrate the optimization process and results. First of all, the desired proppant concentration is set to 1000 kg/m 3 , which is related to the proppant pumping schedule and will be discussed later in the paper. According to the actual situation in the field, the maximum sand ratio is set to 35% and the proppant pumping schedule is divided into eight stages. The maximum proppant concentration during fracture propagation is 700 kg/m 3 . The other treatment parameters will be optimized to satisfy the desired proppant concentration, optimal propped fracture conductivity and dimensions.

Basic Parameters
To validate the proposed method, it was applied to a case study reservoir: that of the Daniudi Gas Field located at the junction of the city of Yulin, Shaanxi Province and the city of Ordos, Inner Mongolia Autonomous Region, China. It contains lower Paleozoic marine carbonate rock as one of the main gas-bearing strata. Its basic parameters are listed in Table 1. The homogeneous rock is presented in this case study for illustration. Nevertheless, the heterogeneous rocks are usually divided into several homogeneous segments [19] and the proposed method can also be used for each segment respectively.

Optimization Method and Results
A proppant volume of 18 m 3 is used to illustrate the optimization process and results. First of all, the desired proppant concentration C s is set to 1000 kg/m 3 , which is related to the proppant pumping schedule and will be discussed later in the paper. According to the actual situation in the field, the maximum sand ratio is set to 35% and the proppant pumping schedule is divided into eight stages. The maximum proppant concentration C smax during fracture propagation is 700 kg/m 3 . The other treatment parameters will be optimized to satisfy the desired proppant concentration, optimal propped fracture conductivity and dimensions. As a matter of fact, all the related treatment parameters can be optimized automatically by an exhausted search method. However, due to the field treatment capability, if certain constraints are set to some parameters [52], then the optimization efficiency can be increased significantly. Next methods demonstrate how the treatment parameter constraints facilitate the optimization procedure and how they affect the optimization result.
Usually, the pad fluid volume is one of the key parameters to be optimized. If the injection rate, fluid rheological parameters and proppant pumping schedule can be selected according to the actual treatment conditions, and then the pad fluid volume needs to be optimized automatically. Taking an injection rate Q of 5 m 3 /min as an example, if the fluid's apparent viscosity µ a is set as 58 mPa·s (for instance, the fluid consistency coefficient K is set to 0.7 Pa·s n , and the flow index n is set to 0.6), and the proppant pumping curve index b is set as 0.63, then the curve for the calculated propped fracture half-length and width versus different pad fluid volumes can be illustrated, as in Figures 4 and 5. For the sake of comparison, the optimal fracture half-length and width are also plotted in the figure as horizontal lines. As shown in the figures, the propped fracture half-length decreases while the propped fracture width increases with increasing pad fluid volume. When the calculated fracture half-length curve intersects the optimal fracture half-length horizontal line and the calculated fracture width curve intersects the optimal fracture width horizontal line, the optimal pad fluid volume can be obtained.
Energies 2018, 11, x FOR PEER REVIEW 9 of 23 As a matter of fact, all the related treatment parameters can be optimized automatically by an exhausted search method. However, due to the field treatment capability, if certain constraints are set to some parameters [52], then the optimization efficiency can be increased significantly. Next methods demonstrate how the treatment parameter constraints facilitate the optimization procedure and how they affect the optimization result.
Usually, the pad fluid volume is one of the key parameters to be optimized. If the injection rate, fluid rheological parameters and proppant pumping schedule can be selected according to the actual treatment conditions, and then the pad fluid volume needs to be optimized automatically. Taking an injection rate Q of 5 m 3 /min as an example, if the fluid's apparent viscosity is set as 58 mPa·s (for instance, the fluid consistency coefficient K is set to 0.7 Pa·s n , and the flow index n is set to 0.6), and the proppant pumping curve index b is set as 0.63, then the curve for the calculated propped fracture half-length and width versus different pad fluid volumes can be illustrated, as in Figures 4 and 5. For the sake of comparison, the optimal fracture half-length and width are also plotted in the figure as horizontal lines. As shown in the figures, the propped fracture half-length decreases while the propped fracture width increases with increasing pad fluid volume. When the calculated fracture half-length curve intersects the optimal fracture half-length horizontal line and the calculated fracture width curve intersects the optimal fracture width horizontal line, the optimal pad fluid volume can be obtained.  As shown from the figures, within a certain range of pad fluid volumes, no optimization result can be acquired. Nevertheless, the design scheme can be adjusted in two ways. The first is to increase the injection rate and the second is to increase the fluid's apparent viscosity. For the first adjusted   As a matter of fact, all the related treatment parameters can be optimized automatically by an exhausted search method. However, due to the field treatment capability, if certain constraints are set to some parameters [52], then the optimization efficiency can be increased significantly. Next methods demonstrate how the treatment parameter constraints facilitate the optimization procedure and how they affect the optimization result.
Usually, the pad fluid volume is one of the key parameters to be optimized. If the injection rate, fluid rheological parameters and proppant pumping schedule can be selected according to the actual treatment conditions, and then the pad fluid volume needs to be optimized automatically. Taking an injection rate Q of 5 m 3 /min as an example, if the fluid's apparent viscosity is set as 58 mPa·s (for instance, the fluid consistency coefficient K is set to 0.7 Pa·s n , and the flow index n is set to 0.6), and the proppant pumping curve index b is set as 0.63, then the curve for the calculated propped fracture half-length and width versus different pad fluid volumes can be illustrated, as in Figures 4 and 5. For the sake of comparison, the optimal fracture half-length and width are also plotted in the figure as horizontal lines. As shown in the figures, the propped fracture half-length decreases while the propped fracture width increases with increasing pad fluid volume. When the calculated fracture half-length curve intersects the optimal fracture half-length horizontal line and the calculated fracture width curve intersects the optimal fracture width horizontal line, the optimal pad fluid volume can be obtained.  As shown from the figures, within a certain range of pad fluid volumes, no optimization result can be acquired. Nevertheless, the design scheme can be adjusted in two ways. The first is to increase the injection rate and the second is to increase the fluid's apparent viscosity. For the first adjusted  As shown from the figures, within a certain range of pad fluid volumes, no optimization result can be acquired. Nevertheless, the design scheme can be adjusted in two ways. The first is to increase the injection rate and the second is to increase the fluid's apparent viscosity. For the first adjusted scheme, if the injection rate is increased to 7 m 3 /min, then new curves are obtained, as illustrated in Figures 6 and 7. As shown from these figures, an optimal pad fluid volume within the given range can be acquired. scheme, if the injection rate is increased to 7 m 3 /min, then new curves are obtained, as illustrated in Figures 6 and 7. As shown from these figures, an optimal pad fluid volume within the given range can be acquired.  For the second adjusted scheme, if the fluid's apparent viscosity is increased to 201 mPa·s (for instance, the fluid consistency coefficient K is set to 0.7 Pa·s n , and the flow index n is set to 0.8), then new curves are as illustrated in Figures 8 and 9. As shown from these figures, the optimal pad fluid volume within the given range can be acquired.
According to the above analysis, the treatment parameters that satisfy the optimal fracture dimensions are not unique. For the above two adjusted designs, if the field equipment can meet the required injection rate, then the first adjusted design will be the better option.
As a matter of fact, more parameters can be optimized. This paper proposes the interval search method for automatic optimization. If we set the injection rate as 7 m 3 /min, the range of other parameters needing to be optimized and the search step lengths are as listed in Table 2. The search range can be tuned according to the actual treatment conditions and the search step length can also be tuned based on the optimization accuracy and the calculation efficiency. Generally, a coarse-to-fine search strategy is used. During this optimization, the search step length for pad fluid volume is initially set to 50, and then tuned to 10 to obtain a more accurate result. The complete optimization results are listed in Table 3.  scheme, if the injection rate is increased to 7 m 3 /min, then new curves are obtained, as illustrated in Figures 6 and 7. As shown from these figures, an optimal pad fluid volume within the given range can be acquired.  For the second adjusted scheme, if the fluid's apparent viscosity is increased to 201 mPa·s (for instance, the fluid consistency coefficient K is set to 0.7 Pa·s n , and the flow index n is set to 0.8), then new curves are as illustrated in Figures 8 and 9. As shown from these figures, the optimal pad fluid volume within the given range can be acquired.
According to the above analysis, the treatment parameters that satisfy the optimal fracture dimensions are not unique. For the above two adjusted designs, if the field equipment can meet the required injection rate, then the first adjusted design will be the better option.
As a matter of fact, more parameters can be optimized. This paper proposes the interval search method for automatic optimization. If we set the injection rate as 7 m 3 /min, the range of other parameters needing to be optimized and the search step lengths are as listed in Table 2. The search range can be tuned according to the actual treatment conditions and the search step length can also be tuned based on the optimization accuracy and the calculation efficiency. Generally, a coarse-to-fine search strategy is used. During this optimization, the search step length for pad fluid volume is initially set to 50, and then tuned to 10 to obtain a more accurate result. The complete optimization results are listed in Table 3.  For the second adjusted scheme, if the fluid's apparent viscosity is increased to 201 mPa·s (for instance, the fluid consistency coefficient K is set to 0.7 Pa·s n , and the flow index n is set to 0.8), then new curves are as illustrated in Figures 8 and 9. As shown from these figures, the optimal pad fluid volume within the given range can be acquired.
According to the above analysis, the treatment parameters that satisfy the optimal fracture dimensions are not unique. For the above two adjusted designs, if the field equipment can meet the required injection rate, then the first adjusted design will be the better option.
As a matter of fact, more parameters can be optimized. This paper proposes the interval search method for automatic optimization. If we set the injection rate as 7 m 3 /min, the range of other parameters needing to be optimized and the search step lengths are as listed in Table 2. The search range can be tuned according to the actual treatment conditions and the search step length can also be tuned based on the optimization accuracy and the calculation efficiency. Generally, a coarse-to-fine search strategy is used. During this optimization, the search step length for pad fluid volume is initially set to 50, and then tuned to 10 to obtain a more accurate result. The complete optimization results are listed in Table 3. Table 2. Ranges of parameters used in optimization.

Variable
Pad Fluid Volume (m 3 )

Proppant Pumping
Curve Index b

General Discussion
For rapid optimization, the injection rate and fluid rheological parameters are preset or tuned within certain ranges according to the actual treatment conditions, and the pad fluid volume needs to be optimized automatically. To this end, it is important to analyze the influences of the injection rate and fluid rheological parameters on the optimization results. Furthermore, the reservoir permeability has a strong influence on the optimization results, which will be also discussed in this section. Although the discussions used the specific data above, they are quite general and very useful for the field optimization design especially using the proposed methods in this paper. Moreover, the findings in this section are consistent well with the established knowledge. They are presented here for completeness.

Injection Rate
Based on the above data, we consider several injection rates: 4, 5, 6, 7, 8 and 9 m 3 /min. Other parameters, including proppant pumping curve index b and fluid rheological parameters are obtained from Table 3 and held constant. Then, the curve describing the relationship between pad fluid volume and injection rate is obtained, as shown in Figure 10. To generate the same fracture with the same fluid and proppant pumping schedule, the required pad fluid volume decreases with increasing injection rate. Note that when the injection rate is small, the required pad fluid volume is excessively large, which is unfeasible in terms of treatment difficulty and cost. Therefore, it is suggested to select a larger injection rate, within the treatment limitations, in order to reduce the amount of pad fluid used.

General Discussion
For rapid optimization, the injection rate and fluid rheological parameters are preset or tuned within certain ranges according to the actual treatment conditions, and the pad fluid volume needs to be optimized automatically. To this end, it is important to analyze the influences of the injection rate and fluid rheological parameters on the optimization results. Furthermore, the reservoir permeability has a strong influence on the optimization results, which will be also discussed in this section. Although the discussions used the specific data above, they are quite general and very useful for the field optimization design especially using the proposed methods in this paper. Moreover, the findings in this section are consistent well with the established knowledge. They are presented here for completeness.

Injection Rate
Based on the above data, we consider several injection rates: 4, 5, 6, 7, 8 and 9 m 3 /min. Other parameters, including proppant pumping curve index b and fluid rheological parameters are obtained from Table 3 and held constant. Then, the curve describing the relationship between pad fluid volume and injection rate is obtained, as shown in Figure 10. To generate the same fracture with the same fluid and proppant pumping schedule, the required pad fluid volume decreases with increasing injection rate. Note that when the injection rate is small, the required pad fluid volume is excessively large, which is unfeasible in terms of treatment difficulty and cost. Therefore, it is suggested to select a larger injection rate, within the treatment limitations, in order to reduce the amount of pad fluid used.

Fluid's Apparent Viscosity
Based on the above data, we can take different fluid's apparent viscosities, such as 20 mPa·s (for instance, if K = 0.46 Pa·s n and n = 0.5), 30 mPa·s (K = 0.5 Pa·s n and n = 0.55), 58 mPa·s (K = 0.7 Pa·s n and n = 0.6), 80 mPa·s (K = 0.71 Pa·s n and n = 0.65), and 100 mPa·s (K = 0.79 Pa·s n and n = 0.67). Other parameters come from Table 3 and are held constant. Then, the curve describing the relationship between the pad fluid volume and the fluid's apparent viscosity is as shown in Figure 11. To generate the same fracture with the same injection rate and proppant pumping schedule, the required pad fluid volume decreases with increasing apparent viscosity of the fluid. Note that when the fluid's apparent viscosity is small, the required pad fluid volume is excessively large, which is

Fluid's Apparent Viscosity
Based on the above data, we can take different fluid's apparent viscosities, such as 20 mPa·s (for instance, if K = 0.46 Pa·s n and n = 0.5), 30 mPa·s (K = 0.5 Pa·s n and n = 0.55), 58 mPa·s (K = 0.7 Pa·s n and n = 0.6), 80 mPa·s (K = 0.71 Pa·s n and n = 0.65), and 100 mPa·s (K = 0.79 Pa·s n and n = 0.67). Other parameters come from Table 3 and are held constant. Then, the curve describing the relationship between the pad fluid volume and the fluid's apparent viscosity is as shown in Figure 11. To generate the same fracture with the same injection rate and proppant pumping schedule, the required pad fluid volume decreases with increasing apparent viscosity of the fluid. Note that when the fluid's apparent viscosity is small, the required pad fluid volume is excessively large, which is unfeasible in terms of treatment difficulty and cost. Therefore, it is concluded that a larger apparent fluid viscosity reduces the amount of pad fluid used. unfeasible in terms of treatment difficulty and cost. Therefore, it is concluded that a larger apparent fluid viscosity reduces the amount of pad fluid used. Figure 11. Curve of the relationship between pad fluid volume and its apparent viscosity.

Reservoir Permeability
Based on the above data, we now consider reservoir permeability values of: 0.001, 0.01, 0.1, 0.2, 0.35, 0.46 and 1 md. Figure 12 shows the optimal propped fracture half-lengths and widths for various reservoir permeabilities. With increases in reservoir permeability, the optimal fracture half-length decreases while the optimal fracture width increases. This is consistent with established knowledge. That is to say, for very low permeability reservoirs, long and narrow fractures should be created. For low-to-medium permeability reservoirs, short and wide fractures should be created. Figure 13 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different reservoir permeabilities. From the figure, we can see that with decreases in reservoir permeability, the maximum dimensionless productivity increases. Despite this, as actual production is proportional to reservoir permeability, for very low-permeability reservoirs, absolute production remains low. Reservoir permeability (md) optimal fracture half-length optimal fracture width Figure 11. Curve of the relationship between pad fluid volume and its apparent viscosity.

Reservoir Permeability
Based on the above data, we now consider reservoir permeability values of: 0.001, 0.01, 0.1, 0.2, 0.35, 0.46 and 1 md. Figure 12 shows the optimal propped fracture half-lengths and widths for various reservoir permeabilities. With increases in reservoir permeability, the optimal fracture half-length decreases while the optimal fracture width increases. This is consistent with established knowledge. That is to say, for very low permeability reservoirs, long and narrow fractures should be created. For low-to-medium permeability reservoirs, short and wide fractures should be created. Figure 13 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different reservoir permeabilities. From the figure, we can see that with decreases in reservoir permeability, the maximum dimensionless productivity increases. Despite this, as actual production is proportional to reservoir permeability, for very low-permeability reservoirs, absolute production remains low.
Energies 2018, 11, x FOR PEER REVIEW 13 of 23 unfeasible in terms of treatment difficulty and cost. Therefore, it is concluded that a larger apparent fluid viscosity reduces the amount of pad fluid used. Figure 11. Curve of the relationship between pad fluid volume and its apparent viscosity.

Reservoir Permeability
Based on the above data, we now consider reservoir permeability values of: 0.001, 0.01, 0.1, 0.2, 0.35, 0.46 and 1 md. Figure 12 shows the optimal propped fracture half-lengths and widths for various reservoir permeabilities. With increases in reservoir permeability, the optimal fracture half-length decreases while the optimal fracture width increases. This is consistent with established knowledge. That is to say, for very low permeability reservoirs, long and narrow fractures should be created. For low-to-medium permeability reservoirs, short and wide fractures should be created. Figure 13 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different reservoir permeabilities. From the figure, we can see that with decreases in reservoir permeability, the maximum dimensionless productivity increases. Despite this, as actual production is proportional to reservoir permeability, for very low-permeability reservoirs, absolute production remains low.  During the optimization for various reservoir permeabilities, other parameters, including the injected proppant volume on the ground, injection rate, and desired proppant concentration are taken from Table 3 and kept invariable. The pad fluid volume and the fluid rheological parameters are optimized automatically. Figure 14 shows the optimized pad fluid volume and the fluid's apparent viscosity for different reservoir permeabilities. The figure shows that with increases in reservoir permeability, the pad fluid volume increases and the fluid's apparent viscosity increases accordingly. Thus, it can be seen that for very low permeability reservoirs, low viscosity fluid should be used, and for low-to-medium permeability reservoirs, high viscosity fluid should be used. This is also consistent with established knowledge. During the optimization for various reservoir permeabilities, other parameters, including the injected proppant volume on the ground, injection rate, and desired proppant concentration are taken from Table 3 and kept invariable. The pad fluid volume and the fluid rheological parameters are optimized automatically. Figure 14 shows the optimized pad fluid volume and the fluid's apparent viscosity for different reservoir permeabilities. The figure shows that with increases in reservoir permeability, the pad fluid volume increases and the fluid's apparent viscosity increases accordingly. Thus, it can be seen that for very low permeability reservoirs, low viscosity fluid should be used, and for low-to-medium permeability reservoirs, high viscosity fluid should be used. This is also consistent with established knowledge.

Special Discussion
Besides the general parameters discussed above, two other important parameters are introduced in this paper for coupling the fracturing treatment optimization and fracture geometry optimization: desired proppant concentration in the fracture and injected proppant volume on the ground. The desired proppant concentration in the fracture is used to solve the optimal fracture half-length and width. It also affects the treatment parameters in the fracture propagation model to reach the optimal fracture dimensions. The determination of it will be discussed in this section.
The injected proppant volume on the ground is also a key factor for fracturing. It determines the optimal fracture dimensions. Generally, there are two methods for optimizing it. The first is to select the maximum proppant volume that can be injected by the field treatment equipment, under the condition that all the treatment parameters are optimized to obtain the maximum productivity index using the proposed method. The second is to determine the injected proppant volume so as to  During the optimization for various reservoir permeabilities, other parameters, including the injected proppant volume on the ground, injection rate, and desired proppant concentration are taken from Table 3 and kept invariable. The pad fluid volume and the fluid rheological parameters are optimized automatically. Figure 14 shows the optimized pad fluid volume and the fluid's apparent viscosity for different reservoir permeabilities. The figure shows that with increases in reservoir permeability, the pad fluid volume increases and the fluid's apparent viscosity increases accordingly. Thus, it can be seen that for very low permeability reservoirs, low viscosity fluid should be used, and for low-to-medium permeability reservoirs, high viscosity fluid should be used. This is also consistent with established knowledge.

Special Discussion
Besides the general parameters discussed above, two other important parameters are introduced in this paper for coupling the fracturing treatment optimization and fracture geometry optimization: desired proppant concentration in the fracture and injected proppant volume on the ground. The desired proppant concentration in the fracture is used to solve the optimal fracture half-length and width. It also affects the treatment parameters in the fracture propagation model to reach the optimal fracture dimensions. The determination of it will be discussed in this section.
The injected proppant volume on the ground is also a key factor for fracturing. It determines the optimal fracture dimensions. Generally, there are two methods for optimizing it. The first is to select the maximum proppant volume that can be injected by the field treatment equipment, under the condition that all the treatment parameters are optimized to obtain the maximum productivity index using the proposed method. The second is to determine the injected proppant volume so as to

Special Discussion
Besides the general parameters discussed above, two other important parameters are introduced in this paper for coupling the fracturing treatment optimization and fracture geometry optimization: desired proppant concentration in the fracture and injected proppant volume on the ground. The desired proppant concentration in the fracture is used to solve the optimal fracture half-length and width. It also affects the treatment parameters in the fracture propagation model to reach the optimal fracture dimensions. The determination of it will be discussed in this section.
The injected proppant volume on the ground is also a key factor for fracturing. It determines the optimal fracture dimensions. Generally, there are two methods for optimizing it. The first is to select the maximum proppant volume that can be injected by the field treatment equipment, under the condition that all the treatment parameters are optimized to obtain the maximum productivity index using the proposed method. The second is to determine the injected proppant volume so as to maximize the NPV. For illustration, the influence of the injected proppant volume on the treatment parameter optimization results will be discussed, and the injected proppant volume used for the treatment will be determined using the first method.

Desired Proppant Concentration in the Fracture
When the desired proppant concentration is preset, it can be reached through optimization of the treatment parameters. Once both the desired proppant concentration and the optimal fracture half-length are satisfied during the optimization process, the fracture width will be satisfied naturally. Usually, the desired proppant concentration can be obtained through the parameter optimization of the proppant pumping curve index b. Based on the above data, we consider several proppant concentrations in the fracture: 800, 900, 1000, 1100 and 1200 kg/m 3 . Other parameters, including the injected proppant volume on the ground, the injection rate, and fluid rheological parameters, are taken from Table 3 and held constant. Figure 15 shows the optimal propped fracture half-lengths and widths with different proppant concentrations. Figure 16 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index with different proppant concentrations. As shown by the figures, with increasing desired proppant concentration, both the optimal fracture half-length and width decrease, as does the maximum dimensionless productivity index.
Energies 2018, 11, x FOR PEER REVIEW 15 of 23 maximize the NPV. For illustration, the influence of the injected proppant volume on the treatment parameter optimization results will be discussed, and the injected proppant volume used for the treatment will be determined using the first method.

Desired Proppant Concentration in the Fracture
When the desired proppant concentration is preset, it can be reached through optimization of the treatment parameters. Once both the desired proppant concentration and the optimal fracture half-length are satisfied during the optimization process, the fracture width will be satisfied naturally. Usually, the desired proppant concentration can be obtained through the parameter optimization of the proppant pumping curve index b. Based on the above data, we consider several proppant concentrations in the fracture: 800, 900, 1000, 1100 and 1200 kg/m 3 . Other parameters, including the injected proppant volume on the ground, the injection rate, and fluid rheological parameters, are taken from Table 3 and held constant. Figure 15 shows the optimal propped fracture half-lengths and widths with different proppant concentrations. Figure 16 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index with different proppant concentrations. As shown by the figures, with increasing desired proppant concentration, both the optimal fracture half-length and width decrease, as does the maximum dimensionless productivity index.   Figure 17 shows optimized proppant pumping curves for various desired proppant concentrations. However, a stepped proppant pumping schedule, rather than a continuous one, is used in actual fracturing treatments. This is represented by the dashed line in Figure 17. The optimal fracture conductivity maximum productivity index maximize the NPV. For illustration, the influence of the injected proppant volume on the treatment parameter optimization results will be discussed, and the injected proppant volume used for the treatment will be determined using the first method.

Desired Proppant Concentration in the Fracture
When the desired proppant concentration is preset, it can be reached through optimization of the treatment parameters. Once both the desired proppant concentration and the optimal fracture half-length are satisfied during the optimization process, the fracture width will be satisfied naturally. Usually, the desired proppant concentration can be obtained through the parameter optimization of the proppant pumping curve index b. Based on the above data, we consider several proppant concentrations in the fracture: 800, 900, 1000, 1100 and 1200 kg/m 3 . Other parameters, including the injected proppant volume on the ground, the injection rate, and fluid rheological parameters, are taken from Table 3 and held constant. Figure 15 shows the optimal propped fracture half-lengths and widths with different proppant concentrations. Figure 16 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index with different proppant concentrations. As shown by the figures, with increasing desired proppant concentration, both the optimal fracture half-length and width decrease, as does the maximum dimensionless productivity index.  optimal fracture conductivity maximum productivity index  Figure 17 shows optimized proppant pumping curves for various desired proppant concentrations. However, a stepped proppant pumping schedule, rather than a continuous one, is used in actual fracturing treatments. This is represented by the dashed line in Figure 17. The optimized proppant pumping curve coefficients a, indexes b, and actual proppant pumping schedules for various desired proppant concentrations are listed in Table 4. As shown in the table, with increases of the desired proppant concentration, the initial sand ratio (the 1st stage sand ratio) increases and the differences between the sand ratios of subsequent stages decrease. This style is rather disadvantageous for use in actual treatments, because the sand-laden fluid with the initial sand ratio stays in the fracture for a relatively long time, and the proppant concentration will continuously increase due to fluid loss. In this situation, a sand plug is likely to occur. In addition, as shown in Figure 16, with increasing desired proppant concentration, the maximum dimensionless productivity index decreases. Therefore, a high desired proppant concentration is not advantageous.  Table 4. As shown in the table, with increases of the desired proppant concentration, the initial sand ratio (the 1st stage sand ratio) increases and the differences between the sand ratios of subsequent stages decrease. This style is rather disadvantageous for use in actual treatments, because the sand-laden fluid with the initial sand ratio stays in the fracture for a relatively long time, and the proppant concentration will continuously increase due to fluid loss. In this situation, a sand plug is likely to occur. In addition, as shown in Figure 16, with increasing desired proppant concentration, the maximum dimensionless productivity index decreases. Therefore, a high desired proppant concentration is not advantageous.  With the decreases in the desired proppant concentration, the initial sand ratio decreases and the differences between the sand ratios of subsequent stages increase. This style is advantageous for avoiding sand plugs. However, from Figure 18 note that the optimal pad fluid volume increases with decreasing desired proppant concentration. In particular, when the desired proppant concentration is relatively small, the pad fluid volume becomes very large, which is unfeasible in terms of treatment difficulty and cost. Based on the above analysis, there is a reasonable range for the desired proppant concentration in the fracture. In this paper, the desired proppant concentration is set to 1000 kg/m 3 , which provides a reference for actual designs.  With the decreases in the desired proppant concentration, the initial sand ratio decreases and the differences between the sand ratios of subsequent stages increase. This style is advantageous for avoiding sand plugs. However, from Figure 18 note that the optimal pad fluid volume increases with decreasing desired proppant concentration. In particular, when the desired proppant concentration is relatively small, the pad fluid volume becomes very large, which is unfeasible in terms of treatment difficulty and cost. Based on the above analysis, there is a reasonable range for the desired proppant concentration in the fracture. In this paper, the desired proppant concentration C s is set to 1000 kg/m 3 , which provides a reference for actual designs.

Injected Proppant Volume on the Ground
Based on the above data, we now consider several injected proppant volumes: 6, 12, 18, 24 and 30 m 3 . Figure 19 shows the optimal fracture half-lengths and widths for different injected proppant volumes. Figure 20 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different injected proppant volumes. The figures show that with increases of injected proppant volume, both the optimal fracture half-length and width increase, and the maximum dimensionless productivity index increases as well.

Injected Proppant Volume on the Ground
Based on the above data, we now consider several injected proppant volumes: 6, 12, 18, 24 and 30 m 3 . Figure 19 shows the optimal fracture half-lengths and widths for different injected proppant volumes. Figure 20 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different injected proppant volumes. The figures show that with increases of injected proppant volume, both the optimal fracture half-length and width increase, and the maximum dimensionless productivity index increases as well.

Injected Proppant Volume on the Ground
Based on the above data, we now consider several injected proppant volumes: 6, 12, 18, 24 and 30 m 3 . Figure 19 shows the optimal fracture half-lengths and widths for different injected proppant volumes. Figure 20 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different injected proppant volumes. The figures show that with increases of injected proppant volume, both the optimal fracture half-length and width increase, and the maximum dimensionless productivity index increases as well.

Injected Proppant Volume on the Ground
Based on the above data, we now consider several injected proppant volumes: 6, 12, 18, 24 and 30 m 3 . Figure 19 shows the optimal fracture half-lengths and widths for different injected proppant volumes. Figure 20 shows the optimal dimensionless fracture conductivity and the maximum dimensionless productivity index for different injected proppant volumes. The figures show that with increases of injected proppant volume, both the optimal fracture half-length and width increase, and the maximum dimensionless productivity index increases as well.  However, with increases in the injected proppant volume, the treatment cost and difficulty will both increase. For different injected proppant volumes, other parameters, including the desired proppant concentration, are taken from Table 3 and held constant. The pad fluid volume, fluid rheological parameters and injection rate are to be optimized next. Figures 21 and 22 show the optimized pad fluid volume, fluid apparent viscosity, and injection rate for different injected proppant volumes. From the figures, with increases of injected proppant volume, the pad fluid volume increases, the fluid apparent viscosity increases and the injection rate also increases. Thus, in spite of the increase in the dimensionless productivity index resulting from the increase in proppant volume, the treatment becomes more and more difficult. Therefore, the injected proppant volume should be controlled within a suitable range according to the treatment conditions in order to use it efficiently. However, with increases in the injected proppant volume, the treatment cost and difficulty will both increase. For different injected proppant volumes, other parameters, including the desired proppant concentration, are taken from Table 3 and held constant. The pad fluid volume, fluid rheological parameters and injection rate are to be optimized next. Figures 21 and 22 show the optimized pad fluid volume, fluid apparent viscosity, and injection rate for different injected proppant volumes. From the figures, with increases of injected proppant volume, the pad fluid volume increases, the fluid apparent viscosity increases and the injection rate also increases. Thus, in spite of the increase in the dimensionless productivity index resulting from the increase in proppant volume, the treatment becomes more and more difficult. Therefore, the injected proppant volume should be controlled within a suitable range according to the treatment conditions in order to use it efficiently.

Conclusions
This paper presented a method to couple the fracturing treatment optimization and fracture geometry optimization in order to maximize the dimensionless productivity index. By introducing a desired proppant concentration in the fracture, the optimal fracture half-length and width can be solved based on UFD method and an iterative approach. Moreover, a rapid semi-analytical fracture propagation model was proposed in combination with an interval search method to optimize fracturing treatment parameters in order to reach the optimal fracture dimensions. Based on the case study and analyses contained in this paper, the following conclusions can be drawn: However, with increases in the injected proppant volume, the treatment cost and difficulty will both increase. For different injected proppant volumes, other parameters, including the desired proppant concentration, are taken from Table 3 and held constant. The pad fluid volume, fluid rheological parameters and injection rate are to be optimized next. Figures 21 and 22 show the optimized pad fluid volume, fluid apparent viscosity, and injection rate for different injected proppant volumes. From the figures, with increases of injected proppant volume, the pad fluid volume increases, the fluid apparent viscosity increases and the injection rate also increases. Thus, in spite of the increase in the dimensionless productivity index resulting from the increase in proppant volume, the treatment becomes more and more difficult. Therefore, the injected proppant volume should be controlled within a suitable range according to the treatment conditions in order to use it efficiently.

Conclusions
This paper presented a method to couple the fracturing treatment optimization and fracture geometry optimization in order to maximize the dimensionless productivity index. By introducing a desired proppant concentration in the fracture, the optimal fracture half-length and width can be solved based on UFD method and an iterative approach. Moreover, a rapid semi-analytical fracture propagation model was proposed in combination with an interval search method to optimize fracturing treatment parameters in order to reach the optimal fracture dimensions. Based on the case study and analyses contained in this paper, the following conclusions can be drawn:

Conclusions
This paper presented a method to couple the fracturing treatment optimization and fracture geometry optimization in order to maximize the dimensionless productivity index. By introducing a desired proppant concentration in the fracture, the optimal fracture half-length and width can be solved based on UFD method and an iterative approach. Moreover, a rapid semi-analytical fracture propagation model was proposed in combination with an interval search method to optimize fracturing treatment parameters in order to reach the optimal fracture dimensions. Based on the case study and analyses contained in this paper, the following conclusions can be drawn: (1) Through the semi-analytical fracture propagation model and the treatment optimization method, the desired proppant concentration in the fracture can be achieved by optimizing the proppant pumping curve index b. The optimal fracture half-length can be achieved by optimizing parameters such as pad fluid volume, injection rate, and fluid rheological parameters. Once both the desired proppant concentration and optimal fracture half-length are achieved, the optimal fracture width is also achieved naturally. (2) In order to obtain the optimal fracture dimensions, the treatment parameters are not unique.
The optimal treatment parameters can be determined according to both the actual treatment conditions in the field and the optimization procedure. The prior empirical knowledge about the actual fracturing treatment can also make the optimization more efficiently.  Using the shape factor, the equivalent proppant number can be written as: The maximum dimensionless productivity index and the optimal dimensionless fracture conductivity in rectangular reservoirs can be described by the following equations. The related constants are from Table A2.