井周应力计算方法及计算装置

摘要

本发明公开了一种井周应力计算方法及计算装置，其中，方法包括：建立坐标系；分别选取物理场并确定相互间的耦合关系；建立同心空心圆柱；定义材料物性；选择线性单元划分网格；施加边界条件；在地应力坐标系中，定义初始物理场；选择瞬态分析类型，以得到输出结果；在井筒坐标系中更新输出结果，输出沿径向变化的地层压力和有效应力分量。该方法能在同时考虑任意地应力、各向异性地层和多物理场耦合的情况下，精确高效地计算斜井的井周应力，可用于井壁稳定性实时分析中的应力计算。

权利要求

1.一种井周应力计算方法，其特征在于，包括以下步骤：

建立井筒坐标系、材料坐标系和地应力坐标系，其中，所述井筒坐标系的轴线与井筒方向同向，所述材料坐标系的坐标轴方向与材料主轴方向相同，所述地应力坐标系的坐标轴方向与地应力方向相同，其中，所述地应力坐标系为总体坐标系；

分别选取物理场并确定相互间的耦合关系，其中，根据所关心的物理过程选择相应的物理场，并在方程视图中通过修改变量值和弱形式来定义物理场之间的耦合关系；

建立同心空心圆柱，其中，所述同心空心圆柱的对称轴与所述井筒坐标系的轴线同向；

在所述材料坐标系中定义材料物性；

选择线性单元划分网格，使圆柱两个端面的面网格相同；

施加边界条件，其中，在圆柱外表面为远场边界条件，且根据钻井液状态确定内表面边界条件，使两个端面对应节点的相应自由度相等；

在所述地应力坐标系中，定义初始物理场；

选择瞬态分析类型，以得到输出结果；以及

在所述井筒坐标系中更新所述输出结果，输出沿径向变化的地层压力和有效应力分量。

2.根据权利要求1所述的井周应力计算方法，其特征在于，在建立所述同心空心圆柱时，圆柱代表地层，内孔代表井眼，内孔直径表示井眼直径，圆柱外径不小于五倍内孔直径，以使圆柱外表面可代表远场。

3.根据权利要求1所述的井周应力计算方法，其特征在于，所述选择线性单元划分网格进一步包括：在所述空心圆柱的一个端面划分四边形面网格，并且由此端面的面网格扫略到另一端面，以使两个端面上的面网格完全相同，扫略过程中在圆柱高度方向只划分单层网格，其中，从内孔轴线向外采用由密到疏的网格过渡方法。

4.根据权利要求1所述的井周应力计算方法，其特征在于，在定义所述初始物理场时，地应力在Biot本构中为总应力。

5.一种井周应力计算装置，其特征在于，包括：

第一建立模块，用于建立井筒坐标系、材料坐标系和地应力坐标系，其中，所述井筒坐标系的轴线与井筒方向同向，所述材料坐标系的坐标轴方向与材料主轴方向相同，所述地应力坐标系的坐标轴方向与地应力方向相同，其中，所述地应力坐标系为总体坐标系；

选取模块，用于分别选取物理场并确定相互间的耦合关系，其中，根据所关心的物理过程选择相应的物理场，并在方程视图中通过修改变量值和弱形式来定义物理场之间的耦合关系；

第二建立模块，用于建立同心空心圆柱，其中，所述同心空心圆柱的对称轴与所述井筒坐标系的轴线同向；

第一定义模块，用于在所述材料坐标系中定义材料物性；

选择模块，用于选择线性单元划分网格，使圆柱两个端面的面网格相同；

施加模块，用于施加边界条件，其中，在圆柱外表面为远场边界条件，且根据钻井液状态确定内表面边界条件，使两个端面对应节点的相应自由度相等；

第二定义模块，用于在所述地应力坐标系中，定义初始物理场；

获取模块，用于选择瞬态分析类型，以得到输出结果；以及

输出模块，用于在所述井筒坐标系中更新所述输出结果，输出沿径向变化的地层压力和有效应力分量。

