水力压裂是开发干热岩地热能的关键技术，通过向高温地层注入高压流体形成人工热储，以实现高效取热^[1-3]。然而，干热岩造储过程中，压裂裂缝受地应力场控制，形态单一且不能定向沟通目标井，严重限制了注采井间的连通性^[4]。因此亟需探索能够克服地应力控制，实现裂缝定向沟通注采井的压裂新方法。针对上述难题，提出径向井引导压裂裂缝定向扩展的新思路，利用径向井应力重构作用引导裂缝扩展方向，最终在注采井间形成定向连通的裂缝网络^[5-9]。针对径向井参数和工程地质参数对裂缝起裂与扩展的影响，国内外学者开展了深入研究^[10-15]。其中，孔祥伟等^[10]发现定面射孔可调控裂缝初始破裂位置。Guo等^[11]通过压裂试验证实，相比射孔压裂，径向井压裂在降低破裂压力和引导裂缝扩展等方面展现出更大的优势。郭天魁等^[12]对各影响因素进行分析，其中方位角、地应力和排量是影响裂缝扩展的3个主要因素。曲占庆等^[14]基于扩展有限元方法，指出不同相位的多分支井能有效提升压裂缝网复杂度。尽管已有大量研究，但考虑高温条件下的径向井压裂研究较少，热应力和径向井耦合作用下裂缝扩展模式仍不明晰。为此，笔者建立热流固耦合径向井压裂干热岩裂缝扩展数值模型。通过与解析解及试验对比，验证模型准确性。在此基础上，开展不同温差下干热岩径向井压裂模拟，与垂直井压裂对比，阐明热应力对裂缝扩展方向的影响机制，并研究不同地应力、压裂液排量和黏度条件下的裂缝扩展模式。

1 干热岩径向井压裂模型

采用三维离散格子法建立干热岩水力压裂裂缝扩展模型，该方法通过合成岩体技术和离散格点理论来模拟压裂过程中裂缝的萌生和扩展^[16]。其中模型网格由一系列随机分布的质量节点和弹簧组成，如图1所示。岩石颗粒之间的接触通过弹簧实现，当弹簧力超过其强度时岩石发生破坏，对应微裂纹的形成。流体节点位于弹簧的中心，相邻的流体节点通过管道相连。随着新裂纹的形成，流体管网不断更新。

1.1 热流固耦合模型控制方程

1.1.1 力学模型

离散格子法中每个节点具备3个平动自由度及3个转动自由度，平动自由度采用线性动量平衡方程和位移速度关系来计算^[17]，即

式中，${\dot{u}}_i^{\left(t\right)}$为分量i（i=1，2，3）在时刻t的速度，m/s；$u_i^{\left(t\right)}$为分量i（i=1，2，3）在时刻t的位移，m；∑F_i为作用在质量为m的节点上所有力分量i的总和，N；Δt为时间步长，s；m为节点的质量，kg。

分量i（i=1，2，3）在时刻t的角度自由度w_i由以下中心差分方程计算^[17]：

式中，∑M（t）_i为作用在转动惯量节点I上所有力矩分量i（i=1，2，3）的总和，N·m；I为转动惯量，kg·m²。

在访问完所有节点后，将对所有弹簧进行扫描。如果弹簧没有断裂，则在时间t处进行以下计算：

式中，下角rel表示相对；下角A和B表示通过弹簧连接的两个节点。

式中，下角N和S分别表示法向和剪切；n_i为单位法向矢量。

然后，更新弹簧的法向力和剪切力^[17]为

式中，F_N和F_i，S分别为弹簧的法向力和切向力，N；k_N和k_S分别为弹簧的法向刚度和剪切刚度，N/m；N为法向速度，m/s；${\dot{u}}_{i，\mathrm{S}}$为切向速度，m/s。

最后，将新的弹簧力添加到相关节点的合力中：

微观弹簧的刚度和抗拉强度与宏观岩石的体积模量和抗拉强度的对应关系^[18]为

式中，F_Nmax和F_Smax分别为弹簧的抗拉强度和抗剪强度，N；α_t和α_s分别为抗拉强度校正系数和抗剪强度校正系数，无量纲；T和C分别为岩石的宏观抗拉强度和宏观抗剪强度，Pa；R为单元尺寸，m；μ′为摩擦系数，无量纲。

当弹簧的法向力超过弹簧的抗拉强度或弹簧的切向力超过弹簧的抗剪强度时，弹簧标记为断裂，此时弹簧力归零，F_N=0，F_i，S=0。

