Evaluating essential features of proppant transport at engineering scales combining field measurements with machine learning algorithms

The behaviours of the particle settlement, stratified flow and inception of settled particles are essential features that determine the proppant transport in low-viscosity fracturing fluids. Although great efforts have been made to characterize these features, limited research work is performed at field scales. To test the laboratory outcomes, we propose a machine-learning-based workflow to evaluate the essential features using the measurements obtained from shale gas fracturing wells. Over 430,000 groups of fracturing data (1 s time interval) are collected and pre-processed to extract the particle settlement, stratified flow and inception features during fracturing operations. The GRU and SVM algorithms, trained by these features, are applied to predict fracturing pressure. Error analysis (the root mean squared error, RMSE) is carried out to compare the contributions of different features to the pressure prediction, based on which the features and the corresponding calculations are evalu- ated. Our result shows that the stratified-flow feature (fracture-level) possesses better interpretations for the proppant transport, in which the Bi-power model helps to produce the best predictions. The settlement and inception features (particle-level) perform better in cases where the pressure fluctuates significantly. The features characterize the state of proppant transport, based on which the development of subsurface fracture is also analyzed. Moreover, our analyses of the remaining errors in the pressure-ascending cases suggest that (1) an introduction of the alternate-injection process, and (2) the improved calculation of proppant transport in highly- filled fractures will be beneficial to both experimental observations and field applications.


Introduction
Hydraulic fracturing has become an important technique to enhance hydrocarbon recovery from unconventional gas resources, aiming to meet the growing demand for clean energy globally. To avoid fracture closure after the dissipation of the hydraulic injecting pressure, proppant injections are essential in the hydraulic fracturing process, and their effectiveness plays an important role in enhancing the stimulated reservoir volume (Barree and Conway, 1994;L. Fan, Thompson and Robinson, 2010;Nassir et al., 2014). To inject thousands of tons of proppant particles down into the induced fractures, a deep well in a shale gas reservoir (>4000 m) is usually operated under a wellhead pressure exceeding 100 MPa (Hou, Chang, Fu, Muhadasi andChen, 2019a, 2019b;Mao et al., 2021). Therefore, how to inject proppant particles under safe operating pressures is very challenging, especially with the application of low-viscosity slickwater (Liang et al., 2016). Proppant transport is, therefore, an essential research topic in hydraulic fracturing engineering (Economides and Nolte, 1989).
The behaviours of the proppant settlement, stratified flow and inception of settled particles in low-viscosity fracturing fluid have been characterized numerically based on experimental tests (Gadde et al., 2004;Wei et al., 2020;Zhao et al., 2019), which are the essential features of proppant transport. The recent trend of the proppant transport research is to bring in more realistic subsurface scenarios by replacing the single smooth-panel fracture with artificial-coarse fracture networks (Manchanda et al., 2020;Raki Sahai and Moghanloo, 2019;Tong and Mohanty, 2016). However, it is still challenging to simulate and characterize the realistic morphology (scale, tortuosity, branches et al.) of fractures numerically (Dahi-Taleghani and Olson, 2011). Moreover, the direct subsurface measurements and observations during fracturing operations are still limited. Many numerical models for calculating the proppant transport, therefore, are usually verified by laboratory experiments (Mack et al., 2014;Patankar et al., 2002;Raki Sahai and Moghanloo, 2019). However, the approaches for examining the numerical observations at a field scale are still in demand.
The links between indoor research and field practice may rely on exploring the key property features that control the proppant transports (Cai et al., 2017). Machine learning (ML) is one of the key techniques to perform the exploration at field-practical scales. Researchers built ML workflows to predict the fracturing pressure and forewarn the sand screen-out in real-time based on fracturing records (Ben et al., 2020;Hu et al., 2020). The simulation data is also applied to learn the pressure pattern and detect sand screen-out (Yu et al., 2020). These effects help to improve the safety of fracturing operations by combining data science and engineering knowledge. It is also important to understand how the key features, being tested in the lab, influent the proppant transport during real fracturing operations, which may, in turn, promote fundamental research.
In this study, we examine and compare essential features of proppant transport and their corresponding calculations by using a new workflow, where we introduce machine learning (ML) algorithms, including Support Vector Regression (Al-Anazi and Gates, 2010;Al-Mudhafar, 2017) and Gated Recurrent Units (L. Hou, Cheng, Wang, Ren and Geng, 2022;Sun et al., 2020). The ML method can process the field measurements for proppant transport evaluation directly without in-depth characteristics of the realistic fracture morphology, which may build a bridge between the fundamental research and field applications. Based on the data-driven approach, our study is aimed to 1) propose a new workflow to estimate the proppant transport at field engineering scales; and 2) better understand the essential features that control the proppant transport, which is valuable for both field engineering and basic research work.

