一种湿度运移和干缩开裂的三维数值模拟方法

1.本发明涉及土木工程、岩土工程、地质工程技术领域，尤其涉及一种湿度运移和干缩开裂的三维数值模拟方法。


背景技术：

2.土体或其他材料中的湿度运移及干缩开裂问题涉及土木工程、岩土工程、地质工程等工程应用问题。随着计算机硬件的发展及计算仿真技术的发展，采用数值模拟方法逐渐成为湿度运移和干缩开裂分析计算的有力工具。目前，求解湿度运移和干缩开裂问题的主要计算方法是有限元和离散元方法两种。其中有限元方法能够较好地模拟湿度运移，但对裂缝的起裂、扩展及交汇等难以处理。离散元方法又分为块体离散元和颗粒离散元。其中块体离散元要求对整个求解域进行完整切割，虽然可以对各种非连续问题进行求解，但难以处理块体内部的裂缝的起裂、扩展和交汇；而颗粒离散元，则因为细观输入参数无法直接试验测得的材料参数对应，需要对输入参数进行繁琐的标定。因此，上述这些方法在模拟土体或其他材料中的湿度运移及干缩开裂问题方面仍然存在很多问题。


技术实现要素：

3.有鉴于此，为解决上述问题，本发明的实施例提供了一种湿度运移和干缩开裂的三维数值模拟方法。
4.本发明的实施例提供一种湿度运移和干缩开裂的三维数值模拟方法，包括以下步骤：
5.s1通过网格划分将求解域划分为若干个多面体实体单元；
6.s2在相邻所述实体单元间可选择地设置公共单元；
7.s3根据所述实体单元内节点之间的湿度传递计算每一所述实体单元内的每一节点的水分质量流量；
8.s4将步骤s2设置的公共单元视为不存在，然后进行湿度运移计算；根据所述实体单元内每一节点的水分质量流量，获得每一节点总的水分质量流量；根据有限差分积分策略，利用每一节点总的水分质量流量更新每一节点的湿度，用每一所述实体单元中的各个节点的湿度来表示该实体单元的湿度分布，完成整个所述求解域的湿度运移的三维数值模拟；
9.s5根据所述实体单元的湿度分布计算收缩应力；
10.s6将收缩应力作为体荷载施加在所述实体单元上计算应力场；
11.s7对于设有公共单元的实体单元，根据应力场判断公共单元是否断裂；对于未设有公共单元的实体单元，根据应力场计算实体单元公共面的拉应力和剪应力，确定公共面两侧的实体单元是否要发生分离；
12.s8将断裂的公共单元删除，要分离的公共面两侧的实体单元进行分离，循环重复步骤s1
‑
s8，直至完成湿度运移和干缩开裂模拟。
13.进一步地，步骤s3包括以下步骤：
14.s3.1根据高斯散度定理，求得所述多面体实体单元内的湿度梯度，湿度梯度为：
[0015][0016]
式中，w为湿度，v是多面体实体单元的面积，是多面体实体单元中节点l所对的面的单位外法向向量，w
l
是多面体实体单元中节点l的温度，s
(l)
是多面体实体单元中节点l所对的面的面积；
[0017]
s3.2设多面体实体单元内的湿度运移与湿度梯度成正比关系，通过湿度梯度得到单位时间内沿着i方向单位截面积的水的质量流量m
i
为：
[0018][0019]
式中，k
ij
为湿度传导系数张量，w为湿度；
[0020]
s3.3根据步骤s3.1和步骤s3.2，多面体实体单元流入节点l的水分质量流量为：
[0021][0022]
进一步地，步骤s4包括以下步骤：
[0023]
s4.1多面体实体单元流入节点l总的水分质量流量为：
[0024][0025]
式中，t
i
表示共用节点l的多面体实体单元；
[0026]
s4.2在下个时间步δt内，节点l的湿度可以根据下式更新为：
[0027][0028]
式中，m
s
是整个材料的干质量。
[0029]
进一步地，步骤s5中，收缩应力的计算公式为：
[0030]
δσ
ij
＝
‑
δ
ij
3kαδw；
[0031]
式中：δσ
ij
为收缩应力，i＝1,2,3，j＝1,2,3，k为体积模量，α为收缩系数，δw为湿度改变量，δ
ij
为置换张量。
[0032]
进一步地，步骤s6中，计算应力场可以采用有限元(fem)、混合有限
‑
离散元(fdem)、离散元(dem)、非连续变形分析(dda)、数值流行(nmm)当中的任意一种。
[0033]
进一步地，步骤s1中，所述多面体实体单元可为n面体(n为自然数且大于或等于4)、三维voronoi中的一种或多种进行网格划分。
[0034]
进一步地，步骤s2中，相邻所述多面体实体单元间的所述公共单元为节理单元、界
面单元、cohesive单元中的一种。
[0035]
本发明的实施例提供的技术方案带来的有益效果是：湿度运移和干缩开裂的三维数值模拟方法，可以非常方便地模拟土体或其他材料中的湿度运移，同时由于裂缝扩展沿着实体单元边界，在模拟裂缝起裂和扩展和交汇时，不用追踪裂缝扩展和网格从划分，可以非常方便地模拟任意复杂裂缝的起裂、扩展和交汇。
附图说明
[0036]
图1是本发明提供的湿度运移和干缩开裂的三维数值模拟方法中的湿度运移计算模型。
具体实施方式
[0037]
为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地描述。
[0038]
本发明的实施例提供一种湿度运移和干缩开裂的三维数值模拟方法，包括以下步骤：
[0039]
s1通过网格划分将求解域划分为若干个多面体实体单元，多面体实体单元可为n面体(n为自然数且大于或等于4)、三维voronoi中的一种或多种进行网格划分。
[0040]
s2在相邻所述实体单元间可选择地设置公共单元(在相邻多面体实体单元间设置或不设置公共单元)，公共单元为节理单元、界面单元、cohesive单元中的一种。
[0041]
s3根据所述实体单元内节点之间的湿度传递计算每一所述实体单元内的每一节点的水分质量流量。
[0042]
具体地，步骤s3包括以下步骤：
[0043]
s3.1根据高斯散度定理，求得所述多面体实体单元内的湿度梯度，湿度梯度为：
[0044][0045]
式中，w为湿度，v是多面体实体单元的面积，是多面体实体单元中节点l所对的面的单位外法向向量，w
l
是多面体实体单元中节点l的温度，s
(l)
是多面体实体单元中节点l所对的面的面积；
[0046]
s3.2设多面体实体单元内的湿度运移与湿度梯度成正比关系，通过湿度梯度得到单位时间内沿着i方向单位截面积的水的质量流量m
i
为：
[0047][0048]
式中，k
ij
为湿度传导系数张量，w为湿度，又称为含水率，定义为材料中水的质量除以整个材料的干质量。
[0049]
s3.3根据步骤s3.1和步骤s3.2，多面体实体单元流入节点l的水分质量流量为：
[0050][0051]
s4将步骤s2设置的公共单元视为不存在，如未设置则无需处理，然后进行湿度运移计算；根据所述实体单元内每一节点的水分质量流量，获得每一节点总的水分质量流量；根据有限差分积分策略，利用每一节点总的水分质量流量更新每一节点的湿度，用每一所述实体单元中的各个节点的湿度来表示该实体单元的湿度分布，完成整个所述求解域的湿度运移的三维数值模拟。
[0052]
具体地，步骤s4包括以下步骤：
[0053]
s4.1多面体实体单元流入节点l总的水分质量流量为：
[0054][0055]
式中，t
i
表示共用节点l的多面体实体单元；
[0056]
s4.2在下个时间步δt内，节点l的湿度可以根据下式更新为：
[0057][0058]
式中，m
s
是整个材料的干质量。
[0059]
s5根据所述实体单元的湿度分布计算收缩应力，式中，收缩应力的计算公式为：
[0060]
δσ
ij
＝
‑
δ
ij
3kαδw；
[0061]
式中：δσ
ij
为收缩应力，i＝1,2,3，j＝1,2,3，k为体积模量，α为收缩系数，δw为湿度改变量，δ
ij
为置换张量。
[0062]
s6将收缩应力作为体荷载施加在所述实体单元上计算应力场，计算应力场可以采用有限元(fem)、混合有限
‑
离散元(fdem)、离散元(dem)、非连续变形分析(dda)、数值流行(nmm)当中的任意一种。
[0063]
s7对于设有公共单元的实体单元，根据应力场判断公共单元是否断裂；对于未设有公共单元的实体单元，根据应力场计算实体单元公共面的拉应力和剪应力，确定公共面两侧的实体单元是否要发生分离；
[0064]
s8将断裂的公共单元删除，要分离的公共面两侧的实体单元进行分离，循环重复步骤s1
‑
s8，直至完成湿度运移和干缩开裂模拟。
[0065]
下面通过一个具体的例子来详细阐述本发明所述的湿度运移和干缩开裂的三维数值模拟方法。
[0066]
本实施例中，以多面体实体单元为四面体单元为例来对求解域进行网格划分。
[0067]
构建图1所示的3d湿度扩散模型，将与节点连接的所有四面体单元质量的1/4等效到该节点，用节点的湿度来代表整个区域的湿度分布，其中位于四面体单元内任意一点的湿度可以根据四面体四个节点的湿度线性插值求得。图1中，节理单元的交汇节点便是四面体和节理单元的公共节点，如图1中的节点1
‑
7，用这些离散节点的湿度来表征整个连续体的湿度场。下面介绍如何基于图1所示的拓扑连接来计算整个连续体的湿度场。
[0068]
以图1中的节点1为例，与节点1连接的四面体单元有δ1236，δ1346，δ1456，δ
1526，δ1237，δ1347，δ1457，δ1527共8个四面体单元，图1中的8为相邻多面体实体单元间设置的公共单元。由于与节点1直接相连的节点2、3、4、5、6、7的湿度可能与节点1的湿度不同，因此在这些区域可能会发生湿度扩散。以与节点1相连的其中一个四面体单元δ1236为例。假设四面体单元内湿度场分布服从线性分布，那么同一四面体单元内任意一点的湿度梯度为常量，可表示为：
[0069][0070]
由高斯散度定理，(1)式可写为：
[0071][0072]
式中，v是四面体单元的面积，是四面体的节点l所对的面的单位外法向向量，w
l
是四面体单元中节点l的温度，s
(l)
是四面体单元中节点l所对的面的面积。
[0073]
设单元内的水量运移与含水率梯度成正比关系，结合(2)式，单位时间内沿着i方向单位截面积的水的质量流量m
i
(i＝1,2,3)为：
[0074][0075]
式中，k
ij
为湿度传导系数张量，w为含水率，定义为土体中水的质量除以土体的干质量。
[0076]
于是，单位时间内，从四面体单元δ1236流入节点1的水分质量流量可通过下式计算：
[0077][0078]
式中，为节点1所对的面的外法向单位向量。
[0079]
这样便可求得四面体单元δ1236中流入节点1的水分质量流量m
δ1236
→1。与上类似，我们可求得其他与节点1直接相连的四面体单元流入节点1的水分质量流量。这样单位时间内，流入节点1的总的水分质量流量可表示为：
[0080][0081]
于是在下个时间步δt内，节点1的湿度可以根据下式更新为：
[0082][0083]
以上便更新节点1在下一个时间步的湿度，其余离散节点湿度的更新可根据上述方法类似求得。据此，我们便可获得整个求解域湿度场的演化。
[0084]
于是根据湿度分布按照如下公式计算收缩应力为：
[0085]
δσ
ij
＝
‑
δ
ij
3kαδw
ꢀꢀꢀꢀꢀ
(7)
[0086]
将该收缩应力作为体荷载施加到四面体单元上，收缩应力引起的等效节点力为：
[0087][0088]
式中，为四面体单元中节点l所对的面的外法向单位向量，s
(l)
为四面体单元节点l对的面的面积。
[0089]
然后，以混合有限
‑
离散元方法为例，计算应力场并判断节理单元或界面单元或cohesive单元是否断裂。
[0090]
最后，根据以上步骤循环，即可以完成土体或其他材料内湿度运移和干缩开裂计算。
[0091]
以上所述的具体实施，仅以四面体单元作为例子，其余以五面体、六面体、三维voronoi单元、其他任意多面体等其中一种或者多种实体单元，计算流程与上述类似，均包含在本发明的保护范围之内。
[0092]
在不冲突的情况下，本文中上述实施例及实施例中的特征可以相互结合。
[0093]
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。