多物理场作用下的混凝土重力坝坝体应力分析方法与流程

1.本发明属于水电工程技术领域，涉及一种多物理场作用下的混凝土重力坝坝体应力分析方法。


背景技术：

2.重力坝以其结构安全可靠、施工便捷，且适用于大多数地形地质环境等优点被广泛运用于水电工程中。重力坝拦水蓄水量较大，但坝体一旦变形破坏，溃决失事便会带来严重后果。混凝土是一种典型的多孔介质，依据多孔介质的传热传质原理对于混凝土湿热耦合传输及变形数值模拟计算，然后进行坝体温度场及应力场的分析，可以掌握坝体应力在不同温度下的演变规律。因而对重力坝进行准确多物理场分析，并采取合理的渗控措施，对大坝的稳定及安全性具有重要作用。重力坝中渗流场与应力场具有耦合效应。渗流通过施加某作用面上的渗透压力和在渗流区域内分布的渗流体积力，影响坝体的应力分布；应力通过改变坝体的体积应变及孔隙率而影响坝体的渗透参数，从而影响坝体渗流场。


技术实现要素：

3.本发明的目的是提供一种多物理场作用下的混凝土重力坝坝体应力分析方法，利用该方法可以掌握在不同物理场作用下大坝坝体应力情况及其变形规律，从而为进一步的评估坝体结构安全提供重要依据。
4.本发明所采用的技术方案是：
5.多物理场作用下的混凝土重力坝坝体应力分析方法，步骤包括：
6.步骤一，简化混凝土重力坝的几何模型与网络划分；
7.步骤二，与网络划分单元几何中构建温度场、渗流场作用下的坝体应力分析模型，对坝体应力分析模型进行计算得到混凝土温度和渗透系数；
8.步骤三，通过应力分析模型，建立渗流场、温度场、应力场耦合分析方程，通过耦合分析得到重力坝的应力分布。
9.本发明的特点还在于：
10.步骤一具体的为：根据重力坝的结构平面图，将坝顶、坝体中间部分、坝底简化成简单几何面；采用分区域剖分的方法，利用自由四面网格对整个坝段二维模型进行剖分。
11.步骤二具体步骤为：坝体应力分析模型包括热传导方程、饱和-非饱和渗流模型。
12.热传导方程为：
[0013][0014]
式中：α为混凝土导热系数，α＝λ/cρ，λ为混凝土导热系数，c为混凝土比热，ρ为混凝土密度，θ为混凝土绝热温升(℃)；
[0015]
上述热传导方程的三个边界约束条件为：
[0016]
第一类边界条件是混凝土的表面温度t为已知函数；
[0017]
第二类边界条件是大坝混凝土表面的热流量为与时间有关的已知函数为：
[0018][0019]
第三类边界条件：大坝混凝土表面的热流量、大坝混凝土的表面温度t以及与混凝土接触的空气温度t
α
之间的差值成正比例关系，β取0.9：
[0020]
通过边界约束条件对热传导方程进行求解，最终求得混凝土温度。
[0021]
饱和-非饱和渗流模型为：
[0022][0023]
式中，ρ为流体密度；cm为容水度；g为重力加速度；se为饱和度；s为储水系数；p为压力；κs为饱和渗透率；η为流体粘度；kr为相对渗透率；d为位置水头；qm为液体源汇项。
[0024]
饱和-非饱和渗流模型的约束条件为：
[0025][0026]
式中，n为边界面外法线方向；z为轴坐标；n0为内向通量；rb为外部电导率；hb＝zb pbρg为外部总水头；zb为外部高程；pb为外部压强；h＝z pρg为总水头。式(6)的混合边界条件在rb＝0时变为第二类边界条件；在rb＝∞时变为第一类边界条件。式(6)中的条件来分割边界条件，pb＝0，zb＝z，出渗面n0＝n，出渗面以上边界n0＝0；p≥0且rb(hb－h)＜0时，rb＝∞，否则rb＝0；
[0027]
通过饱和-非饱和渗流模型的约束条件对饱和-非饱和渗流模型进行求解，最后求得渗透系数。
[0028]
步骤三中渗流场、温度场、应力场耦合分析方程为：
[0029][0030]
其中u为t时刻下的应力，ρ为混凝土密度，s为单位应力，fv为模型受力向量和，ε
p
为t时刻渗透系数，k为水利传导率，μ为动力粘度，p为水头压力，c为混凝土比热，t为混凝土温度，c
p
为水的比热容，q为计算热量，q
ted
为t时刻剩余热能量。在comsol multiphysics软件中对渗流场、温度场、应力场耦合分析方程进行耦合分析后得到应力分布。
[0031]
本发明的有益效果是：
[0032]
(1)由于混凝土的弹模与浆砌石相比差别很大，从弹性力学的观点来看，坝体存在的这种成层异弹模问题，尤其对坝踵的应力影响显著，甚至会因此导致开裂，另外重力坝坝缝之间填充有材料，且坝体之间还有相互作用力，这些均对三维模型的建立和网格剖分造成了很大难度，为此本发明采用了自编的有限元网格剖分方法，很好地处理了这些问题。
[0033]
(2)通过三场耦合、不考虑温度场以及不考虑渗流场三种情况的比较，在施工期中温度场可引起整个重力坝收到较大的水平向拉应力。而坝体第一主应力在不考虑温度场之
后，最大拉应力值减小了5％左右；在不考虑渗流场的情况下，z向最大拉应力减小0.5％左右。在三种情况下相同坝体代表点应力变化趋势基本相同，但应力值有差异，接近实际情况。
附图说明
[0034]
图1是本发明多物理场作用下的混凝土重力坝坝体应力分析方法的流程图；
[0035]
图2实施例大坝截面图；
[0036]
图3单元体温度计算结果图；
[0037]
图4单元体渗流系数计算结果图；
[0038]
图5多物理场作用下的混凝土重力坝应力分布图。
具体实施方式
[0039]
下面结合附图和具体实施方式对本发明进行详细说明。
[0040]
多物理场作用下的混凝土重力坝坝体应力分析方法，步骤包括：
[0041]
步骤一，简化混凝土重力坝的几何模型与网络划分；
[0042]
步骤二，在网络划分单元几何中构建温度场、渗流场作用下的坝体应力分析模型，对坝体应力分析模型分别进行计算，得到渗透系数和混凝土温度；
[0043]
步骤三，通过应力分析模型，建立渗流场、温度场、应力场耦合分析方程，通过耦合分析得到重力坝的应力分布。
[0044]
步骤一具体的为：根据重力坝的结构平面图，将坝顶、坝体中间部分、坝底简化成简单几何面；采用分区域剖分的方法，利用自由四面网格对整个坝段二维模型进行剖分。
[0045]
步骤二中的坝体应力分析模型包括热传导方程、饱和-非饱和渗流模型。
[0046]
步骤二中的的热传导方程为：
[0047][0048]
式中：α为混凝土导热系数，α＝λ/cρ，λ为混凝土导热系数，c为混凝土比热，ρ为混凝土密度，θ为混凝土绝热温升(℃)；
[0049]
上述热传导方程的边界约束条件为：
[0050]
其中第一类边界条件是混凝土的表面温度t为已知函数。
[0051]
第二类边界条件是大坝混凝土表面的热流量为与时间有关的已知函数为：
[0052][0053]
在计算大坝混凝土内部温度时，假设表面绝热，则有
[0054]
在计算大坝混凝土的表面存在与空气相接触的部分，则经过大坝混凝土表面的热流量可表示为：
[0055][0056]
第三类边界条件：大坝混凝土表面的热流量、大坝混凝土的表面温度t以及与混凝
土接触的空气温度t
α
之间的差值成正比例关系，β取0.9：
[0057][0058]
计算中仅考虑大体积混凝土浇筑过程层面散热。若考虑流水养护条件，可将层面散热系数设置为比较大的数值，外界温度可考虑为养护水温和气温的平均值。
[0059]
通过描述饱和—非饱和多孔介质中流体流动的richards方程，来建立饱和—非饱和渗流模型如下：
[0060][0061]
式中，ρ为流体密度；cm为容水度；g为重力加速度；se为饱和度；s为储水系数；p为压力；κs为饱和渗透率；η为流体粘度；kr为相对渗透率；d为位置水头；qm为液体源汇项。
[0062]
饱和-非饱和渗流模型的约束条件为：
[0063][0064]
式中，n为边界面外法线方向；z为轴坐标；n0为内向通量；rb为外部电导率；hb＝zb pbρg为外部总水头；zb为外部高程；pb为外部压强；h＝z pρg为总水头。式(6)的混合边界条件在rb＝0时变为第二类边界条件；在rb＝∞时变为第一类边界条件。式(6)中的条件来分割边界条件，pb＝0，zb＝z，出渗面n0＝n，出渗面以上边界n0＝0；p≥0且rb(hb－h)＜0时，rb＝∞，否则rb＝0；
[0065]
先对坝体作简化假设：
①
坝体结构按线弹性均质材料考虑，水的渗流运动符合达西定律；
②
取一个坝段独立计算，不考虑相邻坝段的影响；
③
材料总的应变是应力引起的应变和水压力导致的应变之和；
④
在基础计算范围内，认为底部固定，坝基为不透水层。这样，得到重力坝二维应力平衡方程：
[0066][0067]
式中，σ
x
、σy、τ
xy
、τ
yx
为应力；x、y为x、y方向的体力。根据有效应力原理及弹性本构关系，改写为位移分量和孔隙水压力表示的平衡方程为：
[0068][0069]
式中，e为坝体混凝土弹性模量23gpa；μ为坝体混凝土泊松比0.167；x0、y0为由初应变引起的等价体力；e为体积应变，通过以下公式获得：
[0070][0071]
公式(8)中为x、y方向的梯度计算因子通过以下公式获得：
[0072][0073]
用总水头h表示孔隙水压力p，根据h与p的关系可得公式(11)。式中，γw为水的容
重：
[0074]-p＝γw(h-y)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)
[0075]
为方便起见将γ
wh
记为h，即总水头的单位与应力单位相同，将式(11)代入式(8)中即为位移与水头表示的本构方程：
[0076][0077]
根据质量守恒定律及达西定律，得到流体的连续型方程：
[0078][0079]
其中渗透系数k(s)：
[0080][0081]
式中，k
x
、ky为x、y向的渗透系数；θ为坝体混凝土的孔隙率；β为水的体积压缩系数，等于水的弹性模量的倒数，一般取1.0；hd为进气吸力；n为形状系数，m＝(l 1)n 2；l为孔隙弯曲度因子，取2。
[0082]
通过上述计算对饱和-非饱和渗流模型进行求解，最后求得渗透系数k(s)。
[0083]
步骤三中在comsol multiphysics软件中，进行渗流场、温度场、应力场耦合分析，偏微分方程组为基础，建立混凝土重力坝断面二维渗流场、温度场、应力场三场耦合分析方程：
[0084][0085]
其中u为t时刻下的应力，ρ为混凝土密度，s为单位应力，fv为模型受力向量和，ε
p
为时刻渗透系数，k为水利传导率，μ为动力粘度，p为水头压力，c为混凝土比热，t为混凝土温度，c
p
为水的比热容，q为计算热量，q
ted
为t时刻剩余热能量。
[0086]
在comsol multiphysics软件中对渗流场、温度场、应力场耦合分析方程进行耦合分析后得到应力分布。
[0087]
实施例
[0088]
本实施例重力坝坝体的实际结构如图2大坝截面图，对图2上坝体进行应力分析：
[0089]
步骤一，将大坝按照实际结构分为两部分：大坝主体和坝顶，其中大坝主体为六面体，横截面为直角梯形，坝顶为长方体。
[0090]
对剖面模型的网络划分
[0091]
采用分区域剖分的方法，利用自由四面体网格对整个坝段二维模型进行剖分。最大单元11.2m，最小单元0.11m，最大单元增长率1.3，曲率因子0.2，狭窄区域分辨率1。针对模型中坝顶区域及挡水面坝底位置，网格进行加密，最大单元1m，最小单元0.01m，最大单元
增长率1.3，曲率因子0.2，狭窄区域分辨率1。因为几何尺寸较小且结构密度较高造成部分应力场和温度梯度分布较大。生成网格，将该模型划分为2043个节点，3888个单元。考虑模型边界条件为渗流边界条件、应力边界条件和位移边界条件三种。在正常蓄水情况下，渗流边界条件包括水头边界、流量边界和混合边界：
[0092]
(1)水头边界为坝体上游面和下游面水位线以下部分；
[0093]
(2)流量边界为坝底与坝基接触面，坝底为不透水边界；
[0094]
(3)混合边界为溢出边界和渗流自由面，位移边界条件为底面完全位移约束，应力边界为考虑自重荷载和渗透压力。
[0095]
由于计算需要，对建好的坝体进行细节调整，即为移除细节，需要移除的细节包括连续相切的顶点、短边、小面、长条面、窄边区域以及薄域。累计共忽略3个顶点，3条边以及一个塌陷域，在参数的选择中，细节大小选择为自动，连续相切容差为5deg。
[0096]
执行步骤二
[0097]
根据公式(1)建立热传导方程，t为已知函数，即t(τ)＝37℃，c为880[j/(kg*k)]；ρ为2300[kg/m^3]；θ为45.38℃。
[0098]
计算出单元体的温度分布图，如图3。根据公式(5)建立渗流方程，ρ为流体密度1000[kg/m^3]；cm为容水度；g为重力加速度；se为饱和度1；s为储水系数0.5；p为压力1mpa；κs为饱和渗透率0.5；η为流体粘度0.6；kr为相对渗透率0.1；d为位置水头10m；qm为水源汇项取1。
[0099]
计算出单元体的渗流系数及方向分布图，如图4。
[0100]
执行步骤三建立渗流场、温度场、应力场耦合分析方程，根据实施例模型给公式(15)代入参数，可得：
[0101][0102]
基于comsol的二维几何的坝体应力变形与分析，代入多物理场方程计算坝体刨面应力分布，到最后加载完成，得到坝体应力较稳定的呈现规律。对水电站混凝土坝模型渗流、温度、应力场的综合研究整体运算结果如图5，整体应力在0.6mpa左右，与设计规范要求的0.5～0.7mpa整体符合，整个重力坝应力分析方法具有实际参考价值，对后续维护和监测提供理论指导。