Methodology
The field measurements of the shale gas fracturing treatments are collected and carefully pre-processed (splitting, trimming, and denoising) for training. The proppant transport features, specifically the velocity ratio and the height of the flowing layers within fractures, are initially calculated by several popular proppant transport models. The calculation outputs consisting of the features relevant to the proppant flow, as long as the other subsurface measurements, are then fed into the machine learning algorithms to predict downhole pressure. The predictions are further analyzed using the control variate method and error analyses to evaluate the proppant transport features and their corresponding calculations.

Data collection and preprocessing
55 stages, including over 430,000 groups, of fracturing measurements (in second) are collected from 10 shale gas wells, which are selected from 5 different platforms in the Sichuan basin, China (Table 1). The field measurements include the geological data (vertical and well depths), clustering data (stage length, cluster number and perforation number), and fracturing data (fluid and proppant types, pump rate, proppant concentration and wellhead pressure). Five of the ten wells are set as the training well (A 1 -E 1 ), of which 50 stages of fracturing data are pre-processed for training the machine learning models. Five testing stages are selected from the remaining five wells (A 2 -E 2 ), defined as testing wells. To constrain the effect of large spatial variation in geological uncertainty and formation properties on the predictions, each training well has its own testing well that is selected from the same platform. For instance, both Well A 1 and Well A 2 (neighbouring wells) are from Platform A, and so forth, as shown in Table 1. This is one of our strategies to eliminate interference factors of pressure variation and promote the influence of proppant transport.
The other strategy during the data preprocessing is to convert the wellhead pressure into the downhole pressure after the perforation hole (Appendix A), defined as the DPP. The conversion can rule out the potential effects of hydrostatic pressure and friction variations, leaving the proppant transport to control the fracture pressure fluctuation (Dontsov and Peirce, 2014;Willingham et al., 1993). Other denoising methods involve trimming the pressure at the beginning (when the fracture is created) and the end (pump-off) of the fracturing operation, repeating predictions and averaging the errors obtained from all the platforms (A -E), as long as applying two different machine learning algorithms.

