Thermal stresses in borehole heat exchangers

The backfilling materials of borehole heat exchangers (BHE), particularly the grout material, must provide a suitable thermal contact and ensure durability to the induced thermal stresses because of the heat loading. In this paper, the thermal stresses that occurred in BHEs because of heat injection or extraction is investigated with an analytical solution of a hollow cylinder model that is adapted for time‐dependent heat loading, the geometry of a BHE, and the thermo‐mechanical properties of surrounding ground conditions. Firstly, the hollow cylinder model is solved with the considered boundary conditions in 2D plane stress. Secondly, the temperature differences at the inner and outer circles of the cylinder are evaluated with the heat line source models for continuous and discontinuous loading to observe the impact of the heat loading schedule. The developed analytical solution for thermal stress investigation is validated with numerical models. It is demonstrated that the analytical solutions agree well with numerical results for two types of BHE configurations (co‐axial and single U‐shaped pipes). Furthermore, the calculated maximum stresses are compared with the tensile strength of grout materials obtained from Brazilian tests. It is predicted that the thermal contraction of the grout, partially constrained by the surrounding rock, generates tensile stresses that may lead to cracking in the BHE. According to the results, the stiffness of rock has a primary role on the developed tensile stresses, and the relationship between the thermal conductivity of the ground and of the grout induces a proportional impact on the magnitude of thermal stresses. Copyright © 2015 John Wiley & Sons, Ltd.     

分布热源作用下单裂隙岩体三维水流-传热的半解析计算方法

针对高放射性核废物地下处置库近场饱和裂隙岩体环境,提出一种由分布热源、饱和单裂隙和两侧无限大岩石构成的三维水流-传热简化模型,建立了控制微分方程和基于拉氏变换域格林函数的积分方程;采用矩形单元把裂隙面域离散化,利用极坐标下的解析方法计算包含奇点的单元积分,利用数值方法计算分布热源和不包含奇点的单元积分,建立拉氏变换域的线性代数方程组,求解后,利用拉氏数值逆变换,计算任意时刻裂隙水和岩石的温度分布。对两个无内热源、流场确定的计算模型进行了计算,与仅考虑岩石沿裂隙面法向一维热传导的解析解进行了对比。计算分析了分布热源作用下饱和单裂隙岩体的三维水流-传热特征及其对裂隙水流速、岩石热传导系数和热源热流集度的敏感度。计算结果表明：与直接采用高斯数值积分相比,提出的解析法奇异积分精度较高;就裂隙水温度而言,单裂隙岩体三维水流-传热半解析计算方法与解析法得到的结果基本一致,但由于半解析计算方法考虑了岩石的三维热传导,使得裂隙水的上游温度较低,而下游温度较高;无分布热源作用时,岩石热传导系数越大,裂隙水温度越低;裂隙水流速越大,裂隙进水温度对裂隙水和岩石温度分布的影响越明显;由于受到裂隙水流动传热的作用,分布热源对裂隙水温度和岩石温度的影响在裂隙水流的下游区域比较显著。     

分布热源作用下裂隙岩体渗流-传热的拉氏变换-格林函数半解析计算方法

以裂隙岩体高放射性核废物地下处置库性能评估为目标,提出了分布热源作用下单裂隙岩体渗流-传热的简化概念模型、控制微分方程和拉氏变换-格林函数半解析法,为进一步采用半解析法计算分布热源作用下多裂隙岩体的渗流-传热问题奠定了基础。针对单裂隙岩体的渗流-传热问题,建立考虑岩石内热源和二维热传导的控制微分方程,利用拉氏变换域微分方程的基本解建立格林函数积分方程,采用解析法处理其中的奇点,通过数值积分和拉氏数值逆变换求解,计算任意时刻裂隙水和岩石的温度分布。通过算例,与基于岩石一维热传导假定的解析解进行了对比,并计算分析了分布热源作用下单裂隙岩体的渗流-传热特征及其对裂隙开度、岩石热传导系数和热流集度的敏感度。算例表明,（1）就裂隙水温度而言,由于考虑了岩石的二维热传导,拉氏变换-格林函数半解析解小于基于岩石一维热传导假定的解析解;（2）裂隙水温度和岩石温度对裂隙开度和热流集度的敏感度较大,对岩石热传导系数的敏感度较小。     