6.根据权利要求5所述的井周应力计算装置，其特征在于，在建立所述同心空心圆柱时，圆柱代表地层，内孔代表井眼，内孔直径表示井眼直径，圆柱外径不小于五倍内孔直径，以使圆柱外表面可代表远场。

7.根据权利要求5所述的井周应力计算装置，其特征在于，所述选择模块还用于在所述空心圆柱的一个端面划分四边形面网格，并且由此端面的面网格扫略到另一端面，以使两个端面上的面网格完全相同，扫略过程中在圆柱高度方向只划分单层网格，其中，从内孔轴线向外采用由密到疏的网格过渡方法。

8.根据权利要求5所述的井周应力计算装置，其特征在于，在定义所述初始物理场时，地应力在Biot本构中为总应力。

说明书

井周应力计算方法及计算装置

技术领域

[0001] 本发明涉及井壁稳定性分析技术领域，特别涉及一种井周应力计算方法及计算装置。

背景技术

[0002] 井壁稳定性是影响钻井效率和成本的关键因素之一。井壁稳定性破坏的两种主要形式是，钻井液密度过大时导致井壁在拉应力作用下破裂，或者钻井液密度过小时导致井
壁在压应力作用下坍塌，因此井壁稳定性分析的核心任务之一是在确定的钻井液密度下计
算井周应力。随着页岩油气等非常规资源开采力度的加大和水平钻井技术的应用日趋广
泛，在井壁稳定性分析中需考虑的因素也越来越多，包括井筒方向的变化、地应力的各向异
性、地层力学性能的各向异性和多物理场耦合(应力、渗流、化学和温度等)等。因此，在井周
应力计算模型中考虑多因素的影响是目前井壁稳定性分析的工程需求和发展趋势。

[0003] 相关技术中，井周应力计算的方法主要有解析法和有限元法。解析法的优势是计算效率高，但大多数实际问题不存在解析解，因此，解析法的应用范围非常有限，一般只在
特定问题中用于验证数值模型。在有限元法方面，二维有限元模型采用平面应变或广义平
面应变假设，不能考虑实际发生的面内翘曲变形，只适用于对部分问题做定性分析；三维有
限元模型为减小地应力在井壁周围产生的边界效应，模型规模一般要做的很大，因此自由
度多，求解时间长，无法满足实时化井壁稳定性分析的要求。井壁稳定性分析专业软件
PBORE-3D等采用伪三维有限元模型计算井周应力，有较高的计算效率和精度，但受软件既
定功能的限制，在多场耦合问题的模拟中缺乏灵活性，而且涉及知识产权问题，性价比较
低。此外，一般的通用商业有限元软件中不提供多物理场耦合的伪三维单元，虽然理论上用
户可通过自己编写单元来解决这一问题，但自定义单元功能要受商业软件中(如在ABAQUS
中)可用自由度等诸多限制，往往难以实现，并且需要较高的专业技能，一般研究者难以胜
任。

发明内容

[0004] 本发明旨在至少在一定程度上解决相关技术中的技术问题之一。

[0005] 为此，本发明的一个目的在于提出一种井周应力计算方法，该方法可以精确高效地计算井周应力。

[0006] 本发明的另一个目的在于提出一种井周应力计算装置。

[0007] 为达到上述目的，本发明一方面实施例提出了一种井周应力计算方法，其包括以下步骤：建立井筒坐标系、材料坐标系和地应力坐标系；分别选取物理场并确定相互间的耦
合关系；建立同心空心圆柱，其中，所述同心空心圆柱的对称轴与所述井筒坐标系的轴线同
向；在所述材料坐标系中定义材料物性；选择线性单元划分网格，使圆柱两个端面的面网格
相同；施加边界条件，其中，在圆柱外表面为远场边界条件，且根据钻井液状态确定内表面
边界条件，使两个端面对应节点的相应自由度相等；在所述地应力坐标系中，定义初始物理
场；选择瞬态分析类型，以得到输出结果；在所述井筒坐标系中更新所述输出结果，输出沿
径向变化的地层压力和有效应力分量。

