Abstract

We explore the evolution of friction and permeability of a propped fracture under shear. We examine the effects of normal stress, proppant thickness, proppant size, and fracture wall texture on the frictional and transport response of proppant packs confined between planar fracture surfaces. The proppant-absent and proppant-filled fractures show different frictional strength. For fractures with proppants, the frictional response is mainly controlled by the normal stress and proppant thickness. The depth of shearing-concurrent striations on fracture surfaces suggests that the magnitude of proppant embedment is controlled by the applied normal stress. Under high normal stress, the reduced friction implies that shear slip is more likely to occur on propped fractures in deeper reservoirs. The increase in the number of proppant layers, from monolayer to triple layers, significantly increases the friction of the propped fracture due to the interlocking of the particles and jamming. Permeability of the propped fracture is mainly controlled by the magnitude of the normal stress, the proppant thickness, and the proppant grain size. Permeability of the propped fracture decreases during shearing due to proppant particle crushing and related clogging. Proppants are prone to crushing if the shear loading evolves concurrently with the normal loading.

1. Introduction

Hydraulic fracturing has been a major well stimulation technique since the 1940s [1]. The process involves the injection of fracturing fluid into a targeted reservoir layer through a wellbore. The high pressure fluid drives the propagation of hydraulic fractures with proppant particles carried with the penetrating fluid along the fracture and into the reservoir formation [2]. Following the injection phase, the fracturing fluid flows back into the wellbore and the created hydraulic fractures will compact due to the release of fluid pressure. The proppant particles, typically made of sand, treated sand, or man-made ceramic materials, however hold the hydraulic fractures open, allowing the propped fractures to act as highly permeable flow paths for the extraction of hydrocarbons.

Over the past decade, massive hydraulic fracturing has been broadly applied for the extraction of tight gas, shale gas, and other unconventional resources [3]. This broad application of horizontal drilling and stimulation, with multiple horizontal wells per pad, multiple fracturing stages per well, and multiple clusters per stage, has resulted in the successful recovery of hydrocarbons from ultralow permeability reservoirs. These include the recovery of unconventional resource from shales and tight sandstones [4]. Hydraulic fracturing in unconventional reservoirs can be significantly different from the hydraulic fracturing of conventional reservoirs. This is due to the presence of preexisting natural fractures or weak planes in the formation that may significantly affect the propagation of the fluid driven fractures. Based on various theories of conventional hydraulic fracturing, a bi-wing type fracture is typically assumed to result from the stimulation of a vertical well [5]. However, based upon field production data, as well as microseismic observations, it is believed that so-called “fracture complexity” may result from interactions between the created hydraulic fracture and preexisting natural fractures. This is especially prevalent in unconventional reservoirs where the contrast in permeabilities between matrix and natural fracture is very high [68]. When a propagating hydraulic fracture intersects a natural fracture, multiple scenarios for the form of the crossing are possible. These include direct crossing, hydraulic fracturing arrested by the natural fracture, crossing with an offset, or even more complicated scenarios when the third dimension is considered [916]. Therefore, fracture branching and the development of complex fracture networks are generally created during multistage hydraulic fracturing (Figure 1(a)) with implications for the state of stress applied on such oblique fractures (Figure 1(b)). Note that the near-wellbore tortuosity [17, 18] is not included in the schematic of Figure 1 but could greatly affect the fracture propagation and cause the increase of injection pressure.

Significant effort has been applied to understand mechanisms involved in creating hydraulic fracture networks in unconventional reservoirs [19]. Seismicity-permeability coupling relationship has been attempted in enhanced geothermal reservoir stimulation [20] and in the caprock for geological sequestration of CO2 [21, 22]. However, less attention has been given to understanding the role of shear deformation (e.g., induced seismicity by hydroshearing) and permeability evolution of the complex fracture system that may evolve during long-term depletion. Injection-induced seismicity has been associated with both the operation of waste water reinjection [23] and hydraulic fracturing [24]. Nevertheless, whether induced seismicity can occur during long-term depletion of unconventional fractured reservoirs is of significant scientific interest but remains poorly understood. In the usual conceptual model of conventional hydraulic fracturing, the created fracture plane is perpendicular to the minimum stress direction; thus, there is no shear stress on the fracture as it propagates and before the depletion stage of reservoir production. After stimulation, additional shear stresses induced during the depletion process should also be negligible for simple fracture geometry. Once depletion begins, considering the Mohr-Coulomb failure criterion, the gradual decrease in pore pressure will increase the effective clamping stress on the fracture plane (Biot coefficient < 1) and further stabilize the fracture planes [25]. Therefore, it can be concluded that the expectation of induced seismicity during conventional hydraulic fracturing is trivial.

