Optimization of Mathematical Function-Shaped Fracture Distribution Patterns for Multi-Stage Fractured Horizontal Wells

: A conventional oil and gas well does not have a natural production capacity, which necessitates a hydraulic fracturing operation. The effectiveness of the fracturing directly impacts the economic beneﬁt of a single well. Among the various parameters, including fracture spacing, fracture width, and conductivity, fracture half-length is one of the main inﬂuencing factors on the productivity of horizontal wells. For conventional homogeneous reservoirs, research mainly focuses on fracture patterns with equal fracture lengths. However, in actual production processes, due to mutual interference and the superimposition of drainage areas between fractures, the production distribution of each fracture is non-uniform. Typical fracture distribution patterns mainly include uniform, staggered, dumbbell, and spindle. While many believe that the dumbbell-shaped fracture distribution pattern has the best effect, there has been no quantitative study on the length of each fracture under the dumbbell-shaped pattern. Based on this, this paper proposes a modeling approach for function-shaped fracture distribution that takes advantage of the high production of edge fractures and the low output of middle fractures in horizontal wells. The inﬂuence of this approach on production capacity is studied. Constant, linear, and parabolic functions are used to establish the relationship between fracture position and fracture half-length, optimizing the fracture distribution function to achieve the best production effect. This method can guide the horizontal well fracture distribution in the block to maximize productivity. The results show that the parabolic function-shaped model is better than the linear function-shaped model and the constant function-shaped model is the least effective. The research presented in this paper offers a new idea for optimizing on-site fracturing plans. It utilizes mathematical expressions to describe the parameters that affect productivity, which provides valuable guidance for designing multi-stage fractured horizontal wells in the ﬁeld. In the future, this research will be extended by exploring the optimal fracture distribution function under different formation conditions.