[0008] 本发明实施例的井周应力计算方法，可以在同时考虑井筒方向变化、地应力各向异性、地层力学性能各向异性和多物理场耦合四个因素的情况下，使用通用有限元平台
COMSOL Multiphysics计算井周应力，且模型中对边界条件的描述准确，并因只采用单层六
面体单元而具有自由度少的特征，模型求解时间短，计算结果精度高，满足对井壁稳定性做
实时分析的要求，并且基于COMSOL Multiphysics的通用功能定义多物理场耦合关系，便于
灵活地定义和修改耦合关系，以及基于通用有限元平台COMSOL Multiphysics的常用功能，
操作简单，并可方便地实现参数化建模，便于工程技术人员掌握和使用。

[0009] 另外，根据本发明上述实施例的井周应力计算方法，还可以具有以下附加的技术特征：

[0010] 进一步地，在本发明的一个实施例中，所述井筒坐标系的轴线与井筒方向同向，所述材料坐标系的坐标轴方向与材料主轴方向相同，所述地应力坐标系的坐标轴方向与地应
力方向相同，其中，所述地应力坐标系为总体坐标系。

[0011] 进一步地，在本发明的一个实施例中，在建立所述同心空心圆柱时，圆柱代表地层，内孔代表井眼，内孔直径表示井眼直径，圆柱外径大于且等于五倍内孔直径，以使圆柱
外表面可代表远场。

[0012] 进一步地，在本发明的一个实施例中，所述选择线性单元划分网格进一步包括：在所述空心圆柱的一个端面划分四边形面网格，并且由此端面的面网格扫略到另一端面，以
使两个端面上的面网格完全相同，扫略过程中在圆柱高度方向只划分单层网格，其中，从内
孔轴线向外采用由密到疏的网格过渡方法。

[0013] 进一步地，在本发明的一个实施例中，在定义所述初始物理场时，地应力在Biot本构中为总应力。

[0014] 为达到上述目的，本发明另一方面实施例提出了一种井周应力计算装置，包括：第一建立模块，用于建立井筒坐标系、材料坐标系和地应力坐标系；选取模块，用于分别选取
物理场并确定相互间的耦合关系；第二建立模块，用于建立同心空心圆柱，其中，所述同心
空心圆柱的对称轴与所述井筒坐标系的轴线同向；第一定义模块，用于在所述材料坐标系
中定义材料物性；选择模块，用于选择线性单元划分网格，使圆柱两个端面的面网格相同；
施加模块，用于施加边界条件，其中，在圆柱外表面为远场边界条件，且根据钻井液状态确
定内表面边界条件，使两个端面对应节点的相应自由度相等；第二定义模块，用于在所述地
应力坐标系中，定义初始物理场；获取模块，用于选择瞬态分析类型，以得到输出结果；输出
模块，用于在所述井筒坐标系中更新所述输出结果，输出沿径向变化的地层压力和有效应
力分量。

[0015] 本发明实施例的井周应力计算装置，可以在同时考虑井筒方向变化、地应力各向异性、地层力学性能各向异性和多物理场耦合四个因素的情况下，使用通用有限元平台
COMSOL Multiphysics计算井周应力，且模型中对边界条件的描述准确，并因只采用单层六
面体单元而具有自由度少的特征，模型求解时间短，计算结果精度高，满足对井壁稳定性做
实时分析的要求，并且基于COMSOL Multiphysics的通用功能定义多物理场耦合关系，便于
灵活地定义和修改耦合关系，以及基于通用有限元平台COMSOL Multiphysics的常用功能，
操作简单，并可方便地实现参数化建模，便于工程技术人员掌握和使用。

[0016] 另外，根据本发明上述实施例的井周应力计算装置还可以具有以下附加的技术特征：