1.1.2 流动模型

裂缝中流体流动是通过流体节点和管道网络来模拟。如图1所示，流体节点位于断裂弹簧的中心，通过管道连接相邻流体节点，多个连通的管道形成流动管网。

沿管道从流体节点A到节点B的流速是基于平面泊素尔流^[17]计算的：

式中，q为流体流量，m³/s；β为修正系数，无量纲；K_r为相对渗透率，无量纲；a为裂缝开度，m；μ为流体黏度，Pa·s；p_A和p_B为节点A和节点B处的流体压力，Pa；ρ_w为流体密度，kg/m³；g为重力加速度，m/s²；z_A和z_B分别为节点A和节点B处的位置水头，m。

K_r为饱和度S的函数，计算公式为

当管道充满流体时（S=1），相对渗透率等于1。

在每个流动时间步长Δt期间内，压力增量Δp计算公式^[17]为

其中

$Q=\sum _i{q}_i.$

式中，Δp为流体压力增量，Pa；q_i为与节点i相连的管道流量，m³/s；Q为连接流体节点所有通道的总流量，m³/s；V为流体节点体积，m³; ${\stackrel{-}{K}}_{\mathrm{F}}$为流体体积模量，Pa；Δt_f为单位流体步长，s。

1.1.3 传热模型

向储层岩石中注入冷流体会引起岩石温度的变化、热致变形以及岩石破裂。温度被储存在晶格弹簧中，基质内的热流直接根据局部弹簧的温度差异计算得出（图1），具体公式^[19]为

式中，$\stackrel{-}{T}$为节点的虚拟温度，℃；η_i为热导，J/（s·℃）；c_V为比定容热容，J/（kg·℃）；Δt_t为热传导时间步长，s。

修正节点间弹簧作用力的法向分量，将热应变计入到晶格模型中：

式中，k_n为弹簧的法向刚度，N/m；α为膨胀系数，与岩石的宏观线性热膨胀系数相同，℃^-1；L为弹簧长度，m；ΔT为温度差，℃。

单位时间和单位面积从岩石到流体的传热q_cv，计算公式为

其中

$h=\frac{k_{\mathrm{r}}}{\beta \sqrt{kt}},k=\frac{k_{\mathrm{r}}}{{\rho }_{\mathrm{r}}{c}_{\mathrm{r}}}.$

式中，h为传热系数，W/（m²·℃）；k_r为岩石导热系数，W/（m·℃）；β为传热参数（1阶无量纲常数，用作校准因子）；t为注入时间，s；k为岩石热扩散率，m²/s；ρ_r为岩石密度，kg/m³；c_r为岩石比热容，J/（kg·℃）。

裂缝中的热平流使用流体单元来存储流体温度。在与流体单元相连接的弹簧处，流体单元中的温度是流体温度T_f，晶格弹簧中的温度是岩石温度T_r，这在裂缝两侧是相同的。流体单元处的新温度由前一个时间步的值计算，使用以下增量关系^[19]：

式中，Q_f为管道体积流量，m³/s；V_f为与流体单元相关联的控制体积m³；T_f，up为位于上游的流体单元处的温度，℃；a_f为流体单元的开度，m。式（22）右侧的第一项是校正项，第二项是流动流体中热平流的温度贡献，最后一项是岩石与流体之间的热交换。

1.1.4 热-流-固耦合过程

本文中使用离散格子法实现了一种紧密耦合的热-流-固耦合模拟。其具体表现为：①流体的流动压力会导致岩石的变形及破坏；②岩石的变形会导致裂缝中流体压力的波动和裂缝开度的变化，进而影响裂缝渗透率；③岩石中的热应变也会引起裂缝开度的变化，进而干扰局部流体压力；④流体在裂缝中流动时受渗透率的影响^[17-19]。

1.2 裂缝扩展模型建立