However, for unconventional reservoirs, induced seismicity may result due to the more complex and oblique fracturing geometry. First, the hydraulic fracture planes are not necessary aligned with the maximum horizontal stress, due to the presence of preexisting natural fractures and related fracture complexity. Second, nonuniform depletion due to heterogeneous permeability fields can cause stress reorientation and additional shear stress on the fracture planes due to poroelastic effects [26, 27]. Third, some field operations, such as the failure of diversion during refracturing and undesired well connection when fracturing a new well (i.e., cross-well communication or frac-hit), can lead to fluid leakage into the preexisting hydraulic fractures which have been under depletion and result in additional slip on preexisting fractures [28]. Fourth, reinjecting fluid into wells under depletion has been implemented in the field to boost the production under some circumstances, for example, if there is a sharp decline of production, or delineated depleted zone around the producing well [29]. Finally, numerous experimental studies show that the interaction between propagating hydraulic fractures and preexisting natural fractures is significantly influenced by differential stress and fracture orientation and frictional strength [9, 12]. As the frictional strength of propped fractures in the first fracturing phase may be altered, it may significantly influence the behavior of fracture propagation in the later refracturing phase. These circumstances make the potential for induced seismicity finite for unconventional resource recovery. Induced seismicity is governed by frictional behavior of the fracture surface contact [21]. Moreover, for both conventional and unconventional hydraulic fracturing, the main purpose is to create propped fractures and to increase the permeability of the formation. Thus, the evolution of permeability of the propped fractures is the key parameter that ultimately affects well production. In the course of production, the permeability of propped hydraulic fractures is expected to decrease due to the impacts of proppant embedment and crushing, concomitant with the gradual increase in the effective clamping stress (Zhang et al., 2015). Previous experimental studies only explored the shear behavior of joints filled with sandy granular materials [30, 31]. However, detailed coupling mechanisms involved in shearing and permeability evolution for propped fractures remain unclear. In this study, we explore the primary and secondary frictional and fluid transport response of propped fractures as the fundamental controlling parameters involved in the reactivation of preexisting fractures (Figure 1(b)). This experimental and analytical work reveals the coevolution of permeability with friction of a propped fracture under shear.

2. Experimental Methods

We explore the evolution of friction and permeability on propped fractures. We first present the procedure of sample preparation for the experiments. Then, we introduce the experimental setup and testing procedure implemented to explore the evolution of permeability and friction of a propped fracture under shearing. Finally, theories and methods to calculate the evolution of friction and permeability are discussed.

2.1. Sample Preparation