[0017] 进一步地，在本发明的一个实施例中，所述井筒坐标系的轴线与井筒方向同向，所述材料坐标系的坐标轴方向与材料主轴方向相同，所述地应力坐标系的坐标轴方向与地应
力方向相同，其中，所述地应力坐标系为总体坐标系。

[0018] 进一步地，在本发明的一个实施例中，在建立所述同心空心圆柱时，圆柱代表地层，内孔代表井眼，内孔直径表示井眼直径，圆柱外径大于且等于五倍内孔直径，以使圆柱
外表面可代表远场。

[0019] 进一步地，在本发明的一个实施例中，所述选择模块还用于在所述空心圆柱的一个端面划分四边形面网格，并且由此端面的面网格扫略到另一端面，以使两个端面上的面
网格完全相同，扫略过程中在圆柱高度方向只划分单层网格，其中，从内孔轴线向外采用由
密到疏的网格过渡方法。

[0020] 进一步地，在本发明的一个实施例中，在定义所述初始物理场时，地应力在Biot本构中为总应力。

[0021] 本发明附加的方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。

附图说明

[0022] 本发明上述的和/或附加的方面和优点从下面结合附图对实施例的描述中将变得明显和容易理解，其中：

[0023] 图1为根据本发明实施例的井周应力计算方法的流程图；

[0024] 图2为根据本发明一个实施例的井周应力计算方法的流程图；

[0025] 图3为根据本发明一个实施例的坐标系示意图；

[0026] 图4为根据本发明一个实施例的几何模型示意图；

[0027] 图5为根据本发明另一个实施例的几何模型示意图；

[0028] 图6为根据本发明一个实施例的模型的有限元网格示意图；

[0029] 图7为根据本发明另一个实施例的模型的有限元网格示意图；

[0030] 图8为根据本发明一个实施例的沿径向变化的地层压力示意图；

[0031] 图9为根据本发明一个实施例的沿径向变化的有效径向应力示意图；

[0032] 图10为根据本发明一个实施例的沿径向变化的有效环向应力示意图；

[0033] 图11为根据本发明一个实施例的杨氏模量变化时沿径向变化的有效轴向应力示意图；

[0034] 图12为根据本发明一个实施例的泊松比变化时沿径向变化的有效轴向应力示意图；

[0035] 图13为根据本发明实施例的井周应力计算装置的结构示意图。

具体实施方式

[0036] 下面详细描述本发明的实施例，所述实施例的示例在附图中示出，其中自始至终相同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附
图描述的实施例是示例性的，旨在用于解释本发明，而不能理解为对本发明的限制。

[0037] 下面参照附图描述根据本发明实施例提出的井周应力计算方法及计算装置，首先将参照附图描述根据本发明实施例提出的井周应力计算方法。

[0038] 图1是本发明实施例的井周应力计算方法的流程图。

[0039] 如图1所示，该井周应力计算方法包括以下步骤：

[0040] 在步骤S101中，建立井筒坐标系、材料坐标系和地应力坐标系。

[0041] 其中，在本发明的一个实施例中，井筒坐标系的轴线与井筒方向同向，材料坐标系的坐标轴方向与材料主轴方向相同，地应力坐标系的坐标轴方向与地应力方向相同，其中，
地应力坐标系为总体坐标系。

[0042] 具体地，如图2所示，步骤S1，建立坐标系：建立圆柱坐标系cy1(井筒坐标系)，其轴线为井筒方向；建立直角坐标系cy2(材料坐标系)，坐标轴方向与材料主轴方向相同；建立
直角坐标系cy3(地应力坐标系)，坐标轴方向与地应力方向相同；将cy3作为总体坐标系。

[0043] 例如，如图3所示，图3是地应力坐标系Oxyz(总体坐标系)和井筒坐标系Ox'y'z'的位置关系(Sh,max是水平最大地应力，Sh,min是水平最小地应力，Sv是垂直地应力，φx为Ox与
Ox’在Oxy平面上投影之间的夹角，φz为Oz与Oz’之间的夹角，φx和φz完全定义了两个坐标
系之间的位置关系)。