图2为干热岩水力压裂几何模型，分别模拟的是常规垂直井和径向井中的压裂过程。模型使用的立方体试样尺寸为0.2 m×0.2 m×0.2 m，在模型中心位置设置垂直井眼以及裸眼段，井眼方向平行于垂向应力，井眼直径为16 mm，长度为150 mm，下部为100 mm的裸眼段。基于前人研究，径向井眼直径和长度分别设置为5和50 mm，方位角（径向井眼轴线方向与最大水平主应力方向之间的夹角）设置为30°。岩石的弹性模量为3.8×10¹⁰ Pa，泊松比为0.24，单轴抗压强度为1.8×10⁸ Pa，抗拉强度为1.1×10⁷ Pa，密度为2650 kg·m^-3，孔隙度为2%，导热系数为3 W·（m·℃）^-1，比热容为950 J·（kg·℃）^-1，使用清水作为压裂液，黏度为1 mPa·s。基于干热岩靶区的现场测试数据开展模拟。中国共和盆地GR1井的水平应力差异系数（K=（σ_H-σ_h）/σ_h）为0.1~0.4，储层温度为218℃^[20-21]。因此在本模拟方案中，岩样初始温度设置为200℃，模型垂向主应力 σ_v为15 MPa，最大水平主应力σ_H为10 MPa，最小水平主应力σ_h为7 MPa，该水平地应力差异系数与真实干热岩储层相符。压裂过程中，注入流体初始温度设置为20℃并以30 mL/min的恒定速率进行连续注入压裂。

通过模拟不同岩石温度（T_r=25、100、200℃）下常规垂直井压裂和径向井压裂裂缝起裂与扩展过程，阐明热应力和径向井应力重构作用对裂缝扩展的定向引导效果。在此基础上，模拟和分析不同地应力差（σ_H-σ_h=0、1.5、3、5、7 MPa，其中σ_H保持恒定）、注入排量（Q=30、90、180、360 mL/min）和压裂液黏度（μ=1、30、50、100、150 mPa·s）条件下的裂缝扩展特征，明确影响定向引导效果的主控因素和最优工程参数范围，为干热岩储层强化裂缝定向沟通提供技术途径和理论基础。

2 模型验证

2.1 水力压裂裂缝扩展验证

通过三维Penny模型裂缝开度随裂缝半长的变化情况，验证离散格子法模拟水力压裂裂缝扩展的准确性。Dontson^[22]提出了断裂韧性控制机制中裂缝开度和裂缝半径的近似解。本文中采用相同的模型设置进行裂缝扩展模拟：建立50 m×50 m×50 m地质模型，垂直井筒位于模型中心，半径为0.2 m，井筒末端预置半径20 m的圆形平面裂缝。流体以速率为0.1 m³/s恒定注入，进行50 s的压裂模拟。图3（a）为水力压裂裂缝开度的分布云图。水力裂缝的开度和半径的近似解^[22-23]为

式中，w_k为裂缝开度，m；R_k为裂缝半径，m；ρ（ρ=r/R）为沿裂缝半径方向正则化坐标；r为计算点与注入点之间距离，m；t为流体注入时间，s；Q₀为流体注入速率，mL/min。

图3（b）为裂缝开度随裂缝半长的变化情况，其中数值解通过模型计算获得。结果表明，数值模拟所得裂缝开度随注入点位置的变化关系与解析解吻合较好，验证了本文中所使用数值模拟方法的准确性。

2.2 热流固耦合压裂模型验证

基于前人开展的高温花岗岩水力压裂物理模拟试验结果，验证热流固多场耦合压裂模型的准确性。参考Yang等^[24]的相关研究，采用相同的岩石物性参数、边界条件和初始条件进行压裂模拟：模型尺寸为0.2 m×0.2 m×0.2 m，模型中心设置长度为130 mm、直径为16 mm的垂直井筒，底部预留60 mm的裸眼段注入清水进行压裂。岩石温度为200℃，注入流体温度为20℃，注入速率为30 mL/min。模型中垂向主应力σ_v为22.4 MPa，最大水平主应力σ_H为15 MPa，最小水平主应力σ_h为11.2 MPa。从图4可以看出，在离散格子模型和水力压裂试验中，主裂缝均沿最大水平主应力方向扩展，整体裂缝形态可描述为“热应力裂缝（井周热损伤区）+主裂缝+分支裂缝”。同时，数值模型和压裂试验的破裂压力相近，分别为29.8和25.8 MPa。结果表明，采用离散格子法建立的热-流-固多场耦合压裂模型与物理模拟压裂试验结果吻合较好，验证了数值模型及其裂缝扩展的准确性。

3 裂缝起裂与扩展数值模拟

3.1 裂缝起裂与扩展特征

