对包括岩石在内的现实介质物理结构的真实描述大致是数值分析工作者心中的诗和远方，或许始终遥远而不可及，也或许像圆周率精确值一样无解。尽管软件研发及其硬件技术的发展一贯高歌猛进，但现实要求对科技升级似乎是一道永远无法破解的魔咒，计算机模拟应用仍然客观且普遍存在计算效率瓶颈问题，在岩土工程数值分析领域同样如此。

以颗粒流方法水压致裂研究为例，本期短文剖析模型力学特点，从计算效率角度讨论现场大尺度模型细观力学参数的研究方法。主要内容：

1.      颗粒流方法压裂分析模型的力学特点；

2.      工程尺度颗粒流压裂分析模型参数标定原理；

3.       工程尺度压裂模型参数研究建议方案。

 

颗粒流方法压裂分析模型的力学特点

颗粒流方法采用刚性颗粒及颗粒间的可变形接触来描述岩体，该细观尺度意义上的表达方式可以更好地反映岩石的物理结构特征，因此特别适用于开展复杂岩石力学特性相关研究，典型如岩石细观破裂萌生直至出现宏观破坏现象的力学过程。

对于颗粒间可变形接触力学性质的定义是颗粒流模型分析中特别重要的环节之一。接触力学性质采用接触模型来描述，自Potyondy于2004年提出BPM（bonded particle model）模型以来，该接触模型被广泛地用于岩石破裂特性研究中，典型如水压致裂分析。到目前为止，Itasca公司已以该模型为基础，分别基于PFC环境和XSite环境研发形成两套水压致裂分析方法，其中XSite是基于PFC方法向水压致裂专题研究的定制化延伸。

Griffith针对脆性介质（玻璃）开展了大量的试验与分析研究，其成果为断裂力学典型了基本理论基础。按断裂力学定义，脆性材料中的破裂按力学机制可分为张型（KI）、剪型（KII）和以上两种机制共同存在的复合型（KIII）三种。在水压致裂法注压过程中，水压作用导致裂隙会沿着致裂所需能量最小的方向即在与最小主应力垂直的平面内起裂并扩展，属于典型的KI型机制。其具体的断裂力学意义是，水压作用导致在岩石原生裂隙尖端部位出现应力集中，当尖端的应力强度因子达到岩石断裂韧度KIC时，岩石产生新破裂、进而随尖端应力调整可进一步扩展，表现为一个动态化的力学过程。由此可见，岩石断裂韧度KIC是水压致裂分析研究需特别考察的重要参数。

断裂韧度KIC影响因素较为复杂，主要包括岩石的结构特征（主要是岩石中原生细观破裂长度）、抗拉强度，此外还可受到应力状态影响。Potyondy最早采用BPM建立规则化颗粒模型开展了岩石破裂研究，认为岩石断裂韧度与模型中的颗粒接触强度和尺度条件满足下式定义：

(1)

式中，为颗粒间法向抗拉强度，是模型的平均颗粒半径。

由式（1）可见，断裂韧度KIC在颗粒流模型中与颗粒大小有关，具有尺寸效应性质，这一现象与断裂力学对断裂韧度的定义一致，即岩石断裂韧度与原生细观破裂的几何长度密切相关。

Potyondy采用BPM模型针对含预制裂隙的岩样破裂条件和特征展开了大量模型分析，同时对比了Contact-Bond模型与BPM两种模型条件下岩石破裂特征的差异，成果见图1。其中，试样长、宽均为W，a、R分别为裂隙半长和颗粒平均半径，参数Ψ为预置裂隙半长a与颗粒直径2R之比，反映了模型对岩石细观结构（或原生破裂）的解译度（精度）。因此，在固定颗粒半径参数R的条件下，通过改变裂隙半长a参数可帮助考察岩石断裂韧度KIC的尺寸效应性质。

 

图1:BPM模型数值试验

图1中，横轴描述了解译精度参数Ψ，纵轴表示经归一化处理（预置裂隙试样抗拉强度与颗粒接触抗拉强度之比）得到的预置裂隙岩样的抗拉强度。模型成果主要揭示以下特点：

1)      当解译度Ψ小于10时，岩石抗拉强度随Ψ降低表现出愈发强烈的尺度（解译度）无关性；或者说，在尺度意义上，岩石抗拉强度的塑性性质越来越明显；

2)       当解译度Ψ超过10时，岩样抗拉强度与解译精度近似线性相关，且曲线斜率近似趋近于-0.5，这与线弹性断裂力学LEFM对预置裂隙岩体的抗拉强度的描述特点基本一致。其指示的断裂力学意义具体在于，当Ψ超过10时，颗粒半径近似代表了岩石中原生细观破裂的长度，岩石断裂韧度条件KIC也因此得到更为准确的反映；