[0044] 在步骤S102中，分别选取物理场并确定相互间的耦合关系；

[0045] 进一步地，如图2所示，步骤S2，选取物理场并定义相互间的耦合关系：根据所关心的物理过程选择相应的物理场，并在方程视图中通过修改变量值和弱形式来定义物理场之
间的耦合关系。

[0046] 在步骤S103中，建立同心空心圆柱，其中，同心空心圆柱的对称轴与井筒坐标系的轴线同向。

[0047] 进一步地，在本发明的一个实施例中，在建立同心空心圆柱时，圆柱代表地层，内孔代表井眼，内孔直径表示井眼直径，圆柱外径大于且等于五倍内孔直径，以使圆柱外表面
可代表远场。

[0048] 进一步地，如图2所示，步骤S3，几何建模：建立同心空心圆柱，其对称轴为cy1的轴线坐标轴；圆柱代表地层，内孔代表井眼，内孔直径d即井眼直径；圆柱外径D≥5d，目的是使
圆柱外表面可代表远场，具体值可以综合考虑计算精度和时间确定；圆柱高度h的选择应优
先保证内孔周围的单元质量。

[0049] 例如，如图4和图5所示，图4是几何模型，图5是图4中框图部分放大后，其中，cy1为井筒坐标系，cy3为地应力坐标系，cy3的一个坐标轴因指向读者而未画出，φx＝30°，φz＝
60°。

[0050] 在步骤S104中，在材料坐标系中定义材料物性。

[0051] 需要说明的是，如图2所示，步骤S4，定义材料参数：在cy2坐标系中定义材料物性，此时材料参数的表达式最为简洁；需要注意的是，在渗流-应力耦合分析中定义Biot本构关
系时，Biot有效应力张量αij在cy2中只有i＝j的元素是非零的。

[0052] 在步骤S105中，选择线性单元划分网格，使圆柱两个端面的面网格相同。

[0053] 进一步地，在本发明的一个实施例中，选择线性单元划分网格进一步包括：在空心圆柱的一个端面划分四边形面网格，并且由此端面的面网格扫略到另一端面，以使两个端
面上的面网格完全相同，扫略过程中在圆柱高度方向只划分单层网格，其中，从内孔轴线向
外采用由密到疏的网格过渡方法。

[0054] 进一步地，如图2所示，步骤S5，划分网格：选择线性单元；首先在空心圆柱的一个端面划分四边形面网格，再由此端面的面网格扫略到另一端面，保证两个端面上的面网格
完全相同，扫略过程中在圆柱高度方向只划分单层网格；考虑到内孔附近有应力集中，从内
孔轴线向外采用由密到疏的网格过渡方法，以同时保证计算精度和效率。例如，如图6和图7
所示所示，图6是模型的有限元网格，图7是图6中框图部分放大后。

[0055] 在步骤S106中，施加边界条件，其中，在圆柱外表面为远场边界条件，且根据钻井液状态确定内表面边界条件，使两个端面对应节点的相应自由度相等。

[0056] 进一步地，如图2所示，步骤S6，施加边界条件：在圆柱外表面上取各自由度值为远场值(如位移为零)；根据钻井液的状态定义空心圆柱内表面上的边界条件，如使其内表面
上的地层压力和温度等分别等于钻井液的相应值；使空心圆柱两个端面对应节点的相应自
由度相等，所谓对应节点即其中一个节点是由另一节点扫略得到的。

[0057] 在步骤S107中，在地应力坐标系中，定义初始物理场。

[0058] 其中，在本发明的一个实施例中，在定义初始物理场时，地应力在Biot本构中为总应力。

[0059] 进一步地，如图2所示，步骤S7，定义初始物理场：在地应力坐标系cy3中，以初始应力 的形式在整个模型上施加地应力 在整个模型上施加初始地层压力p0、温度和化
学势等；在多物理场耦合问题中需注意的是，地应力在Biot本构关系(以及基于Biot本构推
广得到的其他本构关系)中为总应力 因此，定义初始应力场和初始地层压力场等之后
应保证总应力等于地应力，即 和 在cy3中地应力 只有i＝j的元素
非零，但αij在cy3中的表达式需要通过由cy2中的表达式做坐标转换得到。

