蒋中明，邓自源，唐 栋 ,3，万 发，欧阳钰榕

(1.长沙理工大学水利工程学院，湖南 长沙 410114；2.长沙理工大学水沙科学与水灾害防治湖南省重点实验室，湖南 长沙 410114；3.长沙理工大学洞庭湖水环境治理与生态修复湖南省重点实验室，湖南 长沙 410114)

气垫式调压室作为一种性能优越的水锤和涌浪控制设施，相较于常规调压室，在交通、环境以及投资等方面更具有优势[1-2]。电站在变负荷工况下会发生水锤现象，调压室水面上下波动，气室气体因此不断膨胀、收缩，压缩空气温度和压力发生循环升降变化[2-3]。调压室内压缩空气压力和温度的变化不但影响到调压室的水力特性[4]，同时也会改变围岩的应力变形状态，进而影响围岩的稳定性。

目前国内外关于气垫式调压室的研究重点主要集中在调压室的水位变化规律和压力变化过程等水力过渡过程特性方面[5-7]，对水位波动过程中调压室内压缩空气的温度变化较少关注。实际上空气的温度也随着压力的改变而变化，并将影响到围岩和密封层等结构的受力和变形状态[8-11]。为分析压缩空气热力学过程对调压室围岩应力和变形的影响，需要准确获得调压室水位变化条件下压缩空气的压力及温度的动态变化过程，为此，笔者以考虑传热条件下的热力学控制方程为基础，根据气垫式调压室水位波动变化特点，建立适用于气垫式调压室的压缩空气热力学分析模型，并利用所建立的模型计算调压室内的空气压力和温度变化过程，进而将空气压力和温度作为调压室围岩应力变形分析的边界条件，对气垫式调压室围岩在瞬态热传导和压力变化工况条件下的受力特性进行研究。

1 气垫式调压室压缩空气热力学模型

1.1 气垫式调压室水位波动方程

水电机组突然甩负荷工况下，其流量由Q0突变减至0，以调压室初始水位为基准水位，建立气垫式调压室压力变化分析的运动方程、连续性方程、气体状态方程[3]分别如下：

(1)

(2)

(3)

式中：L、f、v分别为引水隧洞的长度、断面面积和断面流速；F为调压室水平断面面积；p0、l0分别为机组在额定工况下，初始时刻t=0 s时气垫式调压室内气体初始压力和气室初始高度；p为任意时刻气垫式调压室内气体压力；α、β分别为引水隧洞水头损失系数和阻抗孔水头损失系数；n为气体多方指数；t为时间；z为水位；g为重力加速度。

气垫式调压室水位波动[12]建议方程如下：

z=Acosφ

(4)

(5)

φ=wt+φ0

(6)

(7)

(8)

式中：ρ为水的密度；v0为初始时刻引水隧洞中的初始流速；A、φ分别为水位偏差值和相位；A0、φ0分别为初始时刻水位偏差值和相位。

1.2 气垫式调压室空气压力及温度计算方程

利用式(1)至式(4)得到调压室内空气压力时，没有考虑空气压缩引起的温度变化对压力的影响。研究表明调压室内的空气温度和压强的变化过程具有高度复杂的相关性[13-16]。为获得调压室内水位变化条件下空气压力和温度计算的新方法，参考压气储能地下储气库内压缩空气热力学过程的分析方法[14]，建立考虑空气与围岩和水体间的热量相互传递影响的气垫式调压室的压缩空气温度和压力计算模型。根据热力学第一定律，调压室内部空气的能量守恒方程为

(9)

单位质量空气内能可用下式表示：

(10)

将式(10)两边乘以m，然后对时间t求导，可得

(11)

利用体积做功关系式dW=-pdV和式(11)，可将调压室内空气能量守恒方程(9)改写为

(12)

空气与接触界面(壁面和水面)之间的对流换热速率计算表达式为

(13)

(14)

式中：hcw、hcr分别为空气与围岩和水体间的热交换系数；Acw、Acr分别为气体与围岩、水间的接触面积；Tsw、Tsr分别为围岩和水的温度。

调压室中的气体状态方程为

pV=ZmRT

(15)

式中：Z为空气压缩因子，假定调压室气体为理想气体，取压缩因子Z=1；R为气体常数。将式(15)左右两边分别对t求导，有

(16)

由(14)和式(16)得压力对时间的导数表达式为

(17)

式(14)和式(17)即为调压室内压缩空气的压力和温度动态变化方程。由此可知，调压室内空气压力变化与调压室内空气占据的体积空间V(可由水位z和调压室截面积计算求得)和温度密切相关。计算体积V所需的z值可以由式(4)得到。利用式(4)计算水位z时所需的压力由初始条件和式(17)提供，而不是利用式(3)计算气室内的气体压力。

在对调压室围岩应力及变形特性进行分析时，需将式(14)和式(17)计算得到的压力和温度值分别作为压力边界和对流换热边界条件直接施加到调压室的洞壁上作为计算边界，从而可实现对围岩的单向热力耦合变形与受力特性分析。

