几何模型、网格与边界条件设计介绍

图6-6和图6-7是模型示意图。

从图6-6看出，洞穴的平面形状类似于足球场，两端是半圆形的，圆心处为两个注水的井口，热水从这里注入地下，从而溶化地下盐岩层，造出空腔，用于储气。图6-7的柱形图给出了沿深度不同地层的材料变化情况。由于结构有XOZ和YOZ两个对称面，计算中可以只取1/4结构进行离散。结构整体网格如图6-8所示，总共采用了13916个节点、11828

图6-7 岩层分布深度柱形图

图6-8 模型网格图

个8节点长方体单元对结构进行离散。整个结构从上到下全长1188m，含有6个不同的岩层构造，分别为：上部弹性层（1025m厚），凝灰岩层（30m厚），页岩层（5m厚），钾盐层（40m厚），含磷岩层（12m厚），钾盐层（76m厚）。结构的长和宽分别为140m和105m。储藏洞室选在1100m深度处的蠕变性能良好的含磷岩层。图6-9给出了施工结束时储藏洞室的结构形状示意图，两端呈半圆形，中间呈矩形，洞室内部净高为12m。洞室中间井孔轴线的下部有一个初始融化产生的深10m、直径为10m的凹坑。为了清楚起见，图6-9中隐去了上部的弹性层，只显示下部的岩层构造。

图6-9 模型内部形状示意图

模型中用于输水的套管尺寸为外径为0.194m、内径为0.172m，套管与岩石之间有摩擦接触的相互作用。上水管和注水管相距80m。图6-9中只给出了注水管的位置，水管在被对称性简化去掉的那一部分上。在溶化开采施工时，围岩融腔表面的液体压力为注入的高压水压力，其值为13MPa；施工结束时这个压力变为储存的液体天然气的压力。由于施工表面上始终作用着常量液压，这个问题的数值模拟相当于一个计算力学中的移动边界问题，数值计算时与一般的洞室开挖如隧道施工等相比有较高的难度。本例结合ABAQUS的∗Internal Surface模型来施加这个水压力。这个水压力的作用有两个：一是溶解岩层，二是支护岩层。计算时在开挖步施加之前，提前在预计开挖到的表面上施加了水压力。边界条件分别为：对4个侧面施加了各面法向的单向位移约束，下底面施加了z向约束，上表面自由。施工的具体模拟步骤为：①施加地应力；②移出套管和水泥环单元形成表面自由的井孔，然后再植入套管和水泥环单元，并在套管单元内表面施加液体压力；③在给定的施工面（内腔面）上施加液体载荷，移出融化开挖掉的单元；④重复执行步③，模拟整个洞穴的溶化开采过程。为了保证计算稳定，整个融化开挖采用了20个开挖步模拟融化开挖过程。

1）施加地应力。初始地应力的构造采用了结合实测结果（施工层竖向地应力分量实测为25MPa）、预加沿竖向分段线性的z向应力分量σ_z、分层取侧应力系数，然后施加z向体积力的方法，采用初始条件中的地应力项构造平衡的初始地应力场。岩盐层的侧压力系数取0.99，其他岩层的侧压力系数取0.85。这样构造的初始应力场比仅采用重力法构造的应力场所对应的位移场和应变场要好得多。在inp文件中的语句为：

∗∗

∗∗在分析步的前面首先施加初始地应力：

∗INITIAL CONDITIONS，TYPE=STRESS，GEOSTATIC

part-1-1.SET-HALITE-3，-27436315.2，0，-25834995.2，76，0.99，0.99，

part-1-1.SET-SYLVINITE-2，-25834995.2，76，-25590387.2，88，0.99，0.99，

part-1-1.SET-HALITE-2，-25590387.2，88，-25442897.2，95，0.99，0.99，

part-1-1.SET-SYLVINITE，-25442897.2，95，-25381745.2，98，0.99，0.99，

part-1-1.SET-HALITE，-25381745.2，98，-24728575.2，129，0.99，0.99，

part-1-1.SET-SHALE，-24728575.2，129，-24602155.2，134，0.85，0.85，

part-1-1.SET-LIMESTONE，-24602155.2，134，-23843635.2，164，0.85，0.85，

part-1-1.SET-CLASTICS，-23843635.2，164，0，1188，0.85，0.85，

∗∗

在第一个分析步中施加相应的体积力、模拟重力，以平衡初始地应力，语句为：

∗∗Name：GRAVITY-1 Type：Gravity

∗Dload

part-1-1.SET-CLASTICS，BZ，-23284.8

part-1-1.SET-LIMESTONE，BZ，-25284

part-1-1.SET-SHALE，BZ，-25284

part-1-1.SET-HALITE，BZ，-21070