[0060] 在步骤S108中，选择瞬态分析类型，以得到输出结果。

[0061] 进一步地，如图2所示，步骤S8，提交分析任务：选择瞬态分析类型，为保证计算精度初始时间增量步宜小，随着时间增长增量步可逐渐增大。

[0062] 在步骤S109中，在井筒坐标系中更新输出结果，输出沿径向变化的地层压力和有效应力分量。

[0063] 进一步地，如图2所示，步骤S9，后处理：在井筒坐标系中计算结果更容易解读，因此，在cy1坐标系中更新计算结果。

[0064] 另外，如图8-12所示，图8是沿径向变化的地层压力，其中，r是距井筒中心的径向距离，R是井筒半径，nν＝ν/ν′，ν是横观各向同性面内的泊松比，ν′是过对称轴平面内的泊
松比，t是时间，d是“天”)；图9是沿径向变化的有效径向应力(图例同图8)；图10是沿径向变
化的有效环向应力(图例同图8)；图11是杨氏模量变化时沿径向变化的有效轴向应力其中，
nE＝E/E′，E为横观各向同性面内的杨氏模量，E′为垂直于横观各向同性面方向的杨氏模
量，其他图例同图8)；图12是泊松比变化时沿径向变化的有效轴向应力(图例同图8)。

[0065] 具体而言，如图2所示，本发明实施例的计算方法包括：S1建立井筒、材料和地应力坐标系；S2选取物理场；S3建立同心空心圆柱几何模型；S4在材料坐标系中定义材料物性；
S5划分网格，保证圆柱两个端面的网格相同；S6施加边界条件，使圆柱两个端面对应节点的
相应自由度相等；S7定义初始物理场，用预定义应力场的方法模拟地应力；S8提交分析；S9
在井筒坐标系中做后处理。本发明可在通用有限元平台COMSOL Multiphysics上方便地实
现，能在同时考虑任意井眼方向、各向异性的地应力和材料性能、多物理场耦合的情况下，
精确高效地计算井周应力，可用于井壁稳定性实时分析中的应力计算。

[0066] 下面以一个具体实施例对本发明实施例的计算方法的原理进行详细描述。

[0067] 在本发明的一个具体实施例中，本发明实施例的计算方法包括以下步骤：

[0068] 1.考虑渗流-应力耦合的斜井井周应力分析。

[0069] (1)建立坐标系：

[0070] 建立圆柱坐标系cy1，其轴线为井筒方向；建立直角坐标系cy2，坐标轴方向与材料主轴方向相同；将总体坐标系(直角坐标系)定义为地应力坐标系cy3； 为30°， 为60°。

[0071] (2)选择物理场和定义材料参数：

[0072] 为将计算结果与解析解对比，采用横观各向同性材料，且其对称轴为井筒轴向；在物理场中选取Poroelasticity模块，并在cy2坐标系中定义和修改材料参数poro.Dij、
poro.Slij、poro.eelvol、poro.evol、poro.S(i，j＝1-6)，建立横观各向同性材料的Biot本
构关系，同时对应地修改相关的弱形式表达式；其中在各向同性条件下，弹性模量E＝
20.6GPa，泊松比v＝0.189，Biot模量M＝15.8GPa，固相体积模量Ks＝48.2GPa，渗透率κ＝
8.64×10-6m2/MPa，通过给定nE和nν可确定E′和ν′，从而完成横观各向同性材料参数的定义。

[0073] (3)建立几何模型：

[0074] 建立同心空心圆柱，其对称轴为cy1的轴线坐标轴，内孔直径d＝0.2m，外径D＝3m，高度h＝0.01m。

[0075] (4)划分网格：