2 模型合理性分析

为验证式(14)和式(17)的合理性，假定气垫式调压室不发生气体渗漏，即m为常数的条件下，以小天都气垫式调压室为研究对象，以蒋中明等[17]建议的差分法对式(14)(17)进行求解及合理性分析。表1为小天都水电站气垫式调压室的运行参数[18]，热学计算参数取值如表2所示。

表1 气垫式调压室计算参数

表2 热学计算参数

图1(a)为利用式(17)与式(3)得到的调压室内气体压力变化过程对比，在同一水位波动条件下，采用两种方法得到的气体压力变化规律基本相同，压力最大值分别为4.57 MPa、4.55 MPa，最小值分别为3.14 MPa、3.11 MPa。虽然两者在数值上存在一定的差异，但偏差较小，均方差仅为0.059 1，且均与杜鹏侠等[18-19]监测得到的压力变化范围3.28～4.44 MPa较为接近，表明笔者提出的调压室内压缩空气压力的计算方法是合理的。造成压力值偏差的原因是在利用式(3)进行求解计算时仅依据气体多方指数的取值来反映压力的变化规律，没有完全考虑气体温度变化对压力的变化的影响；而实际运行过程中压力与温度之间是相互影响的。正是由于笔者提出的压缩空气压力计算模型考虑了温度变化对压力的影响，从而使得两种方法得到的计算结果之间存在一定偏差。

图1 气体状态演化过程计算结果分析

图1(b)为根据式(14)计算得到的空气温度变化过程图，调压室内气体的温度变化趋势与压力变化趋势基本相同，其中空气温度的最高、最低值分别为33.3℃、1.2℃。由图1(a)可知，调压室内空气的峰值压力在初始值3.68 MPa基础上增加了约0.89 MPa，空气的峰值温度在初始值15℃基础上增加了约18.3℃，与Lux[19]对地下洞室压缩空气的温度变化研究结果(即气体压力每增加1MPa气体温度即升高约13℃)基本一致，表明笔者提出的调压室内压缩空气压力和温度变化分析的热力学方程是合理的。图1(b)还表明，在运行过程中气垫式调压室的气体温度变化幅度达到了30℃以上，同时最低温度接近于0℃，存在温度较低的现象，因此在评价调压室结构性能与安全时，不可忽略温度变化对结构应力的不利影响。

3 热力耦合作用下调压室围岩受力特性分析

为分析气垫式调压室围岩在调压过程中的受力及变形演化特性，基于3DEC软件平台，编写调压室内气体压力和温度计算的FISH程序，采用热力耦合数值分析方法对小天都气垫式调压室围岩受力和变形特性进行数值计算。

3.1 数值模型及计算条件

小天都气垫式调压室为圆拱直墙形，尺寸为80 m×16 m×16.5 m，围岩以II～III类为主。在模型建立时，考虑部分主要裂隙的影响，其空间分布如图2(a)所示，裂隙总数为35条，呈椭圆盘状；数值计算模型尺寸为300 m× 200 m×200 m，网格计算单元数量为670 970个，如图2(b)所示。为便于分析，在调压室周围设置若干监测点，位置如图3所示。

图2 三维数值模型

图3 洞室围岩测点分布

计算参数采用小天都电站相关地质勘测资料推荐参数[19-21]；热力学计算参数参考Zhou等[22-23]论文数据，具体取值见表3、表4。围岩的力学本构模型采用摩尔库伦弹塑性本构模型，热力学本构模型采用各向同性热传导模型。模型顶部为压力边界，其余边界均设置位移约束；调压室水位以上部分内表面压力边界和对流换热边界。按自重应力计算初始应力场，按15℃计算初始温度场。

表3 物理力学参数

表4 热学计算参数

工程实践表明：当水电站突然完全甩负荷时，水流涌向调压室，一般按水库最高设计水位计算调压室水位最高涌波水位；当水电站突然增加全负荷时，计算调压室的最低涌波水位。因此，计算工况采用完全甩负荷工况，此时调压室内的水位最高，空气压力和温度最高，对结构影响最为不利。

3.2 调压室围岩温度时空分布特性

图4为洞壁处测点的温度变化过程，其中D组测点位于水下，不受调压室内温度的影响。由图4可知，随着调压室内水位波动变化，围岩温度也随着气体温度的变化而变化；尽管洞壁测点所处的位置不同，但其温度变化规律却基本相同；不同部位测点的温度数值存在一定的差异，其中左侧洞壁A1点的温度相对于其他测点的温度更高，极大值达到16.4℃。其原因是距A1点不远处存在节理裂隙交汇点，而温度在裂隙中传播速度比在壁面中的传播速度更快[24]，故靠近裂隙的A1点呈现出较高的温度。此外，B1和C1测点的温度在工况发生变化后的第5min达到极大值15.95℃左右，且将在该温度状态维持一定的时间。