稀疏不规则裂隙岩体模型三维水流-传热对温度和应力影响的试验研究

选取高放射核废物处置库重要预选场区甘肃北山地区的花岗岩,制作750 mm（宽）×300 mm（厚）×1 000 mm（高）的稀疏不规则裂隙岩体模型,该模型由18块花岗岩和竖向与斜向各两条裂隙组成,在裂隙及岩石内部埋置温度传感器、水压计、直角应变花,并在模型一侧设置局部热源,研究热源温度和裂隙水流速对岩石温度和应力的影响。结果表明,竖裂隙水主要从顶部进水口流向底部出水口,斜裂隙水主要从侧部进水口流向侧部出水口,竖裂隙与斜裂隙在交汇处存在微小流量交换;由于热源处在两条斜裂隙进水口之间,并且斜裂隙长度小于竖裂隙,岩石热传导与斜裂隙水流对岩石温度分布起控制作用,竖裂隙水流对岩石横向热传导起阻滞作用;由于热传导和水流传热的不规则性,上层岩石形成从左向右为主的传热路径,中层和下层岩石形成从上向下为主的传热路径;由于上、下层岩石温度梯度较小,岩石收缩受热拉应力,而中层岩石温度梯度较大,岩石膨胀受热压应力,大主应力的方向大致垂直于斜裂隙面与竖裂隙面的交线,岩石应力增量随斜平面方向的温度梯度增大而增大;热源温度越高,裂隙水流速越低,岩石温度越高、岩石应力越大,系统达到稳态需要的时间越长。     

饱和稀疏裂隙岩体三维水流-传热过程中位移和应力的一种半解析计算方法

针对高放射核废深地质处置库近场环境,建立分布热源作用下饱和裂隙岩体三维水流-传热过程中位移和应力的一种半解析计算方法：采用Goodier热弹性位移势和Laplace变换计算由温度梯度产生的温梯位移和应力;考虑单一裂隙的情况,利用经典弹性力学的Boussinesq解和Cerruti解计算为满足边界条件的约束位移和应力,与温梯位移和应力叠加,可得总体热位移和应力;把裂隙面离散为矩形单元集合,采用极坐标系下的解析法计算包含奇点的单元积分,采用数值法计算与分布热源有关和不含奇点的单元积分。与基于裂隙面法向一维热传导假设的一种解析解对比,结果表明,半解析法与解析法的计算结果基本一致,但由于半解析法考虑岩石的三维热传导,因温度时空分布和演变的不同而导致不同的温梯应力。针对一个假想单裂隙岩体三维水流-传热过程,计算温梯位移和应力、约束位移和应力、总体位移和应力;结果表明,裂隙水流-传热可能对位移和应力的分布和演变有显著影响,距离分布热源较近的岩石因升温膨胀受到约束而出现压应力,而距离分布热源较远的岩石则可能因协调收缩受到约束而出现拉应力。     

Numerical 3D modeling of the effects of strain rate and confining pressure on the compressive behavior of Kuru granite

This paper deals with numerical modeling of the compressive behavior of granite rock under high strain rate dynamic loading and wide range of confining pressure. For this end, a constitutive model based on damage mechanics and viscoplasticity for rock is formulated and implemented in explicit dynamics FEM. The viscoplastic part is based on a simple power law type yield criterion that incorporates the rate-dependency with a linear viscosity term. Moreover, a Rankine type of tensile cut-off is employed. The damage part of the model is formulated with separate scalar damage variables in tension and compression. The model is calibrated for Kuru granite and validated with the experimental data from dynamic compression tests at the strain rate of 600 1/s up to 225 MPa of confining pressure. The numerical simulations demonstrate that, despite the underlying continuum modeling approach, the model captures the correct experimental failure modes, including the transition from single-to-multiple fragmentation, as well as the dynamic compressive strengths at different confining pressures.     

Dynamic Tensile Failure of the Rock Interface Between Tuff and Basalt