3)      BPM和CB模型均可模拟破裂在拉应力作用下的扩展现象。因BPM接触部位考虑了弯矩作用，材料的抗拉强度条件总体低于CB模型结果。

考虑脆性岩石力学性质的代表性取值，如单轴抗压强度、抗拉强度可分别为100MPa、10MPa，断裂韧度条件一般介于0.5~3MPa.m1/2之间。Potyondy、Damjanac等进一步研究指出，为满足对以上参数取值条件的一致性校核标定，颗粒半径需与岩石细观裂隙尺度近似一致，一般可小至厘米（cm）量级。由此可见，岩石断裂韧度的尺寸效应性质决定了颗粒流破裂模型分析应选用的单元大小（颗粒半径），进而影响模型规模和计算效率。该技术特点一般不会对室内试样尺度的模型分析造成显著效率瓶颈，但对现场尺度模型分析而言，通常需首先解决分析效率不足的问题。

 

工程尺度颗粒流压裂分析模型参数标定原理

考虑KI型岩石破裂萌生与扩展性质主要与断裂韧度条件KIC密切相关，抗压强度、抗拉强度等条件属次要因素，针对工程尺度规模的水压致裂BPM颗粒流模型分析，Potyondy和Damjanac建议不对抗拉强度、抗压强度和断裂韧度这些条件均作一致准确的参数标定要求，而是以正确反映断裂韧度条件为主要目标，起到增大颗粒半径，提升分析效率的效果。

参数标定关键环节和流程：

1)      综合工程尺度和计算资源条件选取某一合理的颗粒尺寸，建立BPM现场岩体模型；

2)     以现场岩体的抗压强度、抗拉强度估值为标定依据，校核获得BPM模型参数，特别是接触强度参数（该组参数为中间成果）；

3)        依据以上参数建立含有预置裂隙的模型试样，对接触强度参数中间成果引入折减系数，并考察其对应的岩石断裂韧度条件。当模型获得的断裂韧度与实际值一致时，选择相应的接触强度参数作为工程尺度水压致裂模型分析的最终输入条件；

4)    开展现场尺度压裂模拟计算。

在平面模型条件下，岩石断裂韧度性质KIC可通过Anderson提出的方法开展模型分析进行确定，模型设置和方法为（图2）：

1)     建立长、宽均为W的矩形数值试样（或宽为W，长宽比为2的长方形试样），并输入考虑一系列经比例化（强度）折减得到的BPM模型接触参数；

2)     在试样中部预置长度为W/3的裂隙。注意，裂隙可通过设置接触刚度和强度为0的颗粒来近似表达；

3)     沿预置裂隙垂向施加拉荷载，对试样开展直拉试验，得到试样的抗拉强度。

     据此，按下式计算确定试样的断裂韧度条件：

(2)

式中，为预置裂隙半长；=2/W，对图2所示试样，取=1/3；为：

(3)

 

图2:基于数值模型确定岩石断裂韧度条件

 

工程尺度压裂模型参数研究建议方案

1)      结合现有勘探成果和经验准则确定现场尺度岩体的变形参数（弹性模量、泊松比）和强度参数（抗压强度、抗拉强度、断裂韧度）；

2)      现场尺度岩体颗粒模型细观参数的标定（依据抗压、抗拉强度）：

a)     依据注浆工艺和计算资源合理确定颗粒模型中的颗粒半径，建立现场尺度试样；

b)        开展单轴抗压、直拉强度试验，标定接触变形参数和强度参数的中间成果。

3)        依据断裂韧度进一步修正接触强度参数：

a)        建立图2含预制裂隙的岩体试样；

b)      对接触强度参数中间成果进行折减，通过直拉试验获得裂隙试样抗拉强度，利用式（2）计算得到断裂韧度；

c)       不断修正接触强度参数，若模型断裂韧度和目标值一致时，即得到最终接触强度参数取值结果。

4)      现场尺度岩体模型，开展压裂计算，进行机制研究及注浆压力等因素影响分析。

 

参考文献：

[1]  Potyondy D O, Cundall P A. A bonded-particle model for rock[J]. International Journal of Rock Mechanics & Mining Sciences, 2004, 41(8):págs. 1329-1364.

[2]   Ivan Gil, Branko Damjanac and Fiona Kwok. Analysis of Multiple Hydraulic Fracturing in Horizontal Wellbores[R]. Itasca Consulting Group, 2010.

[3]   Anderson TL (2005) Fracture mechanics: fundaments and applications. 3rd ed. CRC Press.