基本控制方程研究成果

4.3.1.1储层弹塑性力学的基本变量

描述储层岩石的在应力场作用下的变形可以位移矢量｛μ｝表示：

相应的应力状态可用各点的应力矢量｛σ｝来表示：

应变状态可用各点的应力矢量｛ε｝来表示：

处于构造变化中的储层，其应力与应变、应变与位移以及应力与外力之间存在一定关系。

4.3.1.2储层小变形应变描述

在变形较小的情况下，储层岩石的正应变量和剪切应变量与位移矢量的关系如下：

将正应变和剪应变分量合并为应变张量ε，记作：

在变形较小的轴对称情况下应变位移关系如下

对变形较大的情况，为了描述应变张量，可以引入格林应变张量：

4.3.1.3 储层岩石应力-应变关系

1）储层线性弹性本构关系

对于储层的应力状况，应用对称应力张量σ加以描述：

应力张量包含三个正应力和三个剪应力。

线弹性条件下，储层岩石的应力应变关系或本构方程包含初始应力，应变，热效应：

对于各向同性储层岩石，弹性模量矩阵可以写作：

模量矩阵中E为杨氏模量，v为泊松比。对公式（4-66）求逆：

2）储层弹塑性本构关系

塑性尖端扩展模拟均建立在塑性公设的基础上，本节将重点讨论塑性公设及其推导出来的塑性本构关系。

（1）储层岩石材料弹塑性公设

塑性势函数作为一种经典的力学方法，在弹塑性本构关系研究中得到广泛应用。Mises于1928年首次将弹性势函数作为塑性势函数，并进一步提出将塑性势函数的梯度方向作为塑性形变流动方向，据此建立了经典塑性位势理论。

第一，应力循环功理论

Mises建立经典塑性位势理论之后，德鲁克进一步构建了德鲁克塑性公设，创造性得将Mises提出的塑性位势函数和屈服面函数联系在一起。在此基础上提出了，塑性应变增量方向与屈服面正交的正交流动法则，从理论上完善了Mises的经典塑性位势理论。整个弹塑性应力循环包括4个步骤。

①设弹塑性材料处于任意应变历史中的一个平衡状态，在　状态下，处于平衡，初始应变位于加载面内。

②缓慢加载载荷，使应力状态　σij靠近加载面Fσ；应变εij靠近加载面Fε。

③继续加载载荷，应变状态　σij越过加载面Fσ，变为　　 σij+dσij ；εij越过加载面Fε，应变变为 εij+dεij ，此时发生塑性应变。

④卸载载荷，使应力状态变回　；应变应变状态为。

对于稳定材料，在完成整个弹塑性材料应力循环之后，附加外力所做的功大于等于零：

其实，不论储层岩石材料是不是稳定材料，应力循环功都不能小于0，不然外界就可以在这个循环中从储层岩石中吸取能量。由于弹性应变在整个应力循环中做功为0，所以可以将应力循环做功改写为：

由于塑性应变仅仅发生在材料越过加载面之后的一部分，应力循环的其余部分并不产生塑性功，所以可以将以上应力循环功进一步改写为：

由于　dσij非常小，当时，可以将弹塑性应力循环附加力塑性功改写为：

如果应力循环过程的初始应力状态恰好在加载面上，那么则有：

第二，应变循环功理论

但德鲁克对Mises经典塑性位势理论的改造过于严格，而且不满足热力学第二定律，一般适用于储层岩石材料的稳定加载段。尹留辛在德鲁克公设的基础上，基于严格的热力学第二定律，采用应变循环理论。即假设在弹塑性材料的一个应变循环中，外部作用的功是非负的，如果做功是正的，则存在塑性变形，如果应变循环中，外力做功为零，则不存在塑性变形。整个塑性应变循环包括以下4个步骤。

①设弹塑性材料处于任意应变历史中的一个平衡状态，在　状态下，处于平衡，初始应变　位于加载面内。