Features for proppant transport
In general, we divide the proppant transport features into two categories by their scalesparticle level (Fig. 1 a) and fracture level ( Fig. 1  b), including the particle settling velocity (Gadde et al., 2004;Mack et al., 2014;McCabe et al., 1993;Richardson and Zaki, 1954;Yew and Weng, 2014), the critical velocity to restart the settled proppant (also used as the critical turning velocity in complex fractures) (Cao et al., 2006;Hou, Jiang, Li, Zeng andCheng, 2017a, 2017b;Hou, Jiang, Liu, et al., 2017a, 2017bRakshit Sahai, Miskimins, & Olson, 2014), the flowing layer height (H 1 ) (Hou et al., 2019a(Hou et al., , 2019bNovotny, 1977;Patankar et al., 2002;Jing Wang, Joseph, Patankar, Conway and Barree, 2003) and the equilibrium dune level (EDLthe dune height divided by fracture height) (Alotaibi and Miskimins, 2019). Based on the field pumping schedules, those features are further calculated by employing the Velocity, Settling, Bi-power, and EDL models to yield a group of independent variables, which is one of the inputs for ML models (Table 2). Details about the equations and their applications can be found in Appendix B. The representative features, particle settlement (υ settling /υ f ), inception (υ turning /υ f ) and stratified-flow (H 1 and EDL) behaviours, are selected and tabulated in Table 2. The particle settlement and inception are grouped because they are decomposition features of particle movements in vertical and horizontal directions. The selected features in Table 2 control the proppant transport in the low-viscosity fluid, based on which more comprehensive models coupling fracture propagation, fluid leak-off, etc. are derived (Barboza et al., 2021;Isah et al., 2021). Besides, the calculations (Appendix B) for the selected features are analytical, which is more calculational effective to pre-processing our datasets (over 430, 000 groups of measurements) compared with numerical solutions. Furthermore, the models in Table 2 are mainly derived from observations of experimental simulations (Appendix B). By evaluating the calculation outputs at field-practical scales, the experimental techniques may be improved in the aspects of equipment, parameters, methodology, measurements, etc.
During the calculation, the fracture width is the only unknown parameter that is presumably set to a value of 100 × d max (d max is the largest diameter of injected proppant) referring to the result of slant core drilling through a stimulated shale reservoir (Elliott and Gale, 2018). For an alternate pumping schedule (injecting pure fluid and slurry alternatively), the results of the Velocity model are discrete and are all treated as zeros as pure fluid is injected. As shown in Table 2, the fluid typeμ, pump rate -Q, proppant concentration -C and proppant typed) are jointly fed into the ML models for the prediction (dependent variable: downhole pressure after perforation -DPP) and error analyses (the Root Mean Square Error -RMSE). The non-numeric variables (fluid and proppant types) are replaced with the values of fluid viscosity and averaged proppant diameters, respectively. For comparison purposes, the original field measurements alone are directly used to train the ML models to predict the DPP, defined as the unprocessed DPP (Table 2).

Machine learning models and workflow
To constrain the high variance and boost the prediction performance, two different machine learning algorithms are applied for training and predicting. The Support Vector Regression (SVR) model, with a Radial Basis Function (RBF) kernel, is capable of both linear and non-linear regression (Al-Anazi and Gates, 2010), being of memory efficiency, and performing well in various petroleum engineering applications (Goel et al., 2017;Guo et al., 2018). Furthermore, we apply Gated Recurrent Units (GRU) to the same datasets. The GRU is a deep learning algorithm designed for extracting information from time-sequence data. In GRU models, the current state and prediction can be influenced by the preceding state and will affect the following prediction at the next time step as well, making the GRU models appropriate for handling continuous hydraulic fracturing data (Sun et al., 2020;Jinjiang Wang, Yan, Li, Gao and Zhao et al., 2019). According to previous modelling experience (Cho et al., 2014;D. Fan et al., 2021), a three-layer (including the output layer) GRU model is constructed with the activation function of ReLU. The dropout (0.2) layers are applied to avoid the overfitting of the model (Gal and Ghahramani, 2015). The Adam optimizer is used in the model with a learning rate starting at 0.0005 (Kingma and Ba, 2014).
Using Platform A as an example, our workflow for the data processing is shown in Fig. 2. Model i represents one of the four proppant transport models given in Table 2. The reference is the DPP converted directly from the surface pressure records, and Prediction i is the predicted DPP based on Model i. The pressure prediction is made for each platform (A to E), and for each prediction, a new GRU model and SVR model are created and trained respectively. Eventually, the prediction errors for each platform are averaged to evaluate the performance of the selected features and the corresponding models. Fig. 1. Proppant transport features at (a) particle-level (particle settlement and inception); (b) fracture-level (H 1 -the height of the flowing layer) (Hou, Jiang, Liu, et al., 2017a, 2017bPatankar et al., 2002).  Fig. 2. Schematics of the data processing workflow using Platform A as an example.