Introduction
Currently, there is an abundance of unconventional oil and gas resources in China; however, their low porosity, low permeability, and strong heterogeneity [1,2] make them challenging to exploit. Therefore, the technology of multi-stage fractured horizontal wells should be adopted to improve the productivity of single wells [3][4][5][6][7][8].
Both domestic and foreign experts have identified many factors that affect horizontal well productivity, with fracture being the most critical factor [9]. However, many scholars use electric model experiments and numerical simulations to optimize fracture parameters without considering the effect of fracture placement on horizontal well production [10].
In order to improve the fracturing effect, optimizing the reasonable fracture parameters is necessary; therefore, the optimization of fracture distribution patterns has an essential impact on improving productivity after fracturing [11]. Currently, the common fracture distribution patterns mainly include uniform-type, staggered-type, dumbbell-type, and spindle-type [12][13][14].
Since 2002, the multi-stage fracturing technology of horizontal wells [15] has been widely used in companies and oilfield enterprises, with the development of this technology becoming more and more mature [16].
In 2006, Gao Haihong et al. [17] used electrical simulation experiments to optimize the number of fractures, fracture length, horizontal wellbore length, the angle between the horizontal wellbore with fractures, and fracture spacing. However, the experiment was carried out under the condition of an equal fracture half-length.
In 2008, Mayerhofer M J et al. [18]. presented a correlation between the total amount of fracturing fluid and the total length of the fracture network in five Barnett Shale Wells; they proved that the greater the total amount of fracturing fluid, the longer the length of the fracture network.
In 2010, Zhou Linbo et al. [19] proposed a fracturing scheme with unequal fracture half-length to optimize fracture parameters under the overall fracturing development mode of the diamond-shaped reverse nine-spot well pattern. This was an early consideration of the impact of unequal fracture half-length on productivity; however, it was proposed on the premise of injection-production well patterns.
In 2014, Pu Chunsheng et al. [13] studied the fracture distribution pattern and fracture parameters of a fractured horizontal well. They also optimized the impact of spindleshaped, equal-length, and dumbbell-shaped staged fracturing on productivity. In the same year, Du Xianfei et al. [20] used numerical simulation to optimize the fracture parameters and distribution pattern of the fractured horizontal well. Li Zhongxing et al. [21] studied the influence of different fracture distribution patterns and fracturing parameters on the productivity of fractured horizontal wells under different well-row and well-pattern conditions. Ren Long et al. [14] utilized numerical simulation methods and micro-seismic monitoring results to study the influence of different fracture distribution patterns on the productivity of horizontal wells based on single factors. However, these studies were also conducted under the condition of an injection-production well pattern model.
In 2016, He Jun et al. [22] reported that the fracture production of horizontal wells has a high production rate for edge fractures and is lower in the middle. It is suggested that the production of horizontal wells can be increased by adopting a staggered fracture arrangement, prolonging the half-length of edge fractures and appropriately reducing the half-length in the middle. However, the research did not follow the principle of keeping the total fracture length constant when studying unequal half-length fracture distribution patterns.
In 2022, Deng Haiyang et al. [23] presented a new automatic integrated optimization algorithm that utilizes the SPSA nested binary search algorithm to optimize the discrete and continuous fracture parameters. They concluded that the optimal fracture arrangement presents a spindle-shaped distribution. However, it should be noted that they did not establish a functional relationship between fracture half-length and location.
Based on the current research status of fracture parameter optimization at home and abroad, it is commonly believed that the spacing between adjacent fractures is equal and the length of all fractures is the same in the homogeneous reservoir [24]. Although some research has also considered the impact of unequal spacing and unequal fracture length on the productivity of fractured horizontal wells, most of these studies analyze the basis of well networks [13,14,[19][20][21]. In earlier studies on the optimization of unequal fracture halflength in a single well, a dumbbell pattern was found to have the best effect [13]. However, there are too many fracture distribution schemes for the dumbbell patterns, making it difficult to determine the best one. Currently, there has been no quantitative study on each fracture length in the dumbbell model; additionally, no functional relationship between fracture half-length and location has been established. Therefore, this paper proposes a function-shaped fracture distribution pattern to study the effect of equal spacing and unequal length on productivity, taking into account the characteristics of the high production of edge fractures and the low production of middle fractures. The fracture patterns were constructed using a constant, linear function and parabolic function to establish the relationship between fracture half-length and fracture location. Among them, the linear and parabolic function-shaped patterns are variants of the dumbbell-shaped fracture distribution pattern. More kinds of dumbbell-shaped mathematical functions can be studied later. Numerical simulation models of multi-stage fracturing in horizontal wells under different fracture distribution patterns were established to optimize the fracture distribution scheme and function. The research presented in this paper offers a new idea for optimizing on-site fracturing plans. It utilizes mathematical expressions to describe the parameters that affect productivity, which provides valuable guidance for designing multi-stage fractured horizontal wells in the field.
This study mainly consists of six steps: (1) constructing the coordinate system and establishing the relationship between fracture half-length and fracture location; (2) collecting the parameters required for modeling; (3) establishing fracturing numerical models under different function-shaped fracturing schemes and conducting model verification; (4) simulating the cumulative output under different schemes; (5) obtaining the optimal function and the optimal function coefficient by comparing the simulation results; and (6) summarizing to reach a conclusion.

Basic Governing Equation
In this paper, the productivity of fractured horizontal wells under function-shaped fracture distribution patterns was analyzed using the WellWhiz4.0 software. The software utilizes the IMEX black oil model [25] in the CMG classic reservoir simulator as its computing core, which effectively solves the problems of simulating multi-phase flow and non-Darcy effects [26,27]. Additionally, the software uses the local mesh refinement method [28] to characterize fracture properties and simulates dual media, which allows for the analysis and evaluation of the potential stimulation effects of different fracturing plans.
During the process of fluid flow in a dual media reservoir, fluid exchange occurs between the matrix and fracture. The matrix-fracture transfer term is a physical quantity used to describe this fluid exchange, which reflects the ability of the fluid in the matrix to transfer into the fracture. Its value depends on the ratio of matrix permeability to fracture permeability and is affected by the stimulation of the matrix modified by artificial fractures. The larger the ratio of matrix permeability to fracture permeability or the greater the fracture density, the greater the matrix-fracture transfer term.
In the software model, the flow of the matrix with different phases in the doubleporosity or double-permeability reservoir into the fracture model is mainly divided into two situations. The first is the conventional model, which does not consider the gravity effect; the matrix fracture transfer is calculated using either Warren and Roots' method [29] or Gilman and Kazemi's [30] approach. The second is the complete phase separation model, which considers the gravity effect.
(a) Conventional model Normally, the matrix fracture transfer is calculated following either Warren and Roots' method [29] or Gilman and Kazemi's approach [30]. Based on the Warren-Root model, considering the pseudo-steady flow from the matrix into the fracture and the stress sensitivity of fracture, a horizontal well seepage model is established to provide guidance for the analysis of horizontal well pressure distribution.
The matrix-fracture transfer term τ jmf for phase j would be calculated as: where, P g = P o + P cog (2) K r is the relative permeability of phase j; meanwhile, K x , K y , K z are absolute permeability in the x, y, and z directions.
V is the grid block volume, ρ j is the reservoir condition density of phase j; meanwhile, L x , L y, and L z are the fracture spacing in the x, y, and z directions, respectively.
The matrix and fracture blocks are assumed to be at the same depth; thus, no gravity term is included above. Generally, the matrix capillary pressure is much larger than the fracture capillary pressure; the above expression for the transfer term τ jmf cannot correctly model the gravity-drainage process.

(b) Complete phase separation model
The complete phase separation model considers the gravity effect, that is, the gravity drainage process needs to be simulated. This assumes complete gravity segregation of the oil, water, and gas phases as follows: where, h is the height of the matrix element in the direction of gravity (h is normally the fracture height unless the block thickness is less than the fracture height, in which case h is the block thickness). The f subscript represents the fracture property and the m subscript represents the matrix property.
Note that normally ∆γ m =∆γ f and the additional gravity terms just account for the difference in the height of the fluid columns.

Mathematical Model of Function-Shaped Fracture Distribution Patterns
In order to address the characteristics of fractured horizontal wells that exhibit high production at both ends and low production in the middle, the length of the fractures at both ends was increased and the function of fracture location and fracture length of the horizontal wells was established using constant, linear, and parabolic functions. This enabled the research on the optimization of fracture distribution patterns under different schemes to be carried out.
This study employed the mid-point of the horizontal wellbore as the origin, the direction of the wellbore as the X-axis, and the direction of the vertical wellbore (that is, along the direction of the fracture) as the Y-axis, In this way, a rectangular coordinate system was established to study the relationship between the fracture and the fracture half-length. Along the X-axis, the coordinates from the origin to the toe of the horizontal wellbore are positive values and those from the origin to the root of the wellbore are negative values. On the Y-axis, the coordinates from the origin to the fracture tip on one side are negative values and on the other side are positive values. It is understood that the horizontal wellbore length is l m and the total number of fractures is N. When the function between the fracture location and the fracture half-length is a constant, linear function, or a parabolic function (represented by f 1 (x), f 2 (x), and f 3 (x), respectively), the best fracture distribution pattern for maximizing production in fractured horizontal wells is studied, as shown in Figure 1.

Mathematical Model of Function-Shaped Fracture Distribution Patterns
In order to address the characteristics of fractured horizontal wells that exhibit high production at both ends and low production in the middle, the length of the fractures at both ends was increased and the function of fracture location and fracture length of the horizontal wells was established using constant, linear, and parabolic functions. This enabled the research on the optimization of fracture distribution patterns under different schemes to be carried out.
This study employed the mid-point of the horizontal wellbore as the origin, the direction of the wellbore as the X-axis, and the direction of the vertical wellbore (that is, along the direction of the fracture) as the Y-axis, In this way, a rectangular coordinate system was established to study the relationship between the fracture and the fracture half-length. Along the X-axis, the coordinates from the origin to the toe of the horizontal wellbore are positive values and those from the origin to the root of the wellbore are negative values. On the Y-axis, the coordinates from the origin to the fracture tip on one side are negative values and on the other side are positive values. It is understood that the horizontal wellbore length is l m and the total number of fractures is N. When the function between the fracture location and the fracture half-length is a constant, linear function, or a parabolic function (represented by f1(x), f2(x), and f3(x), respectively), the best fracture distribution pattern for maximizing production in fractured horizontal wells is studied, as shown in Figure 1.    The three functions that define the fracture distribution pattern are all symmetrical about the Y-axis. According to the control variable method comparison principle, when comparing and optimizing the three function-shaped fracture distribution patterns, the sum of the fracture half-lengths must be equal. The equation relation follows: When L 1 = L 2 = L 3 , the total fracture half-length in different schemes can be kept constant. Here, L 1 , L 2 , and L 3 are the sum of the fracture half-lengths under the three schemes. The first scheme assumes a constant function-shaped fracture distribution pattern with the fracture half-length m:

Basic Parameters of the Reservoir Production Simulation Model
Before optimizing different fracture distribution patterns, some basic data need to be collected.
The target block is situated in the tight oil block of Jilin Oilfield, which is located in the southern Songliao Basin. This block has great potential for exploration [31]. The reservoir porosity in this block typically ranges from 7% to 12%, with an average of 10%. The permeability in the area generally falls between 0.01 mD and 0.2 mD, with a standard of 0.06 mD. The average thickness of the reservoir is approximately 15 m. This block is characterized by a low-saturation tight reservoir, which has low porosity and ultra-low permeability. The target well X was drilled to a depth of 3245 m and had a horizontal section length of 1300 m. The sandstone length was 1279 m, and the sandstone drilling rate was 98.38%.
Key parameters of the target well, including average water and oil saturation, average permeability, and average porosity, as well as the relative permeability curves (as shown in Figure 2), initial pressure, middle reservoir temperature, and some fracture parameters required for modeling, should be collected. The relevant data can be found in Table 1

Reservoir Production Simulation Model Establishment
Combining the basic data in Table 1, and based on the relationship between the fracture half-length and fracture position function established in Section 3.1, a numerical simulation model for multi-stage fractured horizontal wells under different fracture distribution patterns, such as constant-type, linear function-type, and parabolic function-type, was established. The fracture distribution schemes of the same function were optimized and studied and the optimal fracture distribution functions were optimized. Among them, the established parabolic function-type pattern numerical model is shown in Figure 3. After considering the target reservoir thickness and permeability, the optimal fracture half-length range has been determined to be between 200 m and 300 m. Therefore, the average half-length of a fracture is set between 200 m and 300 m to study the optimal distribution function.

Reservoir Production Simulation Model Establishment
Combining the basic data in Table 1, and based on the relationship between the fracture half-length and fracture position function established in Section 3.1, a numerical simulation model for multi-stage fractured horizontal wells under different fracture distribution patterns, such as constant-type, linear function-type, and parabolic function-type, was established. The fracture distribution schemes of the same function were optimized and studied and the optimal fracture distribution functions were optimized. Among them, the established parabolic function-type pattern numerical model is shown in Figure 3.

Reservoir Production Simulation Model Establishment
Combining the basic data in Table 1, and based on the relationship between the fracture half-length and fracture position function established in Section 3.1, a numerical simulation model for multi-stage fractured horizontal wells under different fracture distribution patterns, such as constant-type, linear function-type, and parabolic function-type, was established. The fracture distribution schemes of the same function were optimized and studied and the optimal fracture distribution functions were optimized. Among them, the established parabolic function-type pattern numerical model is shown in Figure 3.  The left part of Figure 3 shows the horizontal wells having equal fracture spacing but different fracture half-lengths. The X-axis represents the fracture position while the Y-axis represents the fracture half-length. There is a functional relationship between the fracture half-length and its position. When X > 0, the Y-value increases as X increases; whereas, when X < 0, the value of Y gradually decreases as X increases. This means that the fractures in the middle are short and are long at both ends. At the same time, the functions are symmetrical about the Y-axis. The right portion of Figure 3 depicts the fracture model, which comprises two primary sections: the main fracture and the stimulated area. A specified range of high permeability areas is set around the main fractures of the reservoir to simulate the effect of volume fracturing in tight oil reservoirs [32]. The red area indicates the reservoir matrix area that is not stimulated by volume fracturing and the blue area is the stimulated area around the main fracture where the permeability is higher than that of the reservoir matrix [32]. The stimulated area length is primarily provided by field volume fracturing microseismic data. Zhao Chaofeng et al. [33] conducted a series of micro-seismic tests in the Jilin exploration area; the results revealed that the ratio of stimulated area width to length varies from 9% to 29%, with an average of 20%. Based on these findings, a volume fracture model was established. The stimulated area length is approximately equal to the main fracture length and the stimulated area width is roughly 20% of the average fracture length.

Results and Discussion
Based on the numerical simulation model of a multi-stage fractured horizontal well under different fracture distribution patterns established in Section 3.3, we optimized the fracture distribution schemes of the same function and studied the optimal fracture distribution functions.

Model Verification
The target well was fractured for 13 stages with a designed fracture half-length of 200 m. The initial daily production was 8.1 m 3 /d. The cumulative production after one year of production was 1985 m 3 .
Based on this, as well as the data in Table 1, a fracturing model of an equal-stage spacing of 100 m and equal half-length of 200 m was established. The daily production during the initial simulation period was 35 m 3 /d; for the stable period, it was 8.5 m 3 /d. After one year of simulation, the cumulative production was 2150 m 3 .
As the reservoir porosity ranges from 7% to 12% and the permeability ranges from 0.01 to 0.2 mD, the reservoir is extremely heterogeneous. However, this heterogeneity is disregarded throughout the modeling procedure and the average reservoir porosity and permeability values are utilized in the model input. Therefore, there exists some discrepancy between the established homogeneous model and the actual model [24]. Upon comparing the model forecast data with the actual data, it has been discovered that the model forecast daily oil production error rate is approximately 4.9% and the model forecast 1-year cumulative oil production error rate is around 8.3%. These error rates are below 10%, which falls within an acceptable range. As a result, this model can be leveraged to guide the optimization of field schemes without considering the influence of some conditions, such as reservoir lithology and heterogeneity. Further studies will be carried out based on this model.

Optimization of Linear Function-Shaped Fracture Distribution Patterns
The formula of the linear function-shaped fracture distribution pattern is presented in Equation (11) There are a total of 20 plans with an average fracture half-length of 200 m (as depicted in Figure 4). Additionally, there are 20 schemes with an average fracture half-length of 300 m. Among these, Plan 1 is a constant function plan where the fracture half-lengths are equal and Plans 2-20 represent linear function schemes. There are a total of 20 plans with an average fracture half-length of 200 m (as depicted in Figure 4). Additionally, there are 20 schemes with an average fracture half-length of 300 m. Among these, Plan 1 is a constant function plan where the fracture half-lengths are equal and Plans 2-20 represent linear function schemes.   The linear function is symmetrical about the Y-axis; as the coefficient b (10 2 ) increases, the opening of the linear function decreases. Meanwhile, the slope of the linear function increases with regard to the difference in half-length between the fractures at both ends and in the middle of the horizontal well. The figure above shows that as coefficient b increases within a certain range, the cumulative oil production of multi-stage fractured horizontal wells after 5 years also increases. Once the coefficient b reaches a certain value, the cumulative oil production reaches its maximum value. However, with the increase of the b value beyond that point, cumulative oil production gradually decreases due to mutual interference and the super-imposition of the drainage area between fractures. The production of fractures in the middle of the horizontal wellbore is generally lower compared to the production at both ends due to non-uniform production distribution among fractures. Thus, to maximize production while keeping the total fracture half-length constant, the fracture half-length should be increased (the fracture half-length increases to a constant extent) at both ends and reduced in the middle. However, reducing the middle fracture half-length too much would result in a smaller oil drainage area, leading to a decrease in overall production. The figure above shows that when the average fracture halflength is 200 m, the optimal range for the parabolic function coefficient b (10 2 ) is between 40.24 and 46.43. Beyond that range, the trend of cumulative oil production increases at a slower rate. Similarly, when the average fracture half-length is 300 m, the optimal range for the parabolic function coefficient b (10 2 ) is between 46.43 and 52.62. Beyond that range, the trend of cumulative oil production increases at a slower rate as well.

ER REVIEW 10 of 17
The cumulative oil production results of horizontal wells over 5 years under the linear function-shaped fracture distribution patterns are illustrated in Figure 5 below, where the average half-length of fractures is 200 m or 300 m. The average fracture half-length of 300m

Optimization of Parabolic Function-Shaped Fracture Distribution Patterns
The parabolic function-shaped fracture distribution pattern is represented by Equation (12) and in the target well, there are 13 fractures with a fracture spacing of 100 m. The fracture position at the origin, X 7 , has a value of 0 while the corresponding fracture half-length, Y 7 , has a value of e. To establish the fracture distribution plan, two scenarios with average fracture half-lengths of 200 m and 300 m, respectively, are considered. Table 3 shows the specific values of d and e, where d represents the coefficient of the quadratic term of the parabolic function for average fracture half-lengths of 200 m or 300 m and e 1 and e 2 represent the intercepts of the function with the Y-axis.
There are a total of 20 plans with an average fracture half-length of 200 m (as shown in Figure 6 below). Additionally, there are 20 schemes with an average fracture half-length of 300 m. Among these, Plan 1 represents a constant function plan where the fracture half-lengths are equal while Plans 2-20 represent parabolic function schemes. The sum of the fracture half-length in each plan is equal to 2600 m. The production is studied under different schemes to optimize the fracture distribution pattern.        The parabolic function relation is ) ( +e =dx x f 2 3 ; as the coefficient d (10 5 ) increases, the opening of the parabola decreases and the degree of parabola bending increases to the half-length difference between the fractures at both ends and in the middle of the horizontal well. The figure above shows that as the coefficient d increases within a specific range, the cumulative oil production of the multi-stage fractured horizontal wells increases after 5 years. Once the coefficient d reaches a specific value, cumulative oil pro- The parabolic function relation is f 3 (x) = dx 2 + e ; as the coefficient d (10 5 ) increases, the opening of the parabola decreases and the degree of parabola bending increases to the half-length difference between the fractures at both ends and in the middle of the horizontal well. The figure above shows that as the coefficient d increases within a specific range, the cumulative oil production of the multi-stage fractured horizontal wells increases after 5 years. Once the coefficient d reaches a specific value, cumulative oil production hits the maximum. Beyond that, the cumulative oil production gradually decreases. Mutual interference and super-imposition of the drainage area between fractures lead to the nonuniform production distribution of each fracture.
The production of fractures in the middle of the horizontal wellbore is low, while the production is high at both ends. Hence, when the total fracture half-length is constant, increasing the fracture half-length at both ends (the closer to the wellbore, the more significant the increase in fracture half-length) and reducing it in the middle can maximize production. However, after reducing the middle fracture half-length to some extent, the oil drainage at the center of the fracture becomes smaller, decreasing overall production. The figure above indicates that when the average fracture half-length is 200 m, the optimal parabolic function coefficient d (10 5 ) is between 92.86 and 107.14. After that, the increasing trend of cumulative oil production slows down with the rise of the d value. When the average fracture half-length is 300 m, the optimal parabolic function coefficient d (10 5 ) is between 71.43 and 85.71. Again, after that, the increasing trend of cumulative oil production slows down with the rise of the d value. Figure 8 shows the cumulative oil production over 5 years under different fracture distribution patterns, with an average fracture length of 200 m; meanwhile, Figure 9 displays the cumulative oil production over 5 years under varying fracture distribution plans, averaging a fracture half-length of 300 m. The abscissa represents the plan number, and the ordinate denotes the cumulative oil production over 5 years. The green line represents the results of a constant function plan (i.e., equal-length fractures), the blue line represents the linear plan, and the orange line represents the parabolic plan. The sum of all of the fracture half-lengths in each plan is equal to 2600 m. distribution patterns, with an average fracture length of 200 m; meanwhile, Figure 9 dis-plays the cumulative oil production over 5 years under varying fracture distribution plans, averaging a fracture half-length of 300 m. The abscissa represents the plan number, and the ordinate denotes the cumulative oil production over 5 years. The green line represents the results of a constant function plan (i.e., equal-length fractures), the blue line represents the linear plan, and the orange line represents the parabolic plan. The sum of all of the fracture half-lengths in each plan is equal to 2600 m. Based on the figure above, it can be observed that among the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns with an average fracture half-length of 300 m, the best plan for maximizing cumulative output over 5 years is the parabolic function plan, followed by the linear and constant function plans. The most optimal number of parabolic function-shaped fracture distribution patterns falls within the 14-16 range. The optimal range for the coefficient d (10 5  107.14. However, as the number of plans increases, the increase in cumulative oil production gradually diminishes. Based on the figure above, it is evident that among the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns with an average fracture half-length of 300 m, the optimal plan for cumulative output in 5 years is the parabolic function plan, followed by the linear and constant function plans. The ideal number of optimal parabolic function-shaped fracture distribution patterns is roughly within the 11-13 range and the optimal range for the coefficient d (10 5 ) is between 71.43 and 85.71. However, as the number of plans increases, the increase in cumulative oil production gradually decreases.

Pressure Transfer Analysis
According to previous research results, the parabolic function-shaped fracture distribution pattern for fractured horizontal wells has the best effect. Therefore, in this study, we analyze the post-production pressure distribution under this fracture distribution pattern, as shown in Figure 10. This analysis aims to guide the optimization of the field fracturing scheme [34].  Based on the figure above, it can be observed that among the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns with an average fracture half-length of 300 m, the best plan for maximizing cumulative output over 5 years is the parabolic function plan, followed by the linear and constant function plans. The most optimal number of parabolic function-shaped fracture distribution patterns falls within the 14-16 range. The optimal range for the coefficient d (10 5 ) is between 92.86 and 107.14. However, as the number of plans increases, the increase in cumulative oil production gradually diminishes.
Based on the figure above, it is evident that among the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns with an average fracture halflength of 300 m, the optimal plan for cumulative output in 5 years is the parabolic function plan, followed by the linear and constant function plans. The ideal number of optimal parabolic function-shaped fracture distribution patterns is roughly within the 11-13 range Energies 2023, 16, 4987 14 of 16 and the optimal range for the coefficient d (10 5 ) is between 71.43 and 85.71. However, as the number of plans increases, the increase in cumulative oil production gradually decreases.

Pressure Transfer Analysis
According to previous research results, the parabolic function-shaped fracture distribution pattern for fractured horizontal wells has the best effect. Therefore, in this study, we analyze the post-production pressure distribution under this fracture distribution pattern, as shown in Figure 10. This analysis aims to guide the optimization of the field fracturing scheme [34]. 300 m (The optimal range is highlighted in red circle).
Based on the figure above, it is evident that among the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns with an average fracture half-length of 300 m, the optimal plan for cumulative output in 5 years is the parabolic function plan, followed by the linear and constant function plans. The ideal number of optimal parabolic function-shaped fracture distribution patterns is roughly within the 11-13 range and the optimal range for the coefficient d (10 5 ) is between 71.43 and 85.71. However, as the number of plans increases, the increase in cumulative oil production gradually decreases.

Pressure Transfer Analysis
According to previous research results, the parabolic function-shaped fracture distribution pattern for fractured horizontal wells has the best effect. Therefore, in this study, we analyze the post-production pressure distribution under this fracture distribution pattern, as shown in Figure 10. This analysis aims to guide the optimization of the field fracturing scheme [34].  Figure 10a shows that at the initial production stage, only the fluid inside the fracture flows due to the wellbore storage effect, which completely supplies the production currently. In the second stage, illustrated in Figure 10b, fluid-flow mainly occurs on both sides of the fracture, whereas linear flow occurs in the fracture and inner zone. In the third stage, demonstrated in Figure 10c, the pressure interference between the fractures starts and progressively increases over time. Additionally, due to the matrix's extremely low permeability, virtually closed boundaries emerge on both sides of the fractures. As the properties of the inner and outer zones vary, a virtually closed boundary forms between these zones. The boundary along the wellbore's side is mainly in the shape of a parabola, whereas the wellbore root and toe's sides take on an elliptical shape. In the fourth phase, depicted in Figure 10d, as the oil well's production proceeds, the pressure wave propagates outwards in the outer zone; the outer zone now presents a quasi-radial flow, gradually. The boundary shape of the outer zone eventually approaches that of the inner zone [35,36].

Conclusions
Under the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns, the optimal plan for cumulative output in 5 years is the parabolic function plan, followed by the linear and constant functions.
In the parabolic function-shaped fracture distribution patterns with an average fracture half-length of 200 m, the optimal range of coefficient d (10 5 ) is between 92.86 and 107.14. In the parabolic function-shaped fracture distribution patterns with an average  Figure 10a shows that at the initial production stage, only the fluid inside the fracture flows due to the wellbore storage effect, which completely supplies the production currently. In the second stage, illustrated in Figure 10b, fluid-flow mainly occurs on both sides of the fracture, whereas linear flow occurs in the fracture and inner zone. In the third stage, demonstrated in Figure 10c, the pressure interference between the fractures starts and progressively increases over time. Additionally, due to the matrix's extremely low permeability, virtually closed boundaries emerge on both sides of the fractures. As the properties of the inner and outer zones vary, a virtually closed boundary forms between these zones. The boundary along the wellbore's side is mainly in the shape of a parabola, whereas the wellbore root and toe's sides take on an elliptical shape. In the fourth phase, depicted in Figure 10d, as the oil well's production proceeds, the pressure wave propagates outwards in the outer zone; the outer zone now presents a quasi-radial flow, gradually. The boundary shape of the outer zone eventually approaches that of the inner zone [35,36].

Conclusions
Under the three schemes of constant, linear, and parabolic function-shaped fracture distribution patterns, the optimal plan for cumulative output in 5 years is the parabolic function plan, followed by the linear and constant functions.
In the parabolic function-shaped fracture distribution patterns with an average fracture half-length of 200 m, the optimal range of coefficient d (10 5 ) is between 92.86 and 107.14. In the parabolic function-shaped fracture distribution patterns with an average fracture halflength of 300 m, the optimal range of coefficient d (10 5 ) is between 71.43 and 85.71. After that, with the increase in the number of plans, the increase of cumulative oil production gradually decreased.
The research presented in this paper offers a new idea for optimizing on-site fracturing plans. It utilizes mathematical expressions to describe the parameters that affect productivity, which provides valuable guidance for designing multi-stage fractured horizontal wells in the field.
In the future, this research will be extended by exploring the optimal fracture distribution function under different formation conditions.  Data Availability Statement: Datasets related to this article can be found by connecting the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.