The dynamic tensile strength properties of the rock interface and its host rocks sampled from the Baihetan Hydropower Station from Western China were measured using a split Hopkinson pressure bar (SHPB). The results were compared with those for its two host rocks. The dynamic tensile strengths of the two host rocks, tuff and basalt have typical loading rate dependence. However, the dynamic response of the rock interface is much more complicated and at a given loading rate, varies between those of tuff and basalt. To explain the observation, numerical simulation using the discrete element method (DEM) was conducted to determine the detailed tensile failure process of the rock interface. The numerical simulation verifies that the variation of the dynamic tensile strength of the rock interface is a result of the variation of the interface geometry.     

Numerical modeling of dynamic rock fracture with a combined 3D continuum viscodamage‐embedded discontinuity model

This paper deals with numerical modeling of dynamic failure phenomena in rate‐sensitive quasi‐brittle materials, such as rocks, with initial microcrack populations. To this end, a continuum viscodamage‐embedded discontinuity model is developed and tested in full 3D setting. The model describes the pre‐peak nonlinear and rate‐sensitive hardening response of the material behavior, representing the fracture‐process zone creation, by a rate‐dependent continuum damage model. The post‐peak response, involving the macrocrack creation accompanied by exponential softening, is formulated by using an embedded displacement discontinuity model. The finite element implementation of this model relies upon the linear tetrahedral element, which seems appropriate for explicit dynamic analyses involving stress wave propagation. The problems of crack locking and spreading typical of embedded discontinuity models are addressed in this paper. A combination of two remedies, the inclusion of viscosity in the spirit of Wang's viscoplastic consistency approach and introduction of isotropic damaging into the embedded discontinuity model, is shown to be effective in the present explicit dynamics setting. The model performance is illustrated by several numerical simulations. In particular, the dynamic Brazilian disc test and the Kalthoff–Winkler experiment show that the present model provides realistic predictions with the correct failure modes and rate‐dependent tensile strengths of rock at different loading rates. The ability of initial embedded discontinuity populations to model the initial microcrack populations in rocks is also successfully tested. Copyright © 2016 John Wiley & Sons, Ltd.     

Numerical modeling of rock failure under dynamic loading with polygonal elements

Polygonal finite elements are gaining an increasing attention in the computational mechanics literature, but their application in rock mechanics is very rare. This paper deals with numerical modeling of rock failure under dynamic loading based on polygonal finite elements. For this end, a damage-viscoplastic constitutive model for rock based on the Mohr-Coulomb criterion with the Rankine criterion as a tensile cutoff is employed and implemented with the polygonal finite element method. Moreover, the mineral mesostructure or rock is described by randomly mapping groups of polygonal elements representing the constituent minerals into a global mesh and assigning these groups with the corresponding mineral material properties. The performance of the polygonal elements is compared with that of the linear and quadratic triangular and bilinear quadrilateral elements in numerical simulations of controlled shear band formation under uniaxial compression and lateral splitting failure in the dogbone tension test. Numerical simulations of uniaxial tension and compression tests as well as dynamic Brazilian disc test under increasing loading rates demonstrate that the present approach predicts the correct failure modes as well as the dynamic increase in strength of rock.     

裂隙岩体不稳定温度场的复合单元算法研究

基于复合单元法,结合三维热传导-对流方程和“充填模型”,提出了裂隙岩体不稳定温度场的复合单元模型。该模型前处理简便快捷,计算网格生成时无需考虑裂隙的存在,网格剖分不受限制,随后利用复合单元前处理程序,依据裂隙的位置和方位将其自动离散在单元内。对常规热传导-对流方程进行自伴随性调整,应用变分原理,推导出裂隙岩体不稳定温度场的复合单元算法,该算法可分别计算出岩块子单元和裂隙的温度值,且可真实反映裂隙中水流与相邻岩块间的热能量交换规律。将复合单元数值模型计算的不稳定温度场结果与相应的实测数据进行对比分析可知,数值计算结果与实测数据基本一致,验证了裂隙岩体不稳定温度场复合单元算法的可靠性与有效性。算例分析表明,裂隙中水流与相邻岩块间有明显的热传导和热对流作用。     

FDEM-TH3D: A three-dimensional coupled hydrothermal model for fractured rock

