Discrete-Element Modelling of Bedding- and Stress-Controlled Hydraulic Fracturing in Deep Qiongzhusi Shale

The Cambrian Qiongzhusi Formation shale of the Sichuan Basin, buried below 4000 m, is a key deep-shale-gas target whose hydraulic fracturing is jointly controlled by bedding weak planes and high in situ stress. This study separates the two controls quantitatively and establishes a reproducible criterion for selecting favourable bedding attitudes and scaling treatment energy. A fluid–solid coupled discrete-element model was built on PFC2D, combining parallel-bond matrix contacts, smooth-joint bedding planes installed on a discrete fracture network, servo-controlled stress boundaries and a cubic-law fluid network, with micromechanical parameters calibrated against a triaxial test on a core from a 4601.76 m depth. Seven bedding inclinations (0–90°) were simulated under a unified stress state, and nine stress combinations were arranged in a 3 × 3 orthogonal design of basal horizontal stress (88.1–125 MPa) and vertical–horizontal stress difference (0.1–6 MPa) at a fixed 60° bedding, then analysed by quadratic response surface and two-way ANOVA. The peak fracturing pressure follows a U-shaped relation with bedding inclination spanning 110–160 MPa. Low-angle bedding (0–30°) favours lateral extension (ratios: 0.63–0.68) but limits vertical connectivity, whereas high-angle bedding (75–90°) promotes vertical breakthrough with lateral reach falling to 0.08–0.09. Ranked by a composite stimulation score built from four normalised sub-scores—initiation efficiency, lateral reach, vertical reach and an SRV proxy—60° is the optimum (S = 0.76), ahead of 15° (0.72) and 45° (0.67), combining the lowest fracturing pressure of 110 MPa with progressive branching rather than sudden breakthrough. Raising basal σh from 88.1 to 125 MPa lengthens propagation approximately fourfold and shifts the failure mode from tensile splitting to shear–tensile ductile rupture, with the tensile-fracture proportion falling from 1.0 to 0.25–0.35. The two-way ANOVA attributes 82–99% of the variance of every response to σh (p ≤ 0.031), whereas the stress difference contributes at most 8.2% and reaches significance for none (p > 0.17) The bedding optimum is placed on a stated objective function rather than qualitative judgement, and the two stress controls are shown to combine additively rather than multiplicatively within the Qiongzhusi stress envelope. Practically, treatment energy should be scaled with basal σh rather than with the stress difference, and target intervals selected at 45–60° bedding.