教程(二）COMSOL中实现流固耦合理论介绍

一般固体变形控制方程主要由三个方程构成：应力平衡方程、几何变形方程、本构方程，一般以Navier的形式出现,此公式推导涉及到弹性力学的基本理论，在此不再推导。式(1)中εv、u、Fx分别为体应变、位移、x方向体载荷。

化成张量形式(一种表述方法)

如果出现流固耦合或者温度、煤基质变形引起应力耦合，则需要添加额外项。此过程以太沙基提出的有效应力原理为基础，如式(3)。式中α为有效应力系数即Biot系数，p为孔隙压力。太沙基的有效应力方程是针对单孔隙提出的，而对于像煤层这些双重孔隙/裂隙介质的多孔介质而言，需要作出一些修正，如式(4)。式(4)中考虑了基质中孔压与裂隙中孔压对有效应力的影响。对于流固耦合问题，便是讨论有效应力下的变形控制方程，这样便考虑到孔压对固体变形的影响。将式(3)带入到式(2)得到，得到考虑流固耦合的张量形式，如式(5)。

式(5)考虑了孔压对有效应力影响，还可以考虑其他应力对有效应力影响如温度引起的热应力、煤体基质变形引起的应力等。对于多孔介质中流体的流动方程，一般采用达西流动，非饱和流动的理查兹方程，其中达西流动较为简单，一般适用于低速线性流动，如式(6)。固体中的渗透率一般与应力或者应变有关系，此时固体变形将会通过影响孔隙率和渗透率，进而影响流体的流动，流体的流动又导致孔压发生变化，影响固体的有效应力，达到流体和固体之间的双向耦合。

COMSOL中如何实现流固耦合?按照前文推导的公式，选用“固体力学”模块与“达西定律”模块。固体力学模块中线弹性材料中的控制方程便是式(2)，还需要添加一项代表孔压的影响。从式(5)分析可以看到，把孔压项当做体载荷，输入到COMSOL中。Fi为重力引起的体载荷，在需要考虑重力项时，可以把重力项加入到体载荷中，不需要考虑时，即可忽略Fi此项。图1为体积力设置项，选择体载荷。图2是体载荷设置，选择“单位体积的力”，在x,y栏分别输入alpha1*dl.px1与alpha1*dl.px2。dl.px1、dl.px2表示压力在x、y方向的梯度，即压力p对x或y求偏导。设置好体载荷后，然后设置边界载荷和边界条件，这样固体变形控制方程就在COMSOL设置好了。

图1 COMSOL中体载荷与重力栏

图2 体载荷设置

对于达西定律，此物理场设置较为简单。按照流体和基本属性的栏顺序，依次输入。边界设置边界压力，同时设置初始压力。以上设置完成后，选择瞬态求解，把固体力学与达西定律均选上，设置求解时间即可求解。流固耦合问题的难点在于对控制方程的理解以及在COMSOL中的输入，对于其他耦合问题，可以参考流固耦合的设置。