This paper proposes a three-dimensional coupled hydrothermal model for fractured rock based on the finite-discrete element method to simulate fluid flow and heat transport. The 3D coupled hydrothermal model is composed of three main parts: a heat conduction model for the rock matrix, a heat transfer model for the fluid in the fractures (including heat conduction and heat convection), and a heat exchange model between the rock matrix and the fluid in the fractures. Four examples with analytical solutions are provided to verify the model. A heat exchange experiment of circulating water in a cylindrical granite sample with one fracture is simulated. The simulation results agree well with the experimental results. The effects of the fracture aperture, fluid viscosity, and pressure difference on the heat exchange between the fluid and rock are studied. Finally, an application concerned with heat transport and fluid flow in fractured rock is presented. The simulation results indicate that the 3D fully coupled hydrothermal model can capture the fluid flow and temperature evolution of rocks and fluids.     

米尺度裂隙岩体模型水流-传热试验的数值模拟分析

为了研究高放射性核废物地下处置库近场的水流-传热耦合问题,采用国内高放废物地下处置库预选场址--甘肃北山地区的花岗岩石块体,加工组合成米尺度的规则裂隙岩体模型,设置边界热源和裂隙水流,试验模拟裂隙水水流与传热之间的相互作用。作为该室内模型试验的前期理论研究,采用等效孔隙介质数值模型,着重分析了裂隙开度、裂隙流量和热源功率对流场和温度场的影响。在设定条件下,计算分析表明：热传导和裂隙水水流由热源作用初期的不耦合很快转化为耦合;不流动的裂隙水主要表现为热存储和热传导,而流动的裂隙水还引起流动传热和水与岩石之间的对流换热,使岩体温度场明显不同于单纯热传导的情况;如果保持裂隙水流量不变,则裂隙开度的变化对水流-传热影响不大;如果保持裂隙水流速不变,则裂隙开度的变化对水流-传热影响显著;热源功率越大,通过裂隙水的热流量越大,裂隙水压强越大,而当温度超过100 ℃时,裂隙水会因汽化而压强显著增大;加热7 d时,热量的输入和输出几乎相等,裂隙水流带走的热量接近热源供给的热量,模型系统基本达到了热平衡。     

Experimental and Numerical Studies on Determination of Indirect Tensile Strength of Rocks

Indirect tension tests using Brisbane tuff Brazilian disc specimens under standard Brazilian jaws and various loading arcs were performed. The standard Brazilian indirect tensile tests caused catastrophic, crushing failure of the disc specimens, rather than the expected tensile splitting failure initiated by a central crack. This led to an investigation of the fracturing of Brazilian disc specimens and the existing indirect tensile Brazilian test using steel loading arcs with different angles. The results showed that the ultimate failure load increased with increasing loading arc angles. With no international standard for determining indirect tensile strength of rocks under diametral load, numerical modelling and analytical solutions were undertaken. Numerical simulations using RFPA^2D software were conducted with a heterogeneous material model. The results predicted tensile stress in the discs and visually reproduced the progressive fracture process. It was concluded that standard Brazilian jaws cause catastrophic, crushing failure of the disc specimens instead of producing a central splitting crack. The experimental and numerical results showed that 20° and 30° loading arcs result in diametral splitting fractures starting at the disc centre. Moreover, intrinsic material properties (e.g. fracture toughness) may be used to propose the best loading configuration to determine the indirect tensile strength of rocks. Here, by using numerical outcomes and empirical relationships between fracture toughness and tensile strength, the best loading geometry to obtain the most accurate indirect tensile strength of rocks was the 2α?=?30° loading arc.     

A three-dimensional numerical investigation of the fracture of rock specimens containing a pre-existing surface flaw