为阐明高温条件下径向井应力重构对裂缝扩展特征的影响，开展径向井压裂和垂直井压裂数值模拟研究。图5为不同时刻下干热岩压裂裂缝形态的演化情况。①在裂缝起裂初期（t=0.3 s），高温岩石由于低温压裂液的注入诱发热应力，在井筒周围形成热损伤区并产生大量微裂缝；此外，由于径向井的存在，注入流体与高温岩石接触面积更大，导致微裂缝更多；②在裂缝扩展中期（t=0.9 s），径向井压裂裂缝形态受径向井眼控制，表现为沿径向井眼方向定向扩展；而垂直井压裂裂缝形态受热应力影响，形成沿多个方向扩展的多条裂缝；③在裂缝扩展后期（t=1.5 s），径向井压裂裂缝扩展行为从最初的径向井控制阶段转变为地应力主控阶段，导致裂缝方向出现扭转，最终沿最大水平主应力方向扩展。同样，垂直井压裂裂缝受热应力影响逐渐减弱，地应力开始主导裂缝延伸，即扩展方向近似平行于最大水平主应力。图6为岩石温度200℃、地应力σ_v/σ_H/σ_h=15/10/7 MPa不同时刻干热岩的温度场、位移场和应力场分布演化情况。观察到如下变化：①随着压裂时间的增加，注入流体与高温岩石发生热传递，井筒周围产生温度变化，但由于岩石基质渗透性极低，裂缝附近并未出现明显温降；②相较于垂直井压裂，径向井周围拉应力显著增强，位移变化更为明显，且在增压作用下径向井上下两侧形成最大主应力区，促进裂缝沿径向井轴线方向定向起裂和扩展；③裂缝自径向井根部起裂并以增长缝高的形式纵向延伸，当裂缝扩展超出临界点进入原始地应力分布区域后，此时受水平地应力差的影响开始转向，最终裂缝面趋向于沿最大水平主应力方向。综上所述，径向井压裂是基于重构的应力场，注水增压，使径向井井周拉应力优先达到岩石破裂强度，引导压裂裂缝沿径向井扩展，促使裂缝深穿透，最终形成定向裂缝面。

干热岩压裂井筒压力-时间曲线特征如图7所示。径向井压裂和垂直井压裂过程中，井筒压力随流体注入缓慢上升，达到岩石破裂压力后逐渐降低并趋于平衡。由于径向井条件下的井筒体积更大，其压力曲线的增压速率较低。在相同条件下，径向井压裂和垂直井压裂的破裂压力分别为19.9和22.7 MPa。相比之下，径向井压裂更容易促进岩石的起裂与扩展，破裂压力降低了12.3%。这主要是因为冷流体与高温岩石间的巨大温差诱导热应力，热应力表现为井筒周围的拉应力，使与低温流体接触面积更大的径向井压裂更容易使高温岩样达到拉伸破坏状态。

3.2 温差诱导热应力对裂缝扩展的影响

为阐明温差诱导热应力对裂缝扩展模式与应力场重构特征的影响，对比分析不同岩石温度条件下径向井压裂和垂直井压裂的裂缝形态和应力场分布特征（图8）。地应力σ_v/σ_H/σ_h为15/10/7 MPa，压裂液注入排量为30 mL/min，总注入时间为1.5 s。可以发现，热应力对水力裂缝形态具有重要影响。对于径向井压裂，当岩石温度为25℃时，注入流体与岩石温度接近，井筒周围未出现因热应力而引起的微裂缝，裂缝形态仅受径向井影响，主井筒两侧裂缝沿径向井方向定向扩展；当岩石温度为100℃时，注入流体与岩石的温差增大，形成“主裂缝+热应力裂缝（井周热损伤区）”的裂缝形态，主裂缝在扩展过程中受地应力的影响发生转向；当岩石温度为200℃，注入流体与岩石的巨大温差诱导产生更多的微裂缝（井周热损伤区范围增加），裂缝形态受地应力和径向井的共同影响，最终形成“多主裂缝+热应力裂缝”的复杂缝网，其中主裂缝引导能力显著降低，径向井一侧主裂缝定向扩展一段距离后逐渐转向最大水平主应力方向，无径向井一侧主裂缝扩展平行于最大水平主应力方向。对于垂直井压裂，当岩石温度为25℃时，裂缝形态仅受地应力影响，呈现出近似平行于最大水平主应力的双翼裂缝，形态单一；当岩石温度为100℃时，裂缝形态受热应力和地应力共同影响，表现为“主裂缝+热应力裂缝”；当岩石温度为200℃时，热应力的影响增强，井筒周围热破裂程度显著，造成多条主裂缝沿不同方向起裂。综上，对于垂直井压裂，裂缝形态受地应力主导控制，但低温流体与高温岩石之间的温差诱导热应力可以显著改变井周应力分布，促进形成复杂缝网。对于径向井压裂，裂缝形态受径向井主导控制，温差诱导热应力越大，径向井井周应力重构区域越小，导致裂缝定向引导能力显著降低。