part-1-1.SET-SYLVINITE，BZ，-20384

part-1-1.SET-HALITE-2，BZ，-21070

part-1-1.SET-SYLVINITE-2，BZ，-20384

part-1-1.SET-HALITE-3，BZ，-21070

在给定的位移边界约束下，模型中承受初始地应力的岩层的初始位移很小，小于1mm，对计算结果影响甚微，本例中忽略不计。

2）井孔模拟的相关操作语句。涉及井孔的操作包括：①移出套管和水泥环单元；②形成具有自由表面的井孔后再植入套管和水泥环单元；③在套管单元内部表面施加液体压力。详细如下：

∗∗

∗∗STEP：REMOVE

∗∗

∗STEP，NAME=REMOVE

∗STATIC

1.，1.，1E-05，1.

∗MODEL CHANGE，REMOVE

PIPE，CONCRETE

∗END STEP

∗∗

这里把“移出井孔处的套管和水泥环单元”作为一个单独的分析步是为了模拟钻孔过程。钻孔之后，围岩收缩，由此形成的井壁内表面没有内压力。在此之后，实际工程操作是将没有初始应变的套管自由地植入直径稍大的井孔，然后浇注混凝土固井。因此，下面的套管和混凝土保护环相关单元植入时设为strain free状态。

∗∗STEP：Step-add

∗∗

∗STEP，NAME=STEP-ADD

∗STATIC

0.01，1.，1E-05，1.

∗∗

∗∗LOADS(www.xing528.com)

∗∗

∗MODEL CHANGE，ADD=STRAIN FREE

PIPE，CONCRETE

∗∗

∗∗LOADS

∗∗

∗∗Name：load-pressure type：pressure

∗DSLOAD，OP=NEW

PICKEDSURF279，P，1.328E+07

∗∗NAME：SURFFORCE-1 TYPE：PRESSURE

∗DSLOAD，OP=NEW

PART-1-1.PIPE-SURF，HP，1.328E+07，1188.，88.

∗∗Name：gravity-1 type：gravity

∗DLOAD

CONCRETE，BZ，-23520.

∗∗

上述过程中，在激活单元的同时，给套管和混凝土单元施加了重力，最后施加了套管内表面液体压力。

3）在给定的溶解施工面（内腔表面）上施加液体载荷。这个操作比较复杂。首先要使用internal surface命令定义建立在模型内部的所谓内部表面。下面给出了最初的溶解开采初始洞穴Cavern-Center时的内表面操作命令：

∗Elset，elset=CAVERN-CENTER-SURFSS1，internal，instance=PART-1-1

428，429，430，434，435，436，437，438，439，440，441，442，443，444，445，446

447，1494，1495，1496，1497，1498，1499，1500，1501，1502，1503，1504，1505，1506，1507，1508

1509，1510，1511，1512，1513，1514，1515，1516，1517，1518，1519，1520，1521，1522，1523，1524

1525，1526，1527，1528，1529，1530，1531，1532，1533，1534，1535，1536，1537，1538，1539，1540

1541，1542，1543，1544，1545

∗Surface，type=ELEMENT，name=CAVERN-CENTER-SURFS

CAVERN-CENTER-SURFSS1，S1

CAVERN-CENTER-SURFSS2，S2

CAVERN-CENTER-SURFSS3，S3

CAVERN-CENTER-SURFSS5，S5

然后才能在载荷步中向内部表面施加分布面载荷，即溶解开采的液体产生的压力：

∗∗

∗∗LOADS

∗∗

∗∗Name：SURFFORCE-2 Type：Pressure

∗Dsload，op=NEW

CAVERN-CENTER-SURFS，P，1.328e+07

∗∗

这时所在的分析步就要考虑岩石的蠕变效应。相应的分析步定义为：

∗∗

∗Step，name=cavern-mining，amplitude=RAMP，inc=1000，nlgeom=YES，un-symm=YES

∗Visco，cetol=0.1，CREEP=explicit

0.01，100.，1e-15，100.

∗∗

∗MODEL CHANGE，REMOVE

CAVERN-center

∗∗

值得注意的是，由于塑性问题的收敛性有时候不能保证，这时就需要采用多个分析步，使用∗restart语句进行接续运算。因此，结果输出的操作中加入相关的数据操作是必要的：

∗∗OUTPUT REQUESTS

∗∗

∗Restart，write，frequency=1

∗∗

如果上述的frequency=0，那么将不会输出中间结果的二进制文件，计算时所用的硬盘空间会较节省，但是如果一个分析步迭代不收敛的话就得改进数据后重新从0开始再次计算。对于多分析步的岩土工程弹塑性分析而言，预先输出中间结果以备重启动分析使用的做法将节省更多的计算成本。