The experiments are completed on Green River shale (GRS) as an appropriate analogue of a shale reservoir. To provide a contrast in rock texture, Westerly granite (WG) was also used as a reference, because it has been extensively studied and is well suited for comparison due to its homogeneous and isotropic structure. The mineralogical compositions and mechanical properties of GRS and WG are listed in Table 1. The rock samples were first cored to a length of 2 inches and diameter of 1 inch and then carefully saw-cut into two halves, representing a parallel plate model (Figure 2). The planar surfaces were uniformly polished with abrasive powder (#60 Grit carbide) to provide consistent surface roughness for all fracture analogues. To prevent the dislocation of proppant particles during the process of sample reassembly, a very thin layer of washable glue was placed on the fracture surfaces to temporarily fix the proppant particles. The proppant particles were uniformly and tightly placed on the surface of the fractures, forming a monolayer. To evaluate the effect of proppant thickness, samples with double and triple layers of proppants were also assembled and tested. The reassembled split samples, with proppants embedded, were packed within a latex membrane with an initial offset of 8 mm to accommodate the shear offset applied during shearing. To reduce the friction between the outer wall of the sample and the membrane, we used Teflon tape to cover the outer wall of the sample, through which the extra friction by the system can be significantly reduced.

In this experiment, we used commercial ceramic proppants. Three typical proppant sizes, that is, 40/80 mesh (180~425 μm), 30/50 mesh (300~600 μm), and 20/40 mesh (420~840 μm), were used in the experiments of this study. The proppant size is often referred to as the sieve cut and is typically in the range between 8 and 140 mesh (105 μm~2.38 mm). The exact distributions of proppant grain size in this experiment are shown in Figure 3. The 40/80 mesh was used as a baseline for the purpose of comparison, while the 30/50 mesh and 20/40 mesh were utilized to study the effect of proppant size.

2.2. Experimental Setup and Testing Procedure

The experiments were performed in a triaxial testing apparatus that is able to independently apply confining pressure, pore pressure, and shear displacement at prescribed (constant) velocity. The evolution of fracture permeability during the experiments can also be concurrently monitored (Figure 4(a)). The packed sample was then assembled in the cylindrical vessel. The confining stress (normal stress in this configuration) was gradually applied until the desired magnitude was reached. With the desired normal stress applied, deionized water was circulated through the fracture with a constant upstream pressure for 5 mins to dissolve and remove the glue that was used to fix the proppant during sample assembly. Once a steady flow rate was attained, monotonic shearing at constant velocity was applied. The shear velocity was controlled at 3 μm/s and the shearing was stopped after a displacement of 6 mm was reached. All experiments were performed at room temperature. The shear displacement was recorded by LVDT installed at the end of the displacing piston. Before and after experiment, samples were characterized by white light optical profilometry to observe the possible interaction (e.g., embedment or striations) between proppant particles and the shale surface due to shearing. White light profilometry was performed using a Zygo NewView 7300 profilometer with a 10x objective lens with data processed with Mx™ software (Figure 4(b)). Furthermore, to examine the shear-induced damage or crushing of proppants, the proppant particles were scanned and sized both before and after experiment by a laser particle size analyzer.

2.3. Friction and Permeability Calculation

We calculated the coefficient of friction as a function of shear displacement using the ratio of measured shear stress to applied normal stress as and ignoring cohesion. A parallel plate model for the cubic law is typically employed to describe fluid flow within a fracture [32]; however, as the number of proppant layers inside the fracture increases, the pattern of fluid flow experiences a transition from parallel plate flow to porous medium flow (Slichter, 1899; Kozeny, 1927; Carman, 1937; Li, 2017) (Figure 5). In the direct-shear experiment, the proppant particles may be dislocated, deform, and even break in shearing. Nevertheless, we define an equivalent hydraulic permeability based on Darcy’s law as follows: where (Pa·s) is the viscosity of fluid; (m) is the contact length of the fracture surface; (m) is the fracture width; (m3/s) is the measured flow rate; and (Pa) is the differential pressure between the upstream and downstream extent of the fracture.

The permeability of a porous medium can be estimated based on the Kozeny-Carman equation [33] as follows: The permeability is then related to the porosity of the packed bed and the particle diameter .

3. Experimental Results

In this section, the results are presented to highlight key experimental observations related to the effects of normal stress, proppant thickness, proppant size, and rock texture on the friction-permeability relationships of propped fractures.

3.1. Effect of Normal Stress

Figures 6(a)6(c) show the evolution of permeability, normalized permeability, and friction for the propped fracture during shearing under different normal stress, that is, 1 MPa, 3 MPa, and 5 MPa, respectively. A monolayer of proppant is present in all three cases with a corresponding proppant size in the range 0.18 mm to 0.425 mm (40/80 mesh). The total shear displacement is 5 mm. As expected, the initial permeability before shearing decreases with an increase in normal stress due to the combined influence of reduced porosity and fracture closure (Figure 6(a)). Permeability gradually declines during the shearing for all three cases. The normalized permeability indicates that the reduction in permeability is most profound for the case with the highest normal stress. At the end of loading, the permeability for the three cases decreases to be ~70%, 40%, and 20% of the initial values. A plausible mechanism for explaining this phenomenon is that the proppants crush the most during shearing for the case with the highest normal stress, which causes the largest relative decrease in apparent fracture aperture. Particle crushing is apparent from the grain size distributions collected before and after shearing shown in Figure 6(d), where the particle damage increases with normal stress, inferring its influence on permeability evolution. Note that, for the case of 1 MPa normal stress, there is nearly no change in the grain size distribution before and after testing (i.e., no crushing of particles during shearing). Thus the grain size distribution of the initial proppant packing is not shown in the plot.

The frictional resistance of the propped fracture in Figure 6(c) also displays dependence on the normal stress. The coefficient of friction decreases as the normal stress increases, which is consistent with friction-normal stress relationships from previous experimental studies on simulated gouge [37, 38]. The reduction in frictional strength may be attributed to two possible causes. First, at higher normal stresses, the normalized membrane restraint between the sample surface and the membrane is reduced. As the normal stress increases, the coefficient of friction converges to the actual value representing the contact behavior between proppant particles and the fracture surface. Second, higher normal stress will compact the proppant particles and result in the crushing and embedment of proppants, which changes the contact response at the interface between proppant particles and fracture surfaces. Figure 7 explicitly compares the topographies of the fracture surfaces both before and after slip and for different normal stresses. Each experiment is performed with a virgin fracture. A notable feature is that there are an increasing number of postshearing striations on the fracture surface as the normal stress increases. These striations result from the embedment of proppant particles into the fracture surface with grooving owing to the shear loading. To further characterize the striations, the fracture surfaces were scanned with white light profilometry. Figure 8 compares the profiles of fracture surfaces both before and after shearing for the case of 5 MPa normal stress. The dark channel represents a striation with a depth of about 100 μm. As the normal stress increases, the frictional behavior then gradually transits from being governed by the sliding of particles along the fracture surface to along the surface of striations. Since the surface of the striation is less rough than the initial fracture surface, the friction may decrease for the highest normal stress as shown in Figure 6(c).

3.2. Effect of Proppant Thickness

Figure 9 presents the evolution of permeability, normalized permeability, and friction during shearing for different proppant thickness as influenced by the number of proppant layers. Four cases are represented, showing the behavior for bare surfaces (i.e., no proppants embedded) and with monolayer and double and triple layers of proppant. The normal stress for all four cases is 3 MPa and the proppant size is 40/80 mesh. For idealized close-packing of monosized particles filled between two parallel plates, the monolayer structure has the largest porosity (0.3954) while the porosity approaches the minimum magnitude (0.2595) as the number of layers increases and the aggregate conforms to a face-centered cubic (FCC) structure (Figure 10). From this rationale, and with monodisperse particle sizes, the permeability for a monolayer of proppant is expected to be the largest among all the cases. However, Figure 9(a) shows that the initial permeability of a monolayer of proppant is actually smaller than that for double and triple layers. This may be attributed to the permeability decrease due to proppant embedment where the case with a monolayer of proppant is affected the most. In the case of a monolayer, the fewer displacement degrees of freedom offered where the proppant is sandwiched between rigid faces, rather than compacting to the interior, result in a greater embedment.

Permeability gradually decreases during shearing for all four cases; however the mechanisms of permeability reduction between bare surfaces and proppant embedded surfaces are different. For the bare surfaces, the permeability decreases due to the generated wear products [21], while, for proppant embedded fractures, the permeability is reduced due to proppant crushing, embedment, and clogging during shearing. Furthermore, the normalized permeabilities indicate that the relative decrease in permeability increases as the proppant thickness decreases. The case without proppants has the largest permeability drop (i.e., more than 80%), while the case with triple-layer proppants has less than a 20% decrease in permeability. The grain size distributions of proppants, both before and after the experiments (Figure 9(d)), indicate that the case with double-layer proppants has the most profound crushing of the proppant. However, the reason why the case with double-layer proppants exhibits the most particle crushing (relative to both single and triple layers) is not yet clear.

The friction of the propped fracture for these three cases with proppants shows a clear trend that more proppant layers result in an increase in frictional strength. This is possibly due to the increase in interlocking forces and jamming between particles when multiple proppant layers are present during shearing. The case without proppant exhibits a frictional strength similar to that of double-layer proppants but larger than the friction value for a monolayer of proppants.

3.3. Effect of Proppant Size

Figure 11 shows the evolution of permeability, normalized permeability, and friction for the propped fracture during shearing for the three different proppant sizes shown in Figure 3. All three cases are for a monolayer of proppant at a normal stress of 3 MPa. As expected, the initial permeability decreases as the proppant size decreases owing to the smaller initial apparent aperture and pore-throat size.

The normalized permeability for the 20/40 mesh proppant remains near constant during shearing, while that for the 30/50 mesh decreases by ~20% and that for the 40/80 mesh proppant decreases by ~50%. To explain why the smaller proppant is subject to a larger permeability drop during shearing, Figure 11(d) plots the grain size distributions of the proppants both before and after the experiments for each proppant size. The smallest proppant size (40/80 mesh) shows significant particle crushing; however, particle crushing is not clearly identified with either the 30/50 mesh or 20/40 mesh proppants. While particle crushing could explain the severe permeability drop during shearing for the case with 40/80 mesh, it may not resolve the paradox of permeability drop for the 30/50 mesh proppant. Another possible reason for the permeability drop, besides that of particle crushing, is the potential reorganization of particles during shearing and the possibility of particle clogging. Although not verified in the experiments, due to limitations of measurement, clogging is highly possible since the initial particle packing is relatively loose.

The friction of the propped fracture for these three cases is nearly identical which suggests that, given the identical normal stress, though the proppant size varies, the contact state between the proppant surface and the fracture surface is equivalent in each case. The results indicate that friction is mainly governed by the stress state as well as the degree of embedment between proppant particles and fracture surface.

3.4. Effect of Rock Texture

To investigate the effect of texture of the fracture surfaces and to eliminate the role of striation-formation on the response, the experiments were repeated with Westerly granite at 3 MPa normal stress and with a monolayer of proppant (40/80 mesh). Figures 12(a)12(c) show the evolution of permeability, normalized permeability, and friction during shearing for the two different rock textures. The case with granite has a slightly larger initial permeability than that with shale, although both cases show a decrease in permeability during shearing. Granite is stiffer and of higher strength than shale; therefore the proppant embedment is smaller than shale under the same normal stress. This may explain the higher initial permeability for the granite before shearing. The permeability for both granite and shale sandwiching-fractures converges to similar magnitudes at the end of shearing. The normalized permeability indicates that the decrease in permeability for granite is slightly larger than that for shale. This is because there is more particle crushing in the case with granite than with shale, as apparent in the postexperiment reduction in the particle size distribution for granite (Figure 12(d)).

Friction of the propped fracture for granite is also slightly larger than that for shale, which could be attributed to the mineral-particle contact state (i.e., proppant-quartz for WG and proppant-calcite for GRS) and the larger amount of generated small particles as a result of particle crushing.

4. Discussions

The experimental results show that the main factors controlling the frictional behavior of a propped fracture are the normal stress and proppant thickness. High normal stress results in the crushing of proppant particles, reducing the mean size of the proppant PSD (Figure 6(d)). However, this change in size has only limited impact on the frictional response of the encasing fractures (Figure 11(c)). Under high normal stress, the normalized membrane restraint between the sample surface and the membrane is reduced, implying that at higher normal stress this contributes proportionally less to the frictional resistance and yielding a strength closer to the real strength of the assemblage. This strength is smaller than at lower stresses, due to the absence of the spurious shear restraint provided by the membrane. High normal stress also generates striations on the fracture surfaces during shearing, which allows a smoother contact between proppant particles and the fracture surface. As a result, the overall friction reduces as the normal stress increases. The proppant thickness also plays an important role in defining the friction of a propped fracture. Interlocking of particles and jamming are not expected between proppant particles for the monolayer proppant configuration. However, given the saw-cut fracture surface geometry, as the number of proppant layers increases from one to three, the interlocking forces between particles largely increase the friction of the propped fracture. The effects of proppant size and rock texture on the friction of a propped fracture are secondary. In terms of the response of proppants in fractures, the friction of the fracture without proppant is significantly different from that with proppant (Figure 9(c)). For a single layer propped fracture, the reduced friction implies that, for refracturing in a propped fractured reservoir, the hydraulic fracture may be arrested by a propped fracture. It should be noted that the impact of fracture asperities on the friction of a propped fracture is ignored in this study. This particular fracture analogue comprises two flat surfaces with a uniform but minimal roughness (controlled by the PSD of the grinding powder), while, in reality, the resulting fracture surfaces may have considerably higher amplitude roughness.

The permeability of a propped fracture is mainly governed by the normal stress, the proppant thickness, and the proppant size. The normal stress controls the amount of proppant embedment and thus the dilation of the fracture aperture during shearing. High normal stress not only causes the compaction of the fracture and of the proppant bed, but also leads to the crushing of proppant particles that accelerates the fracture closure. Compared to the multilayered specimen, the monolayer case exhibits the smallest initial permeability due to proppant embedment. Although a larger proppant size favors a higher initial permeability of the propped fracture, the in-fracture transport of large diameter proppants during completion of slickwater fracturing is difficult; the issue of proppant size selection is beyond the scope of this study.

Except for the case with a monolayer of 20/40 mesh proppant, the permeability of the propped fracture decreases during the shearing process for all other cases. The comminution-related decline in permeability during shearing dominates over the effect of shear-induced dilation for a fracture with flat surface geometry. However, it is not excluded that permeability enhancement by shear slip of a fracture with high surface roughness would dominate over the proppant-crushing induced permeability decline. Another unexpected conclusion drawn from this study is that significant proppant crushing occurs during shearing even at a normal stress of 5 MPa. Since proppants are typically designed to withstand normal stresses as high as 50 MPa, they become vulnerable if the shear loading evolves concurrently with the normal loading.

It is also worth mentioning that there are some limitations in the experiments. First, the direct observations of proppant crushing and embedment are not feasible while keeping specimens in an in situ stress state. These results are indirectly reflected by the measurements of particle size and surface characteristics after the experiments. Second, the distribution of proppant after shearing cannot be measured; thus crucial information on proppant clogging is missing. The real time proppant clogging status may only be tested via the imaging techniques such as xCT scanning. Last, the normal stress applied by the apparatus is modulated by the strength of the aluminum ring (shown in Figure 4(a)) used to protect the void left by the offset distance between the two fractures. To prevent radial deformation of aluminum ring, the highest normal stress applied in the experiments is required to be less than 6 MPa. Although the normal stress in the experiments can be interpreted as the effective stress applied on the fracture wall, it is still lower than a typical magnitude in the field.

5. Conclusions

In this study, we explore the evolution of friction and permeability of a propped fracture using shearing-concurrent measurements of permeability during constant velocity shearing experiments. We separately examine the effects of normal stress (1 MPa, 3 MPa, and 5 MPa), proppant thickness (monolayer, double-layer, and triple-layer), proppant size (40/80 mesh, 30/50 mesh, and 20/40 mesh), and rock texture (Green River shale and Westerly granite) on the frictional and transport response of proppant packs confined between planar fracture surfaces. The results indicate that proppant-absent and proppant-filled fractures show different frictional strength, implying that proppants could change the friction of natural fractures and influence the potential for shear failure. For fractures with proppants, we observed the following: (1) the frictional response is mainly controlled by the normal stress and proppant thickness. High normal stress results in the crushing of proppant particles although this change in size has almost no impact on the frictional response of the proppant-fracture system. The observed postshearing striations on fracture surfaces suggest that the magnitude of proppant embedment is controlled by the applied normal stress. Moreover, under high normal stress, the reduced friction implies that shear slip is more likely to occur for the propped fractures in deeper reservoirs. With this simple specific fracture configuration (i.e., saw-cut surface), the increase in the number of proppant layers, from monolayer to triple layers, significantly increases the friction of the propped fracture due to the interlocking of the particles and jamming, suggesting that high proppant density during emplacement would help stabilize the fractures during injection. (2) Permeability of the propped fracture is mainly controlled by the magnitude of the normal stress, the proppant thickness, and the proppant size. Permeability of the propped fracture decreases during shearing, which is plausibly due to proppant particle crushing and related clogging. Compared to the multilayered specimen, the monolayer case which has fewer displacement degrees of freedom exhibits the smallest initial permeability due to proppant embedment. Proppants become prone to crushing if the shear loading evolves concurrently with the normal loading. The above combined conclusions suggest the use of high-concentration proppants in the field, which not only provides high hydraulic conductivity for hydrocarbon production, but also may help to mitigate the risk of induced seismicity.

Disclosure

Yi Fang is now working at Institute for Geophysics, The University of Texas at Austin, Austin, TX 78712, USA.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

The authors acknowledge the support provided by the Young 1000 Talent Program of China, Tongji Civil Engineering Peak Discipline Plan (CEPDP), National Natural Science Foundation of China under Grants 41772286 and 51674267, and US Department of Energy (DOE) under Grant DE-FE0023354. The useful discussions with Professor Chris Marone are also greatly appreciated.