Three-dimensional surface crack initiation and propagation in two kinds of heterogeneous rocks were numerically investigated via parallel finite element analysis using a supercomputer. Numerically simulated rock specimens containing a pre-existing flaw were subjected to uniaxial compression until failure. The initiation and propagation of wing cracks, anti-wing cracks, and shell-like cracks were reproduced by numerical simulations. The numerically simulated results demonstrate that the further propagation of wing cracks and shell-like cracks stop due to their wrapping (curving) behavior in three-dimensional spaces, even if the applied loads continue to increase. Furthermore, rock heterogeneity could significantly influence crack propagation patterns and the peak uniaxial compressive strengths of rock specimens. Moreover, anti-wing cracks only appeared in relatively heterogeneous rocks, and the peak uniaxial compressive strengths of the specimens were observed to depend on the inclination of the pre-existing flaw. Finally, the mechanism of surface crack propagation is discussed in the context of numerically simulated anti-plane loading tests, wherein it was identified that Mode III loading (anti-plane loading) does not lead to Mode III fracture in rocks due to their high ratio of uniaxial compressive strength to tensile strength. This finding could explain the lateral growth of an existing flaw in its own plane, which is a phenomenon that has not been observed in laboratory experiments.     

不同水力加载条件下非均质储层缝网扩展规律研究

水力压裂可显著提高页岩气等致密储层岩体的渗透性以增加油气产量,然而受多因素影响,水力压裂形成缝网结构的机理和压裂优化设计一直是研究的焦点和难点。本研究基于渗流-应力-破坏耦合计算模拟方法,对不同水力加载条件下的非均质储层水力压裂过程进行了模拟和对比研究。研究结果表明:水力压裂过程中起始注水压力和增量大小对水力压裂缝网扩展和改造区域形态有着显著的影响。当起始注水压力小于等于模型材料体抗拉强度,并缓慢增压致裂时,压裂过程可近似视为稳态应力-破坏-渗流耦合作用过程的不同阶段,这种情况下仅在压裂井孔周围形成两组对称式的伞状水力裂缝带。当对模型体施加高于模型材料体破裂压力的注水压力时,相当于对压裂孔快速施加高动水压力,水力裂缝沿压裂孔全方位迅速萌生并快速扩展,当注水压力值高于破裂压力一定幅值时,压裂改造可形成围绕压裂井全方位的放射状裂缝网络,使压裂储层得以最大范围改造。在拟静力和拟动力两种加载条件下,不同水岩相互作用机理是造成不同水力加载条件出现不同缝网结构的力学机制,而对于实际的页岩气储层改造,压裂产生围绕压裂井全方位放射状的缝网结构则是一种最优的体积压裂改造。     

Effect of inherent microcrack populations on rock tensile fracture behaviour: numerical study based on embedded discontinuity finite elements

Inherent microcrack populations have a significant effect on the fracture behaviour of natural rocks. The present study addresses this topic in numerical simulations of uniaxial tension and three-point bending tests. For this end, a rock fracture model based on multiple intersecting embedded discontinuity finite elements is developed. The inherent (pre-existing) microcrack populations are represented by pre-embedded randomly oriented discontinuity populations. Crack shielding (through spurious locking) is prevented by allowing a new crack to be introduced, upon violation of the Rankine criterion, in an element with an initial crack unfavourably oriented to the loading direction. Rock heterogeneity is accounted for by random clusters of triangular finite elements representing different minerals of granitic numerical rock. Numerical simulations demonstrate the strength lowering effect of initial microcrack populations. This effect is substantially stronger under uniaxial tension, due to the uniform stress state, than in semicircular three-point bending having a non-uniform stress state with a clear local maximum of tensile stress.

     

Damage‐viscoplastic consistency model for rock fracture in heterogeneous rocks under dynamic loading

This paper presents a damage‐viscoplastic consistency model for numerical simulation of brittle fracture in heterogeneous rocks. The model is based on a combination of the recent viscoplastic consistency model by Wang and the isotropic damage concept with separate damage variables in tension and compression. This approach does not suffer from ill‐posedness, caused by strain softening, of the underlying boundary/initial value problem since viscoplasticity provides the regularization by introducing a length scale effect under dynamic loading conditions. The model uses the Mohr–Coulomb yield criterion with the Rankine criterion as a tensile cut‐off. The damage law in compression is calibrated via the degradation index concept of Fang and Harrison. Thereby, the model is able to capture the brittle‐to‐ductile transition occurring in confined compression at a certain level of confinement. The heterogeneity of rock is accounted for by the statistical approach based on the Weibull distribution. Numerical simulations of confined compression test in plane strain conditions demonstrate a good agreement with the experiments at both the material point and structural levels as the fracture modes are realistically predicted. Copyright © 2009 John Wiley & Sons, Ltd.     