[0076] 将空心圆柱的一个端面四等分(即每份90°)，在其1/4面划分四边形结构化网格，从内孔轴线向外采用由密到疏的网格过渡方法，再将网格分别拷贝至剩下的3/4面；将已划
分网格的端面作为源面，另一端面作为目标面，划分扫略网格，注意在圆柱高度方向只划分
单层网格；选择线性单元。

[0077] (5)施加边界条件：

[0078] 在圆柱外表面位移分量ui＝0(i＝1-3)，地层压力p＝9.8MPa，空心圆柱内表面上压力和地层压力均为0MPa；令圆柱两个端面对应节点的位移分量和地层压力分别相等，即
psn＝pdn，其中sn代表源面上的某节点，dn代表在目标面上由sn经扫略
得到的节点，这一边界条件可通过在COMSOL Multiphysics中定义线性外推(Linear 
Extrusion)关系方便地实现。

[0079] (6)施加初始地层压力和地应力：

[0080] 在总体坐标系cy3中，对整个模型定义初始地层压力p0＝9.8MPa，并在整个模型上施加初始应力场，施加后的最终效果是初始总应力场等于地应力场，其中水平最大地应力
为29MPa，水平最小地应力为20MPa，垂直地应力为25MPa。

[0081] (7)提交分析任务：

[0082] 选择瞬态分析类型，初始时间步0.01s，最大时间步长100s，总时间步长为86400s，指定在0.001d(86.4s)和1d(86400s)输出计算结果。

[0083] (8)后处理：

[0084] 在cy1坐标系中更新计算结果，输出沿径向变化的地层压力和有效应力分量。

[0085] 2.计算结果的验证

[0086] 将本发明实施例的计算结果与解析解对比，本发明实施例的计算结果在所有算例中均与解析解符合良好。

[0087] 根据本发明实施例的井周应力计算方法，可以在同时考虑井筒方向变化、地应力各向异性、地层力学性能各向异性和多物理场耦合四个因素的情况下，使用通用有限元平
台COMSOL Multiphysics计算井周应力，且模型中对边界条件的描述准确，并因只采用单层
六面体单元而具有自由度少的特征，模型求解时间短，计算结果精度高，满足对井壁稳定性
做实时分析的要求，并且基于COMSOL Multiphysics的通用功能定义多物理场耦合关系，便
于灵活地定义和修改耦合关系，以及基于通用有限元平台COMSOL Multiphysics的常用功
能，操作简单，并可方便地实现参数化建模，便于工程技术人员掌握和使用。

[0088] 其次参照附图描述根据本发明实施例提出的井周应力计算装置。

[0089] 图13是本发明实施例的井周应力计算装置的结构示意图。

[0090] 如图13所示，该井周应力计算装置100包括：第一建立模块10、选取模块20、第二建立模块30、第一定义模块40、选择模块50、施加模块60、第二定义模块70、获取模块80和输出
模块90。

[0091] 第一建立模块10用于建立井筒坐标系、材料坐标系和地应力坐标系。选取模块20用于分别选取物理场并确定相互间的耦合关系。第二建立模块30用于建立同心空心圆柱，
其中，同心空心圆柱的对称轴与井筒坐标系的轴线同向。第一定义模块40用于在材料坐标
系中定义材料物性。选择模块50用于选择线性单元划分网格，使圆柱两个端面的面网格相
同。施加模块60用于施加边界条件，其中，在圆柱外表面为远场边界条件，且根据钻井液状
态确定内表面边界条件，使两个端面对应节点的相应自由度相等。第二定义模块70用于在
地应力坐标系中，定义初始物理场。获取模块80用于选择瞬态分析类型，以得到输出结果。
输出模块90用于在井筒坐标系中更新输出结果，输出沿径向变化的地层压力和有效应力分
量。本发明实施例的计算装置100能在同时考虑任意井眼方向、各向异性的地应力和材料性
能、多物理场耦合的情况下，精确高效地计算井周应力，可用于井壁稳定性实时分析中的应
力计算。