Evaluation of proppant transport features at engineering scales
The pressure predictions by GRU for Well A 2 at Platform A are plotted in Fig. 3, using as an example of the evaluation results. In general, the predicted DPP curves corresponding to the proppant transport models match the reference DPP curve much better than the unprocessed curve, indicating that the introduction of the proppant transport features improves the prediction accuracy and reduces the variance (Fig. 3). The proppant concentration during the fracturing operation is also presented (green solid line in Fig. 3), demonstrating that the proppant-injection-induced pressure fluctuations influence the variation in pressure predominantly. To denoise the pressure variation induced by the injection, the pressure data ranging between the time of 2000 s (the period of the fracture initiation and propagation) and the time of 8000 s (the period of the fracture closure and fluid diffusion after pump-off) are only used for error analyses (the region between two vertical grey dash lines in Fig. 3).
Four more ML predictions are carried out for Platforms B − E. The errors based on different proppant transport models are averaged to compare the performances of the GRU and SVM algorithms, as shown in Fig. 4. Generally, the two algorithms exert close performances according to the generated errors. The GRU-based workflow produces smaller errors for Wells A 2 , B 2 and C 2 . Therefore, the predictions based on the GRU model are selected for further investigations. Besides, the errors for Wells A 2 and B 2 are significantly smaller than in the rest of the cases. We divided the cases into two groups in the later analysesthe small error group (Wells A 2 and B 2 ) and the large error group (Wells C 2 , D 2 , and E 2 ).
The detailed errors between the reference and the ML prediction for each testing well (A 2 − E 2 ) based on the GRU algorithm are summarized in Table 3. By comparing the averaged RMSE, we find that the DPP predictions are enhanced by introducing the Velocity, Settling, and Bipower models, in which the stratified-flow feature performs better than the settlement and inception features. The Bi-power model helps yield the best DPP predictions, followed by the Settling model. The introduction of the EDL model promotes low RMSEs for Wells A 2 , B 2 , and C 2 , whereas leads to large prediction errors for Wells D 2 and E 2 . The performance of the Velocity model is probably limited by the simplification employed under the pure-fluid condition. However, exceptions are observed for Wells D 2 and E 2 , where the Velocity model helps to produce smaller errors than the Settling model does.

Error variance in different cases
According to Table 3, we plot the DDP curves predicted by the GRU algorithm for further investigation, as shown in Figs. 5 and 6. Each dash curve corresponds to a different model used to calculate the input features, including the Velocity (particle-level feature), Settling and Bipower (fracture-level feature) models. As the different calculations are applied for the pure-fluid and slurry injections (Appendix B) when preprocessing the input data, the predicted curves derived from the alternative injection schedule are relatively discrete. The predictions, therefore, fluctuate around the reference at a frequency following the oftenness of the injection switching.
For the small error group (Fig. 5), a relatively flat trend of DPP along with a constant pumping rate can be found throughout the injection treatment in Wells A 2 and B 2 . The effect of proppant transport on pressure variation is moderate, suggesting that the fracture volume may be sufficient for the current proppant injection rate. Besides, the predicted pressure curve is unsmooth and exerts vertical climbing and jumping between slugs (the alternate from pure fracturing fluid into proppant slurry).
In contrast, an ascending pressure trend (red solid curves in Fig. 6) can be found in the large error group (Wells C 2 , D 2 and E 2 ). Compared with the proppant concentrations in Wells A 2 and B 2 (~20%), the proppant concentrations for Well C 2 , D 2 , and E 2 are all under 10%, indicating that their DPPs are relatively sensitive to the proppant transport. In Well C 2 , D 2 , and E 2 , the proppant particles are likely driven into fractures possessing insufficient volume, where the continuously injected proppant may accumulate, and then block the flowing pass (as shown in Fig. 7), resulting in a gradual increase in flowing friction, reflected by the ascending operation pressure (Zhang et al., 2017). Besides, the performance of the Velocity model is unexpected and even better than the Settling and Bi-power models, as shown in Fig. 6 (b) and (c). Integrating the severe fluctuations of fracturing pressure into account, the underground fracture in Fig. 6 may be more complex than that in Fig. 5. The proppant may be transported in fracture networks. Therefore, the Velocity model, calculating the critical condition of proppant turning from the main fracture into the minor fracture, produces better predictions.