②缓慢加载载荷，使应变εij靠近加载面Fε；应力状态σ　ij靠近加载面。

③继续加载载荷，εij越过加载面Fε，应变变为　εij+ dεij ，此时发生塑性应变；应变状态σ　ij变为 σij+ dσij。

④卸载载荷，使应变状态变回，应力状态此时变为，产生塑性残余量　，塑性应变量　，其中D为弹性模量矩阵。

根据热力学第二定律，完成弹塑性应变循环后，外部做功要大于等于零，可记做：

对比弹塑性储层岩石材料的应力循环和应变循环公设可以看出，德鲁克公设的假设条件过于严苛因而适用范围为稳定加载段，而伊留辛公设的假设条件较宽泛，在稳定加载段和非稳定加载段均适合。

同时观察应力循环功图和应变循环功图，可以看出可以利用，即塑性残余量　来表示。应力循环功理论采用塑性残余量　可以重新表述为：

对比应力循环功和应变循环功，可以看出应变循环功比应力循环功多了图中下部的一个小三角形的面积即。所以一般有应变循环功大于应力循环功。同样将应变循环功写成代数形式有：

图50 弹塑性循环功非负公设

（a：弹塑性应力循环稳定性公设；b：应变循环通用性公设）

（2）重叠空间中屈服面和塑性势的分析

根据应力循环功理论，在应力空间中考虑问题，任意一点的应力状态可以用矢量表示，同时任意一点的应变状态也可以用矢量表示。如果将应力空间和塑性应变空间重合，在这个重叠坐标空间中进行分析，可以得到塑性应变增量和屈服面之间的关系。为了分析塑性应变的增量方向与应力空间中屈服面的关系，将塑性应变空间的原点放在屈服面上，同时让塑性应变的坐标和应力空间的应力坐标 重合。假定初始应力在应力空间中的坐标为　，其在应力空间中的矢量表示为 。当应力加载至应力空间屈服面时，应力在应力空间中的坐标为　σij ，此时应力在应力空间中的坐标为。则应力的增量，在应变空间中采用矢量可以表示为　。根据应力循环功理论　，可知　。 由于我们采用重叠空间对应力和塑性应变进行分析，可以得到如下形式的应力循环 功理论表示。其中为重叠空间中塑性应变的矢量表示，根据内积定律改写应力循环功理论为。

其中θ为重叠空间中应力增量和塑性增量的夹角，由上述该关系式可知，在重叠空间中应力增量和塑性增量的夹角小于等于0。由于该结果对于应力空间屈 服面以内的任意初始应力均成立。在重叠空间屈服面上任意取一点　ijσ，做切平面T （σij ），同时取重叠空间中塑性应变增量矢量为，则以　σij为原点，在切平面T （σij）周围做无限小圆周γ上取到　，可知在此圆周上均有。在γ上，如 果，存　在　一　点　令，则取与公共面P内的反向延长线上一点　，令，则与应力循环功理论相矛盾，因此在γ上，均有，即塑性应变增量垂直于重叠 空间中的应力屈服面。

图51 塑性应变增量矢量和应力空间屈服面切面垂直示意图

根据应力循环功理论，在重叠空间中任意初始应力　与应力屈服面上　ijσ点之间的向量与　σij点之处的塑性应变增量之间的夹角θ均应当不大于90度。同时根据上述推导，可知在重叠空间中塑性应变增量矢量平行于重叠空间中的应力屈服面上　σij点之处法向矢量。所以与　σij点之处法向矢量夹角应当不大于90度。所以所有初始应力点　均应在屈服面Fσ以内，并且也在屈服面上点　σij的切平面T （σij）内侧。这说明应力循环功理论要求，重叠空间中 应力屈服面的几何形状必须在屈服面上任意点切平面的内侧。这点将有助于后面分析SEFE方法中所使用的屈服面的形态特征。