A 2D coupled hydro-thermal model for the combined finite-discrete element method

Based on the combined finite-discrete element method (FDEM), a two-dimensional coupled hydro-thermal model is proposed. This model can simulate fluid flow and heat transfer in rock masses with arbitrary complex fracture networks. The model consists of three parts: a heat conduction model of the rock matrix, a heat-transfer model of the fluid in the fracture (including the heat conduction and convection of fluid), and a heat exchange model between the fluid and rock at the fracture surface. Three examples with analytical solutions are given to verify the correctness of the coupled model. Finally, the coupled model is applied to hydro-thermal coupling simulations of a rock mass with a fracture network. The temperature field evolution, the effect of thermal conductivity of the rock matrix thermal conductivity and the fracture aperture on the outlet temperature are studied. The coupled model presented in this paper will enable the application of FDEM to study rock rupture driven by the effect of hydro-thermo-mechanical coupling in geomaterials such as in geothermal systems, petroleum engineering, environmental engineering and nuclear waste geological storage.

     

Numerical Modeling of Pore Pressure Influence on Fracture Evolution in Brittle Heterogeneous Rocks

Rock is a heterogeneous geological material. When rock is subjected to internal hydraulic pressure and external mechanical loading, the fluid flow properties will be altered by closing, opening, or other interaction of pre-existing weaknesses or by induced new fractures. Meanwhile, the pore pressure can influence the fracture behavior on both a local and global scale. A finite element model that can consider the coupled effects of seepage, damage and stress field in heterogeneous rock is described. First, two series of numerical tests in relatively homogeneous and heterogeneous rocks were performed to investigate the influence of pore pressure magnitude and gradient on initiation and propagation of tensile fractures. Second, to examine the initiation of hydraulic fractures and their subsequent propagation, a series of numerical simulations of the behavior of two injection holes inside a saturated rock mass are carried out. The rock is subjected to different initial in situ stress ratios and to an internal injection (pore) pressure at the two injection holes. Numerically, simulated results indicate that tensile fracture is strongly influenced by both pore pressure magnitude and pore pressure gradient. In addition, the heterogeneity of rock, the initial in situ stress ratio (K), the distance between two injection holes, and the difference of the pore pressure in the two injection holes all play important roles in the initiation and propagation of hydraulic fractures. At relatively close spacing and when the two principal stresses are of similar magnitude, the proximity of adjacent injection holes can cause fracturing to occur in a direction perpendicular to the maximum principal stress.     

Finite element simulations of 3D planar hydraulic fracture propagation using a coupled hydro-mechanical interface element

Two-dimensional hydraulic fracturing simulations using the cohesive zone model (CZM) can be readily found in the literature; however, to our knowledge, verified 3D cohesive zone modeling is not available. We present the development of a 3D fully coupled hydro-mechanical finite element method (FEM) model (with parallel computation framework) and its application to hydraulic fracturing. A special zero-thickness interface element based on the CZM is developed for modeling fracture propagation and fluid flow. A local traction-separation law with strain softening is used to capture tensile cracking. The model is verified by considering penny-shaped hydraulic fracture and plain strain Kristianovich‑Geertsma‑de Klerk hydraulic fracture (in 3D) in the viscosity- and toughness-dominated regimes. Good agreement between numerical results and analytical solutions has been achieved. The model is used to investigate the influence of rock and fluid properties on hydraulic fracturing. Lower stiffness tip cohesive elements tend to yield a larger elastic deformation around the fracture tips before the tensile strength is reached, generating a larger fracture length and lower fracture pressure compared with higher stiffness elements. It is found that the energy release rate has almost no influence on hydraulic fracturing in the viscosity-dominated regime because the energy spent in creating new fractures is too small when compared with the total input energy. For the toughness-dominated regime, the released energy during fracturing should be accurately captured; relatively large tensile strength should be used in order to match numerical results to the asymptotic analytical solutions. It requires smaller elements when compared with those used in the viscosity-dominated regime.