[0092] 进一步地，在本发明的一个实施例中，井筒坐标系的轴线与井筒方向同向，材料坐标系的坐标轴方向与材料主轴方向相同，地应力坐标系的坐标轴方向与地应力方向相同，
其中，地应力坐标系为总体坐标系。

[0093] 进一步地，在本发明的一个实施例中，在建立同心空心圆柱时，圆柱代表地层，内孔代表井眼，内孔直径表示井眼直径，圆柱外径大于且等于五倍内孔直径，以使圆柱外表面
可代表远场。

[0094] 进一步地，在本发明的一个实施例中，选择模块还用于在空心圆柱的一个端面划分四边形面网格，并且由此端面的面网格扫略到另一端面，以使两个端面上的面网格完全
相同，扫略过程中在圆柱高度方向只划分单层网格，其中，从内孔轴线向外采用由密到疏的
网格过渡方法。

[0095] 进一步地，在本发明的一个实施例中，在定义初始物理场时，地应力在Biot本构中为总应力。

[0096] 需要说明的是，前述对井周应力计算方法实施例的解释说明也适用于该实施例的井周应力计算装置，此处不再赘述。

[0097] 根据本发明实施例的井周应力计算装置，可以在同时考虑井筒方向变化、地应力各向异性、地层力学性能各向异性和多物理场耦合四个因素的情况下，使用通用有限元平
台COMSOL Multiphysics计算井周应力，且模型中对边界条件的描述准确，并因只采用单层
六面体单元而具有自由度少的特征，模型求解时间短，计算结果精度高，满足对井壁稳定性
做实时分析的要求，并且基于COMSOL Multiphysics的通用功能定义多物理场耦合关系，便
于灵活地定义和修改耦合关系，以及基于通用有限元平台COMSOL Multiphysics的常用功
能，操作简单，并可方便地实现参数化建模，便于工程技术人员掌握和使用。

[0098] 在本发明的描述中，需要理解的是，术语“中心”、“纵向”、“横向”、“长度”、“宽度”、“厚度”、“上”、“下”、“前”、“后”、“左”、“右”、“竖直”、“水平”、“顶”、“底”“内”、“外”、“顺时针”、“逆时针”、“轴向”、“径向”、“周向”等指示的方位或位置关系为基于附图所示的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。

[0099] 此外，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示或者
隐含地包括至少一个该特征。在本发明的描述中，“多个”的含义是至少两个，例如两个，三
个等，除非另有明确具体的限定。

[0100] 在本发明中，除非另有明确的规定和限定，术语“安装”、“相连”、“连接”、“固定”等术语应做广义理解，例如，可以是固定连接，也可以是可拆卸连接，或成一体；可以是机械连接，也可以是电连接；可以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内
部的连通或两个元件的相互作用关系，除非另有明确的限定。对于本领域的普通技术人员
而言，可以根据具体情况理解上述术语在本发明中的具体含义。

[0101] 在本发明中，除非另有明确的规定和限定，第一特征在第二特征“上”或“下”可以是第一和第二特征直接接触，或第一和第二特征通过中间媒介间接接触。而且，第一特征在
第二特征“之上”、“上方”和“上面”可是第一特征在第二特征正上方或斜上方，或仅仅表示
第一特征水平高度高于第二特征。第一特征在第二特征“之下”、“下方”和“下面”可以是第
一特征在第二特征正下方或斜下方，或仅仅表示第一特征水平高度小于第二特征。

[0102] 在本说明书的描述中，参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征、结构、材料或者特
点包含于本发明的至少一个实施例或示例中。在本说明书中，对上述术语的示意性表述不
必须针对的是相同的实施例或示例。而且，描述的具体特征、结构、材料或者特点可以在任
一个或多个实施例或示例中以合适的方式结合。此外，在不相互矛盾的情况下，本领域的技
术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结
合和组合。

[0103] 尽管上面已经示出和描述了本发明的实施例，可以理解的是，上述实施例是示例性的，不能理解为对本发明的限制，本领域的普通技术人员在本发明的范围内可以对上述
实施例进行变化、修改、替换和变型。