Discussion
Our comparison among pressure predictions derived from different proppant transport features exhibits that the stratified-flow feature (H 1 , calculated by the Bi-power and Settling models) can improve the pressure prediction considerably. The unstable predictions produced by the EDL model may be attributed to the limited application range of the empirical equations used in the model (Appendix B). The Velocity model characterizes the particle flowing feature during the slurry injection but fails to take into account the scenario when the pure fluid is injected, which results in a restrained improvement in the pressure prediction.
However, relatively large prediction errors exist for cases of Wells C 2 , D 2 , and E 2 (Fig. 6), which are non-negligible and likely attributed to the following factors: (i) Effect of the injection alternation -The proppant transport models are featured by taking into account the accumulation of the proppant dune being in an equilibrium state. Hence, the prediction curves in Figs. 5 and 6 are relatively discrete under an alternate injection schedule. However, the time interval (around 3-5 min) between the injecting alternation may be insufficient to allow the proppant dune to reach the equilibrium state (Yew and Weng, 2014), thus resulting in the discrete pattern of the predicted curve and vertical pressure variations between slugs. It is likely that the transition state of the proppant dune between two alternating injections influences the pressure substantially, also contributing to the errors for the cases with pressure-ascending trends (Fig. 6). This feature describing the transition state, however, is not reflected by any model we evaluate in this study. (ii) Fracture propagation during proppant injection -The introduction of the velocity feature (describing the critical flow condition that drives the proppant to turn into branching fractures) enhances the prediction performance for Wells D 2 and E 2 (Fig. 6), implying that more complex fracture networks may be generated. The amplitudes of the pressure fluctuations shown in Fig. 6 are broadly larger than those observed in Fig. 5, which may be attributed to the development of branching or minor fractures. The random fracture propagation may cause unexpected pressure variation and thus extra prediction errors (Fig. 6). (iii) Proppant transport in highly-filled fractures -According to the discrepancies between reference and predicted DDP curves, the largest errors emerge at the beginning and end of proppant injections (Fig. 6). Initially, the fracture is underdeveloped with limited volumes. At the end of operations, a large volume of proppant has been injected into the fractures. The similarity of these two conditions is that the fracture is highly filled due to the relative volumes of fractures and proppant. However, few relevant research works can be found during our literature review.
The highly-filled-fracture operating condition may be critical for pressure-sensitive cases and the sand screen-out, thus deserving more studies.
Therefore, we suggest investigating further (1) the evolution of proppant dune based on a staged pumping schedule, and (2) a better assessment of proppant transport in fracture networks and highly-filled Table 3 Summary of the RMSE based on each proppant transport model. fractures.

Conclusions
In this study, we propose a machine-learning-based (GRU and SVM) workflow to process field measurements collected from shale gas fracturing to assess the essential proppant transport features and their corresponding calculations at field-practical scales. The features of particle settlement, stratified flow and inception are evaluated indirectly via algorithm training, fracturing pressure prediction and error analysis. The new workflow paves a path to potentially establish a link between laboratory work and engineering practices. The feature analysis improves the awareness of underground proppant transport in engineering scales, which may provide a complement to numerical and experimental simulations. The main conclusions are generalized as follows: (1) The Bi-power model produces the most accurate prediction of fracturing pressures, followed by the Settling model (stratifiedflow feature) and Velocity model (settlement and inception features). The introduction of the essential features enhances the pressure predictions for the cases where relatively flat trends of pressure evolution (under constant pump rates) are present, implying that the proppant is likely injected into a sufficient volume of fracture, comparable to the conditions simulated by the lab research.
(2) For the cases, where an ascending trend of pressure is shown throughout the proppant injections, all features and calculations bring in relatively large prediction errors. However, the Velocity model (characterizing the critical flow velocity that drives proppant to turn into branching fractures) helps to yield less prediction error in these cases, indicating that the proppant may be transported into un-fully-developed fracture networks. The underground fracture may be more complex in the pressureascending cases according to the more severely fluctuated pressure that may be induced by the generations of branching or minor fractures. (3) The existing errors in the pressure-ascending cases can be improved by enhancing the accuracy of feature calculating models, where the alternate injection schedule and the random propagation of fracture are still missing. Based on the feature tests at field scales, we suggest that the evolution of proppant dune during an alternate pumping schedule may play a critical  role in pressure evolutions. Moreover, the proppant transport in fracture networks and highly-filled fractures should be defined more accurately for the pressure-sensitive operations to mitigate the operating risk and improve the proppant injection.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Data availability
The authors do not have permission to share data.

Acknowledgement
This research has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No. 846775.