图4 30 min内调压室各组测点洞壁处温度变化过程

图5为调压5 min末期时调压室横剖面围岩的温度分布。随着气体温度的升高，调压室内壁直接与气体进行热交换，导致围岩内壁的温度升高至16.0℃左右。由于节理裂隙的影响，围岩温度分布并不均匀，局部区域特别是节理裂隙交汇点处的最高温度升高到16.64℃，见图中N′、N处。同时，在距洞壁一定范围0.5 m内围岩的温度变化梯度较小，与内壁温度基本保持一致，但超过0.5 m后，热能向围岩深处传递的阻力越来越大，围岩的温度变化梯度先增大后减小，温度逐渐降低至与远处围岩温度一致。由于水的存在，调压室底部围岩没有与气体接触进行热交换，故其温度受到气体温度变化的影响较弱。

图5 围岩的温度分布

3.3 围岩应力变形特性

图6给出了调压室工况改变后第5 min横剖面上的第一主应力和第三主应力分布。由图6(a)可知，调压室内壁附近区域的第一主应力基本维持在-6～-4 MPa，且在调压室顶部左侧节理附近及底部拐角区域出现了点状区域的应力集中的情况，局部甚至达到-10 MPa，但其分布范围较小。洞室顶部右侧节理处出现了高达1.26 MPa的拉应力。调压室边墙处的第三主应力均表现为压应力，其压应力最大值为-2.38 MPa。由于上部围岩及底部围岩受开挖卸荷作用的影响，调压室顶拱和底板附近的围岩均出现了一定程度的拉应力，且在洞室顶部节理处存在拉应力集中的情况，其最大值达到1.26 MPa。结果表明洞室围岩的节理裂隙对围岩的应力状态有较大的影响，在进行地下洞室设计时应考虑节理裂隙对洞室稳定性的影响。

图6 围岩应力分布

图7给出了气垫式调压室工况改变后5min时各组测点的位移与距洞壁距离的关系。调压室围岩的位移变化的总体规律为：距离洞壁越远，位移越小；其中底部围岩(D组)位移相较于其他部位的位移值更大，而顶部围岩(B组)的位移值则较小；两侧洞身围岩的位移值大小基本一致。

图7 测点位移与洞壁距离的关系

3.4 测点温度、应力及变形的变化过程

图8为B组测点中的部分测点温度变化过程。随着调压室内压缩空气温度和压力的循环波动变化，围岩的温度也呈现出相应的波动变化特征。调压室洞壁附近区域(0.5 m范围内)的温度受到调压室空气温度的影响较为显著，如B1、B3点；距离洞壁较远区域的温度受到的影响较小，温度上升较为缓慢；当距离达到5 m后，由于距离远、热阻大等原因导致围岩温度基本不受空气温度的影响，如B9点。总体上看，调压室内压缩空气温度的显著变化对洞周围岩的温度绝对升温幅度影响较小。

图8 围岩测点的温度变化过程

图9～11分别给出了调压室B组部分测点在考虑和不考虑温度荷载条件下的位移、应力变化过程的对比，图中符号M表示仅考虑力学荷载的计算成果，TM表示考虑温度荷载和力学荷载共同作用的计算成果。由图9可知，两种情况下测点处围岩的位移值基本一致，考虑温度影响条件下的围岩位移仅比不考虑温度时的略大0.01 mm，表明压缩空气温度对围岩的位移值影响不大。图10和图11表明靠近洞壁附近的围岩应力受温度荷载影响比较显著。在考虑温度荷载条件下，围岩内壁附近第一主应力较不考虑温度荷载时增加了-0.14 MPa，为压应力。相较于不考虑温度荷载的情况，考虑温度荷载条件下的洞壁附近围岩第三主应力在调压过程中先降后增的变化趋势，其第三主应力最大差值较不考虑温度荷载时的应力值减少约0.14 MPa，为拉应力。第三主应力在4.5 min时达到极小值后又呈现出波动起伏且上升的变化趋势，在30 min后达到相对稳定，并接近于不考虑温度荷载下的应力值。

图9 典型测点位移变化过程对比

图10 典型测点第一主应力变化过程对比

图11 典型测点第三主应力变化过程对比

4 结 论

a.气垫式调压室运行工况发生改变时，气体温度可能出现大幅度的变化现象；温度的大幅度变化将对调压室结构的应力应变状态产生不利影响，故分析调压室密封性能和围岩稳定性时应考虑包含气体温度和压力的热力耦合效应影响。

b.气体温度的变化是引起围岩温度变化的直接原因，围岩各区域温度分布呈现出显著的不均匀性，节理裂隙处的温度相对其他区域更高；围岩温度影响范围一般在距内壁0.5 m范围内。

c.相较于不考虑温度荷载的影响工况，温度荷载对围岩的位移值影响不大，但对围岩应力(特别是第三主应力)有较大的影响。