由于塑性应变增量垂直于重叠空间中的应力屈服面，所以在重叠空间中塑性应变增量所包含的3个塑性主应变增量均与重叠空间中应力屈服面的法向保持平行。由矢量分析的基本理论，可知这将导致塑性应变增量的3个分量与应力屈服面法向3个分量差一个乘数因子λ。也就说说塑性应变增量和屈服面法向对应分量有一个统一的比例，这就将塑性应变增量三个分量联系起来。但是必须明确的 是，这个比例的存在，是由于存在一个唯一的应力屈服面Fσ，同时是在应力循环 功理论假设下，才得到的塑性应变增量矢量方向与屈服面垂直的结论。如果在重 叠空间中不存在一个唯一的应力屈服面Fσ，也就是说主应力空间中3个分量有更为复杂的屈服特性，无法在主应力空间中找到一个唯一的描述面Fσ去表示屈服面， 则没有该结论。

图52 重叠空间中应力屈服面的形态特征分析示意图

（3）对公设的批判

根据上面的分析，“塑性应变增量矢量和应力空间屈服面切面垂直”和“重叠空间中应力屈服面位于其上任意点切平面的内侧”两个结论均假定了3个前提条件：第一，应力循环功理论正确；第二，存在唯一塑性屈服面；第三，塑性势面和屈服重合。有了这3个前提假设，才推导出来的前面两个关于屈服面特征的两个结论。目前这3个假设也是被普遍认同和使用的，但是本节将进一步分析这3个假设的合理性。

首先讨论第一个假设，即应力循环功理论正确。由前面的推导可知，应力循 环功理论的一般表述为。如果没有初始应力项，则应力循环功 代表真实应力作用功，但其减去出事应力项后，则表述的是一个应力功的差值， 而不是真正存在的物理功。应力循环功理论和初始应力项　有着密切的关系，前述推导也大量应用了的任意性来论述屈服面的特征。

对于第二个假设，即存在唯一塑性屈服面。这个假设存在一个由全量形式的 主应力表示的屈服面Fσ。Fσ只与主应力大小有关，而与应力增量无关的性质，目前正受到质疑。在一些岩土实验中，确实观测到了Fσ与应力增量有关的情况。所 以这个假设在某些情况下并不成立。

对于第三个假设，即塑性势面和屈服重合。如果塑性势面和屈服面并不重合。塑性势面的法线方向为塑性应变增量的矢量方向。假定塑性应变增量的矢量方向在重合应力空间中为X正方向，则塑性势能面垂直于X轴。则应力屈服面上任一 点　ijσ可做一条垂直于X轴的塑性势能面，且有一个应力屈服面与塑性势能面相 交。应力屈服面内侧区域被塑性势能面分割为两个部分，其中屈服面内侧位于塑性势能面左侧的区域均满足应力循环功理论，而屈服面内侧位于塑性势能面右侧的区域则与应力循环功理论相左。所以应力循环功理论又假定了塑性势能面和屈服面必须重合，否则在应力屈服面内只有同时满足塑性势能面条件的区域才满足应力循环功理论。