3.3 水平地应力差对裂缝扩展的影响

为探究地应力对裂缝定向扩展特征的影响，对比分析不同水平地应力差（σ_H保持恒定）条件下径向井压裂裂缝形态及应力场分布特征（图9）。岩石温度为200℃，压裂液注入排量为30 mL/min，总注入时间为1.5 s。结果表明，当水平地应力差为0 MPa（应力差异系数为0）时，径向井压裂裂缝形态复杂，其中井周岩石受热应力作用产生大量微裂缝，多条主裂缝沿径向井方向定向起裂并继续扩展，同时伴有一条随机方向的主裂缝；当水平地应力差为1.5 MPa（应力差异系数为0.2）时，主裂缝数目明显减少，仅诱导产生一条沿径向井方向持续扩展的主裂缝；当水平地应力差增大到5和7 MPa（应力差异系数大于1）时，地应力影响占据主导地位，主裂缝在径向井根部起裂后迅速转向最大主应力方向。因此地应力差是影响径向井压裂裂缝扩展模式的关键参数。在较低的地应力差（应力差异系数为0~0.4）条件下，径向井井周应力重构区域更明显，有利于实现更优的裂缝定向引导效果。

3.4 压裂液注入排量及黏度对裂缝扩展的影响

压裂液注入排量和黏度是储层改造中重要的可控施工参数，决定着压后缝网的形态特征。基于热流固耦合压裂模型，分析不同排量条件下的径向井压裂裂缝形态及应力场分布特征（图10）。地应力σ_v/σ_H/σ_h为15/10/7 MPa，岩石温度为200℃，保证流体总注入量一致。可以发现，随着压裂液注入排量的增大，径向井井周应力重构区域逐渐增加，利于引导水力裂缝沿径向井方向扩展。同时，较高的注入排量促使井周产生更多不受地应力控制的裂缝，主裂缝数量增加，裂缝网络形态更为复杂。其主要原因在于，高排量导致压裂液在井筒内快速憋压，促使井周多条裂缝沿不同方向同时起裂。值得说明的是，在径向井压裂条件下，提高注入排量可以增强径向井对裂缝扩展的引导能力，但其改善幅度有限，裂缝数量和复杂程度的增加更为显著。

图11为不同压裂液黏度条件下径向井压裂裂缝形态及应力场分布特征。可以发现，压裂液注入排量和黏度对裂缝扩展模式具有相似影响。随着压裂液黏度的增大，井周产生更多不受地应力控制的裂缝，形成复杂的裂缝网络；但由于压裂液黏度大，造成裂缝扩展速度慢，裂缝几何形态相对较小。综合考虑，30~50 mPa·s的中等黏度压裂液既有利于增强径向井对裂缝的引导能力，又有助于形成复杂的裂缝网络。因此为实现干热岩储层高效改造，建议压裂现场施工优选大排量、中等黏度的压裂液，以最大程度提升裂缝定向扩展距离和储层改造复杂程度。

4 结论

（1）对于垂直井压裂，裂缝形态受原地应力场主导控制，主要沿最大主应力方向扩展。而对于径向井压裂，裂缝形态受径向井井周重构应力场影响，表现为沿径向井眼方向定向起裂；随着裂缝进一步扩展，径向井影响减弱，裂缝扩展路径受原地应力场影响增强，导致裂缝转向，近似沿最大水平主应力方向扩展。

（2）径向井应力重构和热应力耦合作用下裂缝扩展模式：压裂液与岩石温差越大，热应力越强，利于形成“主裂缝+热应力裂缝”的复杂缝网，但径向井对裂缝的引导作用明显减弱。

（3）水平地应力差较小（应力差异系数为0~0.4）时，径向井井周应力重构区域大，利于形成多条定向扩展的主裂缝；随着水平地应力差增大（应力差异系数大于1），径向井引导裂缝定向扩展能力显著降低，主裂缝在径向井根部起裂后迅速转向最大主应力方向。

（4）增大压裂液注入排量和黏度会增强径向井对裂缝的引导作用，且裂缝复杂程度显著增大，在实际工程中，建议优选大排量、中等黏度的压裂液，以最大程度提升裂缝定向扩展距离和储层改造复杂程度。

参考文献