Appendix A. Calculations used for pressure conversion
The original wellhead pressure (or surface pressure) is converted into downhole pressure to remove noises of flowing resistance and hydrostatic pressure variation caused by changes of pump rate, fluid type, proppant type and concentration. P downhole = P wellhead + P statics − P pipeloss − P perforation (I) The hydrostatic pressure (P statics ) is calculated by the vertical well depth (h v ) The friction loss of the wellbore (P pipeloss ) is estimated by the Darcy-Weisbach equation where L is the wellbore length from its wellhead to the fracturing stage, m; υ s is the flowing rate of slurry in the wellbore, m/s; μ s is the slurry viscosity, Pa⋅s. The slurry viscosity is calculated by (Dontsov and Peirce, 2014) μ s = μ f where C m is the maximum proppant concentration and is assigned a value of 0.585. The pressure drop through the perforation hole is estimated based on the hydraulic and perforation parameters (Willingham et al., 1993) where d h is the diameter of the perforation hole, m; C p is the coefficient of discharge and is 0.6-0.95 for slurry; n is the number of the opening perforation hole. According to the mini-fracturing test, around half of the designed perforation holes will be opened.

Velocity model
The slickwater, widely used for massive hydraulic fracturing, could only suspend the proppant for minutes during the fracturing operation (Yew and Weng, 2014). Therefore, the proppant settling velocity, one of the most fundamental issues, is given by (Mack et al., 2014;McCabe et al., 1993) where υ settling is the proppant settling velocity, m/s; ρ p and ρ f are densities of proppant and fracturing fluid, respectively, kg/m 3 ; μ f is the fluid viscosity, Pa⋅s; d is the averaged diameter of proppant, m. Settling in a fracture, proppant is slowed down by fracture walls and interactions between particles, which can be modified by (Gadde et al., 2004;Richardson and Zaki, 1954) υ ′ settling = υ settling ( 2.37C 2 − 3.08C + 1 ) where C is the volume fraction of proppant; w is the fracture width, m. For a complex fracture network, the dragging of the carrying fluid is one of the most important motivations to drive the proppant turn into branch fractures (Rakshit Sahai et al., 2014). A minimum flowing rate of the slurry (υ turning ) is required for proppant turning, and is estimated by the proppant restarting condition (Cao et al., 2006;Hou, Jiang, Li, et al., 2017a, 2017b The particle movements could reflect the proppant transporting by evaluating the ratio between particle and fluid velocities (Hou, Jiang, Li, et al., 2017a, 2017b

Settling model
The proppant is tending to form an equilibrium dune in low-viscosity fluids under constant injection conditions (Hou et al., 2019a(Hou et al., , 2019b, as shown in Fig. 1 (b). The height of the flowing layer above the dune is a core parameter that evaluates the proppant transport. It could be estimated by the settling model expressed as (Novotny, 1977) H 1 = 16.67Q wυ eq (5) where H 1 is the height of the flowing layer, m; Q is the pump rate, m 3 /s; υ eq is the flowing rate when the particle settling and restarting reach equilibrium, and is calculated by where ρ s is the density of the slurry, m 3 /s; φ is the porosity of the proppant dune.
The slurry and pure fluid are usually injected alternately for shale gas fracturing. For the slurry condition, H 1 is estimated by Eq. (5). When pure fluid is injected, H 1 is assumed to be consistent with prior H 1 after the last slurry was injected. This process is carried out manually during the data preprocessing.

Bi-power model
The Bi-power law correlations are proposed to directly calculate the height of the flowing layer (Jing Wang et al., 2003), which is defined as where R f , R p , R G and λ are calculated by where Q f is the pump rate of fracturing fluid, m 3 /s; Q p is the pump rate of proppant, m 3 /s. There is a special condition when pure fluid is injected to push the injected proppant deeply into the fracture. The pure fluid may rebalance the proppant dune, which could be calculated by (Patankar et al., 2002) H 1 =

EDL model
A similar empirical power-law formula for dune height has been derived based on a series of sand-carrying experiments (Alotaibi and Miskimins, 2019). The equilibrium dune level (EDL) is defined as the dune height divided by fracture height, which is proposed as EDL = − 0.003496d − 0.3277 [ υ eq (0.7749 − 0.1955 ) C ] + 0.9901d − 0.02667 (10)