进一步讨论应变循环功理论。应变循环功理论的一般表述为。根据应力循环功理论的分析，在应变空间中，同样可以对应变屈服面进行分析。根据应变循环功理论，在应变空间中考虑问题，任意一点的应变状态可以用矢量表示，同时任意一点的应力状态也可以用矢量表示。如果将应变空间和塑性残余应力空间重合，在这个重叠坐标空间中进行分析，可以得到塑性残余应力增量和应变屈服面之间的关系。为了分析塑性残余应力增量方向与应变空间中屈服面的关系，将塑性残余应力的原点放在屈服面上，同时让塑性残余应力的坐标和应变空间的应变坐标重合。假定初始应变在应变空间中的 坐标为，其在应变空间中的矢量表示为。当应变加载至应变空间屈服面时，应变在应变空间中的坐标为εij，此时应变在应变空间中的坐标为。则应变的增量，在应变空间中采用矢量可以表示为。根据应变循环功理论。由于我们采用重叠空间对应变和塑 性残余应力进行分析，可以得到如下形式的应变循环功理论表示为重叠空间中塑性残余应力增量的矢量表示，根据内积定律改写应变循环功理论为，其中θ为重叠空间中 应变增量和塑性残余应力增量的夹角。由上述该关系式可知，在重叠空间中应变增量和塑性残余应力增量的夹角小于等于90度。由于该结果对于应变空间屈服面 以内的任意初始应变　均成立。在重叠空间屈服面上任意取一点ijε，做切平面T （εij ），同时取重叠空间中塑性残余应力增量矢量为，则以εij为原点，在切平面T （εij）周围做无限小圆周γ上取到　，可知在此圆周上均有。在γ上，如 果，存在一点　令，则取公共面P内　的反向延长线上一点，则与应变循环功理论相矛盾，因此在γ上，均有，即塑性残余应力增量矢量垂直于重叠空间中的应变屈服面。根据应变循环功理论，在重叠空间中任意初始应变与应变屈服面上εij点之间的向量与εij点之处的塑性残余应力增量之 间的夹角θ均应当不大于90度。同时根据上述推导，可知在重叠空间中塑性应力 增量矢量平行于重叠空间中的应变屈服面上　εij点之处法向矢量。所以与εij点之处法向矢量夹角应当不大于90度。所以所有初始应变点均应在屈服面Fε以内，并且也在屈服面上点εij的切平面T （εij）内侧。这说明应变 循环功理论要求，重叠空间中应变屈服面的几何形状必须在屈服面上任意点切平面的内侧。

（4）SEFE方法弹塑性本构方程

前面几小节对SEFE建模理论中所使用的线弹性假设和塑性假设进行了论述，但仍然没有给出SEFE建模方法使用的弹塑性本构方程的具体形式。本节将继续讨论SEFE建模理论，并推导分析SEFE建模方法所使用的弹塑性本构方程的增量形式。

仿照线弹性本构方程的形式，可以写出弹塑性增量本构方程：

其中[]σ和[]ε为应力张量　ijσ和应变张量ijε的矢量形式：

弹塑性增量本构方程中的　epD为弹塑性模量矩阵。基于上一节推导的塑性流动法则，可以将　epD写作：

其中，E和G为杨氏模量和剪切模量，ν为泊松比。和 H为有效应力和塑性模量。一般H需要从应力-应变测试数据获取。

4.3.1.4 储层模拟的初始应变

储层模拟的初始应变可以表述为如下形式：

4.3.1.5虚功原理

模型边界条件采用准静态加载方式，将岩层应变控制方程，写成虚功原理的形式：

其中w为总能量，对于线弹性材料w写作：

其中　μσ为应力张量，με为应变张量，wp为地层孔隙水压力，力，是节点位移和节点作用力。 为边界线应

4.3.1.6 增量拉格朗日断层接触解法

图53中Master和Slave分别代表断层的两侧，g为两侧间距。

图53 断层接触算法图解

利用罚函数　np，定义断层两侧压力npT：

其中　np定义为非常大的数字。由该公式可知，如果一侧断层切入另一侧断层，则两盘压力Tnp会很快变为负值。这样为了保持平衡，断层两侧的间隙距离会变小， 直到平衡为止。

断层两盘切向摩擦Ttp：

这里Tttrial的计算公式为：

式中，pt是罚因子；和xm为断层两盘的相对位移；为xm上步的数值，并且　为当前计算步的滑移矢量。

Tterit的计算公式为：

其中μ是摩擦系数；cohe是内聚滑动阻力；　Tmaxt是最大摩擦牵引力。

摩擦模型采用库仑摩擦模型，摩擦系数μ计算公式为：

μd为动态的摩擦系数；vs为滑移速度；　μo为静态摩擦系数；β为衰变系数。