Deep marine shale gas reservoirs are extremely rich in the Sichuan basin in China. However, due to the in situ conditions with high temperature and high pressure (HTHP), in particular reservoir pressure being usually much higher than the test pressure, it is difficult to accurately clarify the adsorption behavior, as seepage theory plays an important role in shale gas reserves evaluation. Therefore, three kinds of sorbent, including illite, quartz and kerogen, and two simulation methods, containing the grand canonical ensemble Monte Carlo method and molecular dynamics method, are synthetically used to determine the methane adsorption behavior under HTHP. The results show that both absolute adsorption and excess adsorption decrease with the increase of temperature. When the pressure increases, the absolute adsorption increases quickly and then slowly, and the excess adsorption first increases and then decreases. The superposition of wall potential energy is strongest in a circular hole, second in a square hole, and weakest in a narrow slit. The effect of pore size increases with the decrease of the pore diameter. Under HTHP, multi-layer adsorption can occur in shale, but the timing and number of layers are related to the sorbent type.

China has abundant shale gas resources with good exploration value and development potential, making it a recent hotspot for exploration and development. It is widely agreed that large-scale hydraulic fracturing is essential for reservoir enhancement in shale formations. However, the evolution of fractures during hydraulic fracturing is highly complex, necessitating research on the influence of various factors on the vertical propagation of hydraulic fractures. Based on geological and engineering parameters from the Luzhou block in southern Sichuan, this study employed the finite element method (FEM) and the cohesive element method to establish a coupled fluid-solid model for the vertical propagation of hydraulic fractures. Numerical simulations were conducted to investigate the interaction between hydraulic fractures and natural weak planes, clarifying the mechanisms involved. This study elucidates how different rock and natural weak plane properties affect the vertical propagation of hydraulic fractures and draws diagrams illustrating these interactions. The research indicated three fracture distribution patterns after the intersection of hydraulic fractures with natural weak planes: passive fractures, ‘I’-shaped fractures, and crossing fractures. The main fractures in these patterns exhibit initial damage and damage evolution characterized by tensile failure. Specifically, in passive fractures, the initial damage and damage evolution of natural weak planes manifest as shear failure. In ‘I’-shaped fractures, the initial damage in natural weak planes is characterized by shear failure, with damage evolution showing tensile failure. Crossing fractures show minimal damage in the weak planes. Under conditions of high natural weak plane cohesive strength, high Young’s modulus, low interlayer rock cohesive strength, high vertical stress difference, low interlayer stress difference, and high intersection angles, crossing fractures tend to form. Conversely, conditions of low natural weak plane cohesive strength, low Young’s modulus, high interlayer rock cohesive strength, low vertical stress difference, high interlayer stress difference, and low intersection angles favor the formation of ‘I’-shaped fractures. Passive fractures form under conditions of low natural weak plane cohesive strength and high vertical stress difference. This study found that Poisson’s ratio has a minimal effect on the vertical expansion of hydraulic fractures under the studied conditions, with natural weak plane strength being the primary control factor for fracture patterns. These findings enhance the theoretical foundation for the vertical propagation of hydraulic fractures in deep shale formations, facilitating the development and implementation of strategies for enhancing production in shale reservoirs with natural weak planes and better optimizing production in different types of shale reservoirs.

Multi-stage hydraulic fracturing is a widely used technology in the development of shale oil and gas reservoirs that creates artificial fractures and forms fracture networks that enhance fluid flow within reservoirs. A well-designed fracture network can significantly enhance the production performance of oil and gas wells, thereby improving the recovery of shale oil and gas reservoirs. To achieve this, understanding the gas production performance of individual artificial fractures is crucial, as it provides valuable insights for refining subsequent fracturing designs, ultimately leading to an optimized fracture network design. At present, numerical simulations are commonly used to study the production performance of individual artificial fractures by modeling the production process of the entire shale oil and gas reservoir. However, due to the heterogeneity of reservoirs and the presence of numerous natural fractures, traditional numerical simulations require high-resolution grids to model the production process, making them computationally expensive and time-consuming. To address this issue, in this work, based on the concept of diffusive time of flight (DTOF), the authors propose a fast simulation method to efficiently simulate the production behavior of individual artificial fractures throughout the shale oil and gas reservoir production process. The DTOF can be obtained by solving the Eikonal equation using the fast marching method (FMM), which is then used to calculate the drainage volume of individual artificial fractures. Subsequently, a geometric approximation of the drainage volume is used to efficiently compute the production rates of individual artificial fractures. Unlike traditional numerical simulations, this method uses a single non-iterative calculation to determine the drainage volume of individual artificial fractures, followed by a geometric approximation to compute the production rates. This eliminates the need for high-resolution grids, significantly reducing computational cost and time, which allows the proposed method to provide faster simulations compared to traditional numerical methods while maintaining sufficient accuracy.