流固耦合基础知识（精选7篇）
流固耦合基础知识 篇1
流固耦合求解体会
对大多数应用，流体和结构网格在流固耦合界面上是不同的，因为流体和结构网格有不同的分辨要求。界面耦合法支持通过不同网格界面的载荷传递。界面耦合法对整个搜索使用线性搜索算法来定位属于一定的有限元网格。
流体—结构耦合分析有以下5个主要步骤。
1）设置流体和结构分析
进行流体—结构耦合分析，首先需要创建流体模型及流体区域的有限元网格；然后，施加合适的边界条件，设置流体分析的材料属性，选择合适的求解选项，创建结构模型及固体区域的有限元网格，施加固体区域合适的边界条件，设置合适的求解选项。
2）标记流体和结构界面
这一步骤是界面耦合法的特色之一，对流体—结构界面加以标记，在该界面发生载荷传递。用 SF 族命令（SF,SFA,SFE,SFL）及 FSIN 在流体结构界面作两次标记，一次是在界面的流体侧，一次是在界面的固体侧。载荷传递将发生在有相同界面号的流体和结构界面处。
3）指定流体结构耦合求解选项
在这一步骤中要设置流体和结构分析的求解顺序，指定是静力还是动力分析，还必须指定载荷传递的插值方法。然后还要根据实际情况设置以下选项：终止时间、载荷步时间及时间增量，交错迭代的最大次数，收敛值，输出频率，等等。
4）获得求解
求解使用 SOLVE 命令。也可以从上一个时间步或结果文件中最后一次的收敛结果重新启动流体结构耦合分析。但须注意，流体结构耦合分析不支持结构分析中的多重启动功能。
5）对结果进行后处理
对流体结构耦合分析进行后处理，必须重新读入数据库。可以用标准的 ANSYS POST1 命令和 POST26 命令来查看结果，但不支持同时对流体和固体文件的后处理。一次只能读入流体或结构的结果文件并进行后处理。流体分析的结果存放在 Jobname.RFL 文件中，而固体分析的结果文件存放在 Jobname.RST 或 Jobname.RTH 文件中。
流固耦合基础知识 篇2
微流体系统中大都要求有微泵。其中,薄膜往复振动式无阀微泵是以扩散/收缩单元为阀门,通过泵膜的振动驱动流体,其制造简单,可以驱动一些非均相的流体,是众多研究的焦点。根据驱动原理不同,无阀微泵可划分为电磁、静电、压电、形状记忆合金等多种类型。其中,静电无阀微泵具有较低功耗,易于与IC工艺集成兼容等优点,而得到广泛关注。
由于微泵的尺度很小,薄膜往复振动式无阀微泵的绝大部分物理量难以通过实验测量,特别是瞬态量。因此,在微泵的研制阶段,利用数值方法对动态性能进行仿真、预测显得很重要。它不仅可降低费用,而且能更好地了解微泵的工作原理以及可能的潜在缺陷。但因微泵中柔性泵膜、电驱动和流体相互耦合,加之流体在扩散/收缩单元中的不同流动方向表现出的不同压力特性,增加了对整个泵的仿真难度和复杂性。早期的研究采用低阶集总参数模型[1,2]和等效电网络模型[3]等简化方法,忽略了参数的空间分布特性。目前大部分的研究没有考虑泵的流固耦合效应,也无法建立电驱动信号量与流体动力学量(比如流量)之间的关系。或是研究驱动器的动态特性[4,5,6,7];或是简化流体场,根据扩散/收缩单元的压力损失系数,将泵膜的运动与泵腔内外压力以及扩散/收缩单元流量联系起来,研究泵膜的动态特性[8,9,10,11];或是将驱动器的位移输出,作为流体场的移动壁面问题来研究泵流体的动力学特性[12,13,14]。实际上,泵膜的动态特性与流体黏滞损失的非线性和不稳定性等是相互影响的。为了理解整个微泵的动态特性,模型必须能有效地描述这些耦合效应[15,16]。本文以静电无阀微泵为对象,建立静电-结构-流体全耦合的3维模型,利用数值方法,仿真并分析各个变量的非线性动态特性。
1 模型的建立
1.1 控制方程、边界条件以及初始条件
图1所示为静电无阀微泵的结构示意图。对于泵送流体,采用任意拉格朗日-欧拉(ALE)描述比较方便。某一时刻,在参考系中
式中,
流体域Ωf的不可压缩流动动量守恒定律以及连续方程可以表示为
式中,ρf、f、σ分别为流体密度、体积力向量和Cauchy应力张量;[0,T]为所考察的时间间隔。
牛顿流体的本构方程为
σ=-pI+2μ∇s
式中,I为单位张量;μ为流体的黏度;∇s
在泵膜流固耦合界面Γf-s处,边界条件如下:
式中,
式(5)为流固界面Γf-s处的无滑移条件;在流固界面Γf-s上,参考系的速度
对于泵膜,可以利用标准的Lagrangian描述建立其运动学方程,即
式中,ρs为固体的密度。
泵膜流固耦合边界Γf-s满足:
边界条件式(9)、式(10)分别对应于式(5)和式(7)。
如果没有包含自由电荷,则描述静电场的Poisson方程可以表示为
∇·ε∇ϕ=0 (11)
E=(E1,E2,E3)=-∇ϕ
式中,ϕ为静电势;ε为介电常数张量;E为静电场电场强度。
Maxwell静电应力张量TE由下式计算:
式中,FE为静电载荷。
各个场变量的初始条件都设置为零。为了计算易于收敛,将一个很小的初值输入到模型中。
1.2 模型参数及材料属性
图1所示静电微泵的泵膜厚度为10μm,泵腔深度为100μm,泵膜半径为2000μm;扩散/收缩单元的几何结构相同,长度为1000μm,最小端面宽度尺寸为80μm,扩张角为10°;连接的进出口管道长度皆为2000μm,宽度为200μm。流体介质为去离子水,物性参数如表1所示。
假设层流流动,考虑对称性,在CFD-ACE中建立微泵的一半模型。利用结构化网格划分流体,单元数为52 880个,泵膜划分为4974个壳单元。使用基于压力修正的SIMPLE-C 算法对流动以及能量进行顺序积分,获得流场的解。流体域的对流和扩散采用一阶迎风格式,时间积分采用了Crank-Nicolson格式。使用基于有限体积法的计算格式求解静电场,这样可以处理不同介质介电常数问题。利用一阶壳单元来表达泵膜,并进行大变形和接触等几何非线性分析。考虑到静电、结构以及流动耦合求解的需要,利用网格重划分工具,使用预测-校正方法对局部变形进行连续修正,将流动速度和结构速度联系起来,实现流固耦合分析。流体和结构求解器的耦合迭代频率设置为fc=1,使得流体场和结构场信息能及时得到交换和更新。当量纲一残差小于10-4时,认为计算是收敛的。施加的周期性驱动电压V=200(1-cos(20πt))(泵腔底电压为零),在Pentium-3.0GHz的PC机中,完成3.5个周期的求解,共350个子步数,求解时间约为59h。
2 结果及讨论
2.1 泵膜的动力学特性
图2所示为计算得到的泵膜中心点垂直方向上的位移随时间变化曲线,图中最大变形量为28.433μm。考虑流固耦合效应后,在流体阻尼作用下,泵膜无法回复至初始零位置,0.653μm为新的平衡位置。稳态情况下,泵膜位移的变化滞后于输入电压信号的变化,一个周期的滞后时间约为0.01s。
2.2 流场的瞬态特性
实际上,泵流体动态特性更为设计者所关注,为了更好地诠释泵膜振动与流体特性之间的耦合关系,用图3所示的泵膜流固耦合界面上中心点处的垂直方向位移与泵腔内流体压力来显示其对应关系,从图中可以明显地看出,流体压力的变化与泵膜的位移存在着密切关系。在静电力的作用下,泵膜下拉变形的速度大于回弹的速度,所以在流体压出泵腔的过程中其压力变化较大,最高压力为945.1Pa;而在流体吸入泵腔过程中,其压力变化相对较小,最大负压为-558.0Pa;整个泵送过程中,流体压力呈现非对称性变化。图3还表明,泵初始工作点从零点开始,经过第一个周期的瞬态过程后,第二个周期的轨迹与第三个周期的轨迹已经重合,表明泵已经进入稳定工作状态。流体压力的波动周期同样滞后于驱动电压的变化(约0.01s)。
图4所示为进出口单元处的瞬态体积流量变化曲线,从图中可以看出,在泵膜下拉(往固定电极)运动阶段(泵送阶段),出口单元的瞬态流量大于进口单元的瞬态流量;而在泵膜回弹阶段(泵吸入阶段),进口单元的流量则大于出口单元的流量。最大的瞬态流量差发生在泵腔内压力达到极大值时刻,而在泵腔内压力与进出口压力相差不大的时间段,进出口单元的流量相差不大。表明要使扩散/收缩单元发挥“整流”作用,其两端的压力差必须达到一定值,并且压力差要尽可能地大。对稳态情况下的泵净流量进行积分计算,得到的净输出流量为10.764nL/s。
图5所示为微泵横截面方向上流体压力的分布情况,由图可知泵腔内的压力几乎分布一致,说明在平面布置的无阀微泵中,采用Reynolds方程来描述泵腔内流体的瞬态动力学特性具有一定的可行性。另外,微泵所连接的进出口管道对微泵的工作性能也有影响。
考虑了流固耦合特性后,模型不仅能够仿真泵流场的动态特性,还能够仿真泵膜的应力特性变化。图6所示为t=28ms时刻的泵腔内流体压力分布与泵膜应力分布情况,以及这一时刻的扩散/收缩单元内流体的速度场分布。此时,泵腔内流体压力达到最大值945.3Pa,泵膜中心挠度为13.20μm;泵膜边缘处VonMises应力为28.39MPa。微泵腔内的流速较小,扩张管/收缩管颈部的平均流速达到最大,计算得到的平均流速为0.0867m/s,雷诺数为7.71,远小于宏观条件下通常认定的临界雷诺数(2000左右),这也证实了前面假设层流流动模式的正确性。另外,图6b和图6c表明,收缩/扩散管内未出现流体与固体壁面分离流动现象。
3 结论
(1)扩散/收缩单元的最大流量差都发生在泵腔内压力与进出口压力相差最大的时刻,表明要使得扩散单元发挥“整流”作用,其单元两端的压力差必须达到一定值。
(2)在仿真的每个时刻,泵腔内流体压力分布几乎一致,表明用Reynolds方程来描述泵腔内流体动力学特性具有可行性。
(3)仿真得到的最大雷诺数远小于通常认定的临界雷诺数,表明微泵流体具有层流特性。
(4)泵腔内流体动态特性与泵膜的运动存在着密切关系,如果忽略了流固耦合效应,简单地将微泵的驱动器与腔内的流场分离处理,必然歪曲了泵流体场的动态特性。泵膜变形、流体压力和进出口流速等响应与输入电压信号之间存在着滞后现象。对微泵的3D全耦合仿真,能够获得驱动电信号与流体流量、泵送压力等输出变量之间的直接关系,有利于从整体意义上实现微泵的优化设计。
流固耦合基础知识 篇3
摘要：电化学放电加工（Electrochemieal Discharge Machining，ECDM）为近年来较受关注的非传统制造技术之一，既可以加工非金属材料，也可以加工金属材料。以有限元软件对工具电极高速旋转电化学放电加工加工状态进行了仿真分析，研究了工具电极旋转速度及工具电极与微小孔之间的偏心程度对工具电极变形及工具电极应力的影响。计算结果表明，随着工具电极旋转速度增加，工具电极最大变形及工具电极最大应力增加；随着工具电极与微小孔偏心量的增加，工具电极最大变形及工具电极最大应力增加。
关键词：液固耦合；工具电极高速旋转；电化学放电加工
电化学放电加工（E1ectrochemical DischargeMachining，EcDM）为近年来较受关注的非传统制造技术之一，既可以加工非金属材料，也可以加工金属材料。使用高速旋转的工具电极进行电化学放电微孔加工，有利于加工产物从微小的加工间隙内排出。在工具电极高速旋转进行电化学放电加工微小孔的过程中，工具电极由于复杂的受力状况，可以视为一个极其复杂的动力学系统。由于工具电极的旋转作用，加工微小孔时，微小孔内充满了具有一定粘度的工作介质，旋转向下运动的工具电极使微小孔内的工作介质受到动压力，并反作用于高速旋转的工件电极上，从而影响工具电极的运动状态。同时，在微小孔加工过程中，工具电极可能处于偏心工作状态，此时，微小孔内的工作介质对工具电极的作用力大小和方向都会发生改变，工具电极的工作状态反过来又会影响加工间隙内的流场状态，从而对加工过程产生影响。由于微小孔电化学放电加工的特殊性，无法直接使用各种测量设备对加工间隙内的状况进行测量；另外，由于电化学放电加工机理的复杂性，无法求得加工过的解析解，目前的情况下，只能使用模拟仿真的方法，定性的进行分析。
以前受到条件的限制，一般只进行基于单独的流场分析或者是基于结构进行分析，这样的分析过程不能反应多物理场综合作用的效果，以及各物理场之间的相互影响状况，仿真分析的结果具有一定的局限性。随着计算机软硬件技术的发展，多物理场耦合计算理论的进一步完善，针对多物理场耦合问题的计算成为可能。流固耦合分析在其他工程问题有不少的研究成果，但在微细电化学放电高速旋转工具电极应用方面却鲜有文献报道。本文对高速旋转的工具电极流固耦合系统的动力特性进行数值分析，分析工具电极旋转速度、工具电极偏心度及微孔深度对加工过程的影响。
1.计算模型
1.1几何模型及网格划分
因为实际的加工过程很复杂，难以用理论方法进行全面分析，因此对加工过程进行合理的模型简化是非常必要的。通过简化分析，定性的认识工具电极旋转时微小孔内流固耦合特性。本文的几何模型包含微小孔和工具电极两部分，工具电极直径为400μm，当工具电极与微小孔同心时，取微小孔直径为450μm，侧壁加工间隙20μm，端面加工间隙25μm，孔深度为2000μm；当工具电极与微小孔偏心时，偏心距e取值为10微米和40微米进行计算，最小侧壁间隙为20μm，端面加工件25μm，微小孔深度为2000μm，如图1所示。图2a所示为加工的几何模型；图2b为工具电极网格；图2c为流体域网格划分；图2d为流体域网格局部放大图。流体域网格以四面体单元为主，合理的选用结构化六面体单元、楔体单元或者锥体单元，对微孔内壁圆柱面设置膨胀层，进行网格细分。
1.2计算方法
单向流固耦合分析是指耦合交界面处的数据传递是单向的，一般是把CFD分析计算的记过（如力、温度和对流载荷）传递给固体结构进行分析，但没有把固体结构分析结果传递给流体分析过程，只有流体分析对固体结构分析有影响，而结构分析的变形等结果非常小，以至于对流体分析的影响可以忽略不计本文分析的高速旋转的工具电极在工作中受到流体介质和离心力的影响，会导致工具电极发生变形，未考虑变形的工具电极反过来影响加工间隙中的流场，因此选用单向流固耦合分析。
2计算结果与分析
2.1工具电极旋转速度与工具电极最大变形的关系
由图3～图5可知，工具电极最大变形发生在电极端部，夹持部分变形最小。从图6可知，当工具电极旋转速度增加时，工具电极最大变形量增加；随着工具电极与微小孔之间的偏心量增加，工具电极的最大变形量增加；当工具电极与微小孔之间同心时，工具电极的最大变形量非常小，接近于零；当工具电极与微小孔之间偏心量较小而且工具旋转速度相对较低时，工具电极的最大变形量也相对非常小，但旋转速度较大时，工具电极最大变形量增加。
2.2旋转速度与工具电极应力的关系
由图7～图9可知，工具电极最大应力发生在电极夹持部分，电极端部应力形最小。从图10可知，当工具电极旋转速度增加时，工具电极最大应力增加；随着工具电极与微小孔之间的偏心量增加，工具电极的最大应力增加；当工具电极与微小孔之间同心时，工具电极受到的最大应力非常小，接近于零；当工具电极与微小孔之间偏心量较小时，工具电极的最大变形量也非常小，但总体呈现随着工具电极旋转速度的增加而增加。
2.3结果分析
进行电化学放电加工微小孔时，由于工具电极旋转作用，加工间隙内的工作液也随着工具电极旋转而旋转，由于离心力的作用被甩向微小孔内壁的四周，可以认为是随工具电极一起旋转的近似圆柱式的液体层，该液体层对微小孔内壁产生一定的作用力，对工具电极的受力状况及运动形式会产生影响。电化学放电过程中，工具电极可能会处于偏心状态下工作，随着偏心量的增加，作用在工具电极上的弯曲力会增加，因此导致工具电极最大变形和最大应力的增加。
3结论
1）由液固耦合计算可知，随着工具电极旋转速度增加，工具电极的最大变形量和最大应力增加。
2）随着工具电极与微小孔之间的偏心量增加，工具电极的最大变形量和最大应力增加。
流固耦合基础知识 篇4
流体动力学是研究流体平衡的条件及压强分布、流体运动规律、以及流体与固体之间的相互作用等, 研究结果对分析管道的振动及影响因素有重要意义。本文针对新疆某石化公司的10-K-302C离心式甲烷制冷压缩机自开机以来润滑油管线振动较大的问题, 通过对管内流体流动状态进行模拟分析, 得出了流体耦合前后动力特性的变化及管道振动的原因。
2 双向流固耦合分析原理
流固耦合要遵循质量守恒定律、动量守恒定律、能量守恒, 所以在流固耦合交界面处, 应满足流体域固体应力 (σ) 、位移 (d) 、温度 (T) 、热流量 (q) 等变量的相等或守恒, 即满足下面四个方程:
式中:df、ds分别为液体、固体的边界位移。σf、σs分别为液体、固体应力。
3 流体和管道的计算模型
就10-K-302C离心式甲烷制冷压缩机装置的润滑油管线位移较大现象, 通过分析润滑油耦合前后的动力学特性, 找出流体运动特性, 对寻找该管道振动原因有重要指导作用。出口管道的管路图如图1:选取润滑油在弯管中心轴线处的1、、2、、3、、4、点, 及在出口处5、为观测点。
用Workbench有限元分析软件对管道结构进行建模, 并利用workbench自带的填充功能生成润滑油模型。将实际的润滑油简化为不可压缩、无粘性、平均密度和压力在整个流域中均一的流体。分析时设定管道内润滑油的材料属性:质量864Kg/m3, 动力粘度为4.33cp, 比定压热容为2063J/ (Kg.K) 。建立润滑油模型, 为了更好地描述边界层处参数变化情况, 在网格划分时设定膨胀层。
4 流体耦合前后分析对比
4.1 流体耦合前后的压力和速度对比
润滑油耦合前后的速度流线图2和图3, 可以看出耦合前润滑油的最大流速达到了5.168m/s, 耦合后的最大速度为4.555m/s, 较耦合前变小, 但耦合前后流体速度均在弯头处较大,
耦合前后润滑油与管道接触壁面的压力云图4和图5。绝对压力均在入口处较大, 弯头处较其连接处的直管压力较大。耦合前润滑油壁面的最大绝对压力为772KPa, 最小绝对压力为759.9KPa, 压力波动值为1.58%, 压力波动较小。流固耦合后接触壁面的压力大小和分布与耦合前几乎相同。
4.2 耦合前后流体观测点随时间的变化对比
观察观测点在耦合前后轴向速度和压力随时间变化曲线 (图6、7) 。从图6、7可以看出, 润滑油在耦合前后各个弯头处的测试点的速度相差较大, 同时达到稳定的时间不同;耦合前测试点3、4、达到稳定时间较长, 出口处的稳定速度较入口变大;耦合后2、、3、、4、速度波动较大且要达到稳定的时间比较长;耦合前后弯管处的1、、2、、3、、4、的压力同时达到了最值;耦合后除了初始阶段压力值波动较大外, 稳定时间和稳定值都一样。
4.3 耦合前后流体中心轴线处的压力速度对比
耦合前后润滑油从入口到出口的中心轴线处的压力和速度变化情况如图8和9。可以看出耦合前后的流体压力和速度变化情况基本一致, 且都是在管道弯头处压力和速度变化都达到大。
5 结论
(1) 车间提供的出口管路的测试压力波动范围在660-680KPa之间, 软件模拟分析的流体的压力波动范围在660-672KPa之间, 与利用cfx计算结果的差值为1.17%, 说明了可以用双向流固耦合有限元模拟方法进行管道流体的分析研究。
(2) 在润滑油速度为1.93m/s的情况下流体耦合前后的流体压力变化很小, 说明流体的压力能在流固耦合后变化很小, 对管道产生的影响较小;耦合前, 润滑油最大流速为5.168m/s, 耦合后的最大速度为4.555m/s, 流体速度较耦合前减小了11.86%, 由耦合时能量守恒可知, 流体的这部分动能传递给了管道。
(3) 耦合时流体的速度波动比较大, 特别是观测点2、3、4处;耦合时流体速度达到稳定的时间较耦合前达到稳定的时间较长, 这样对管道的影响也会变得较大。耦合后最大速度比设定速度1.93m/s大了136%, 且在弯头处变化最大, 这样会对管道产生很大的冲击力, 造成管道的振动。
参考文献
[1]谢龙汉.ANSYS CFX流体分析及仿真[M].北京:电子工业出版, 2012, 1-20
[2]赵兴艳.CFD方法在流体机械设计中的应用[J].流体机械, 2000, 28 (3) ;22-25
流固耦合基础知识 篇5
研究人体上呼吸道气流流场特性，对探索气溶胶在人体上呼吸道内的扩散及沉积模式和分析大气环境状况等因素对人体健康的影响有非常重要的意义。近年来，国内外研究人员采用数值仿真的方法对人体上呼吸道气流流场特性进行了大量有意义的研究工作。徐新喜等[1,2,3]通过建立三维的口喉模型，数值模拟了循环呼吸条件下气流组织特性，进行了初步的实验研究;Wang等[4]采用数值方法，模拟了人体上呼吸道内的气流组织形式;Nithiarasu等[5]也对人体上呼吸道中的层流和湍流等形式进行了仿真研究;Wlofgang等[6]对考虑壁面弹性的人体前四级支气管内气流运动进行了数值模拟研究;X.G.Cui等[7]对理想的人体口喉模型内的气流组织进行了精细的大涡模拟，指出喉部射流将和二次涡流在射流边界处产生相互作用，且喉部射流尾部将发生严重变形。但是目前研究也存在不足，一是采用的模型多数为简化模型，与真实人体上呼吸道结构差异较大;二是大多数数值模拟模型采用低Reynolds数k-ω模型[8,9]和标准k-ε模型[10]来描述气体流动情况，缺乏更高精度的仿真计算。同时，在以往的研究中，通常假设呼吸道壁面是刚性的。流固耦合力学的发展和成熟能够克服这一假设存在的缺陷，使数值仿真的理论分析与呼吸道内的实际气流情况更为相符，从而开拓研究呼吸道内气流运动的思路和方法。基于此，本文运用流固耦合力学的方法，对循环呼吸模式下人体口喉模型内气流和呼吸道壁的相互耦合作用进行了数值模拟，系统地分析了流固耦合作用下口喉模型内的气流运动特点。
2 数值仿真模型和方法
2.1 仿真模型
参考ARLA(aerosol research laboratory of alberta)模型[11]和Stapleton K.W模型[12]的尺寸，文中建立和采用如图1所示的口喉仿真模型。
2.2控制方程
2.2.1流体控制方程
文中采用流体湍流模型为Menter’s SST两方程模型，即耦合了LangtryMenter的转捩模型。该模型基于间歇因子和转捩动量厚度雷诺数2个输运方程，前者与SST湍流模型耦合，后者用来捕捉非当地湍流强度的影响，其张量表达式[13,14]为:
间歇因子输运方程:
式中生成Pγ项的表达式为：
其中，S为应变率张量，Flength用以调整转捩区长度，Fonset用以启动Pγ，包含了“当地转捩”判据，当该比值大于1时，Pγ启动即γ开始增长。
转捩动量厚度雷诺数输运方程：
式中的源项定义如下：
2.2.2固体（壁面）结构控制方程
本文避免了传统的刚性壁面假设，而是将口喉模型壁面材料假设为线弹性材料，这样与实际情况更为接近。参照文献[15]，设定壁面的杨氏弹性模量为9 MPa，泊松比为0.4，为均匀各项同性材料，即服从胡克定律的薄壁弹性模型。
气管壁面方程可表达为：
式中：σs是应力张量，ρs是密度，as是加速度。
在流固耦合交界面上，应当满足以下条件：
式中，d为位移，n为边界法向，下标s、f分别表示固体和流体。
2.3 边界条件和计算方法
假设口腔入口气流流量随时间按正弦函数形式周期性地变化，且气流在入口边界是匀速的，方向为入口平面的法向方向，呼吸流量为30 L/min，呼吸频率为15循环/min。支气管出口采用压力出口边界条件，压力为0 N，其余量梯度为0。在一个循环周期内，选取呼吸吸气和呼气状态下的峰值T3和T7的仿真结果分别进行分析，其中T3=1 s,T7=3 s，如图2所示。
3 结果与讨论
3.1 吸气时刻仿真结果分析
如图3所示，循环吸气过程中呼吸流量达到最大值T=1 s时，口喉模型中不同切面内气流轴向速度分布，其中A-A’为冠状面，其余为水平面。图3(a)显示的是口喉模型的正中矢状面(将口喉模型分为对称的左、右两部分的纵切面)的气流流场分布结果。观察发现，气流在口喉模型内的咽部和喉部，形成了2个速度增长点。在喉部，气流速度达到全场最大值12.33 m/s，而在咽部气流发生分离现象，在咽部外壁形成分离区，如图3(d)所示。咽部外壁分离区的产生，将使气溶胶在此部位回旋、滞留，导致其在咽外壁的沉积几率增加。图3(e)、(f)显示，湍流喷射现象出现在了声门位置，接着气流在声门下方的气管内出现流动分离，这主要是声门的形态结构所致。图3(f)显示声门下方近气管内壁气流速度较高，而图3(g)结果表明，与声门距离增加时气管内壁气流速度仍高于外壁，但气流速度差，并呈现出逐渐减小的趋势。气流速度高的部位，其内部气溶胶颗粒的惯性碰撞作用多，因此将造成多数气溶胶颗粒在气管内壁沉积。
3.2 呼气时刻仿真结果分析
如图4所示，为循环呼气过程中口喉模型内呼吸流量达到最大值T=1 s时刻不同切面内的轴向速度分布图，其中图4(a)显示的是口喉模型的正中矢状面(将口喉模型分为对称的左、右两部分的纵切面)的气流流场分布结果。从图中可以看出，在y=0截面内，在咽部和喉部由于受到截面收缩和气流运动方向改变的影响，分别在咽部前壁、喉部后壁和声门后部形成高速区，气流速度的最大值为13 m/s，如图4(c)～(e)所示。气流速度的增加，容易造成气溶胶颗粒在咽部和喉部的沉积。图4(b)显示，与吸气阶段不同，在呼气阶段在口腔内部发生流动分离现象，形成分离区，分离区的产生使得气溶胶颗粒在口腔中的沉积几率增加。在呼气阶段，由于呼气气流与吸气阶段相比气流流向发生改变，在咽部和喉部未出现气流分离现象，而且气管内气流流速分布比较均匀，在气管不同截面内速度梯度不明显，如图4(f)、(g)所示。如图5所示，为吸气阶段和呼气阶段各个时刻流场内的最大速度曲线。从图中可以看出，呼气阶段气流最大速度略高于吸气阶段，对气道的冲击力也将更强一些。分析原因，一方面是由于气流流向不同，另一方面则与气道的膨胀与收缩有关。
4 结论
流固耦合基础知识 篇6
高杰、郑群等[3]认为间隙流动的产生机理可以分为“攻角驱动型”和“压力驱动型”,而“压力驱动型”则是造成间隙泄漏的主要因素。不同大小的间隙,所产生间隙涡的强度、位置以及对主流的影响都不相同[4,5]。目前普遍认为,较小的叶尖间隙能够改善涡轮性能并能降低噪声; 过大的叶尖间隙会使涡轮性能下降,并可能引起旋转失速、喘振等非定常流动现 象[6]。J. Tallman和B. Lakshminarayana[7]认为泄漏涡主要是由通过间隙的流体所引起的,如果减少叶尖间隙的高度,泄漏涡将变小,间隙和泄漏涡的气动热力损失也减小。王大磊,朴英[8]通过对某跨声速高压涡轮流场模拟研究,发现间隙高度每增加1% 相对叶高,间隙泄漏流量占总流量的比例约增加2. 1% ,造成性能损失线性增加。
间隙大小对涡轮性能有着很大影响,但在涡轮高温、高压、高转速工作条件下,叶顶间隙流动也会受到一系列条件的影响[9]。岂兴明和朴英等[10]利用有限元分析软件发现涡轮叶片位移对叶顶间隙变化影响不大,而温度和转速才是叶片及轮盘径向位移的主要因素。流固热耦合方法是研究涡轮叶片温度场和进行静力分析的重要途径[11,12]。吴正人和王松岭等[13]对离心风机使用流固耦合方法进行研究,发现叶轮叶片在流固耦合条件与非流固耦合条件下所受应力具有很大区别。本文分析了不同间隙高度对涡轮性能的影响,并研究了在流、固、热耦合条件下动叶叶片及轮盘径向变形量对间隙尺寸的影响。
1计算模型以及方法
研究对象为一个3级涡轮,入口处存在入口导叶igv,在单流道中,进口导叶排由一个大叶片另加三个小叶片组成,由燃烧室流入的高温高压燃气在经过导叶之后依次经过动叶R1,静叶S1,动叶R2, 静叶S2,动叶R3流出通道。
图1所示为此研究对象的网格,在网格划分时对叶顶和叶根处进行加密,网格总数约为260 × 104。涡轮中三排动叶转速为8 001 r/min,在涡轮叶顶间隙设置时动叶和静叶间隙保持一致,叶顶间隙大小分别取0 mm、0. 2 mm、0. 4 mm、0. 6 mm、0. 8 mm、1. 0 mm、2. 0 mm。
对于涡轮在流固耦合下的静力分析,使用ANSYS中的WORK BENCH平台进行流体、固体和热耦合计算。用CFX进行气动计算,计算时入口边界条件设定总温、总压,并根据进口导叶的气流角给定速度方向,出口处设置静压参数使用k-epsilon湍流模型进行求解。此外用热模块进行热计算,用结构模块进行受力分析,在workbench中用热模块进行计算时,使用传热系数相较于使用热流量能更好贴合实验结果[14],故而在热应力分析时使用传热系数而非热流量参数。
对涡轮动叶叶片、轮盘以及机匣都使用流固耦合方法进行结构特性分析。叶片结构网格总数为37 × 104左右,轮盘网格数约为91 × 104,机匣网格数为11. 4 × 104。图2为动叶轮盘三维实体图。图2从左至右分别是动叶R1、R2、R3所对应轮盘。图3为机匣三维截图。
2计算结果分析
从气动性能及结构特性两方面分析涡轮的特性,气动性能的分析可以得出不同叶顶间隙高度对涡轮性能影响的变化,对结构特性的简要分析可以看到实际运转条件下,流、固、热多重因素作用下叶片、轮盘以及机匣径向变形量的变化。
2.1气动性能分析
当叶顶间隙高度增加时,对叶顶区域的流场分布影响变大。图4是不同叶尖间隙下出口处的绝对Ma数分布,由图4上可以看出,轮毂处的Ma数小于靠近机匣处,没有间隙时Ma数的分布与有间隙时有明显不同,随着叶顶间隙的增大,叶顶处Ma数分布越加复杂,意味着流场越加复杂。如图中红圈所示,随着叶顶间隙增大,由于受到更强烈的间隙泄漏流的影响,在出口处靠近叶顶区域会形成低速区。 随着间隙增大,低速区域扩大。尤其在叶顶间隙达到2 mm时,此低速区有明显的向相邻叶片通道移动的趋势。
为了说明叶顶间隙变化对涡轮性能的影响,对涡轮各动叶排进出口总压的损失进行了对比。因此定义了总压损失系数 ξ[15],其公式如下所示:
其中Pt1和Pt2分别表示进口和出口的平均总压,Ps2表示涡轮出口的平均静压,在此都使用质量平均。
由于动叶转子叶顶泄漏流对涡轮总压损失的影响远大于静叶[17],因此在图5中给出了各动叶排在不同叶顶间隙高度时总压损失系数沿叶高方向的分布。由于叶顶间隙以及轮毂和机匣端壁面的影响, 各动叶叶排的总压损失系数在中间叶高处达到最小,由中间叶高处向叶顶及叶根处逐渐增大,尤其在叶顶处总压损失系数达到最大。对于任一级动叶, 随之间隙高度的增加,总压损失系数随之增大。从而影响涡轮的总体性能。
叶顶间隙泄漏流与主流掺混产生泄漏涡,这是影响涡轮叶片性能损失的最主要原因。本文针对R1叶片排,重点研究了不同间隙高度下泄漏涡的变化。图6是R1叶片通道三维Ma数分布云图,用9组垂直于叶片弦长的切面来显示叶片通道Ma数的分布。在叶顶区域靠近吸力面部分,产生了一个从叶片前缘至尾缘面积逐渐增大的低速区域,而此区域中心处正是泄漏涡的涡核,低速区域的发展代表了泄漏涡的流动轨迹。
对图7中不同间隙高度下所产生的泄漏涡进行对比。在叶顶间隙较小时,如图7( a) 及7( b) ,在靠近叶顶区域,从叶片前缘流入紧贴吸力面的流体形成泄漏涡的涡核,之后通过叶顶间隙的泄漏流卷入到此涡中。而随着叶顶间隙增大,当在1 mm叶顶间隙高度时,如图7( c) ,通过靠近叶片前缘处的叶顶间隙的流体与贴近吸力面的来流一起形成泄漏涡核。而当间隙高度为2 mm时,如图7( d) 所示,泄漏涡的涡核只由通过上游弦长处间隙泄漏流构成。
如图7所示,随着间隙高度的增加,在叶片通道中靠近吸力面的低速区域随之扩大,间隙泄漏流与主流掺混卷起的泄漏涡的大小及强度亦随之增加,从而对主流通道的堵塞加剧。在叶顶间隙达到2 mm时, 泄漏涡甚至将扩散至整个叶片通道。这与之前得出的总压损失随间隙增大而增大的结论一致。
2.2动叶叶片结构特性分析
使用流固热耦合方法对叶片结构特性进行分析时发现,在不同的叶顶间隙下,叶片的变形量差别不足0. 01 mm,故而只取叶顶间隙为1 mm的叶片做结构特性分析。由于动叶在极高转速的工作条件下受到非常大的离心作用力,因而动叶变形量远高于静叶,故而本文只对动叶的结构特性做出分析。
涡轮叶片的主要材料参数为: 密度8 480 kg / m3,叶片杨氏模量1. 255 × 1011,泊松比为0. 3,热膨胀系数13. 53 × 10- 6℃- 1,导热率为16. 79 W / ( m·℃ ) 。在流固耦合计算时,对叶片根部施加固定约束,即限制其位移方向为沿叶高方向。对叶片施加温度、压力载荷,加载的压力、温度量来源于流场计算结果。
图8为涡轮三排动叶在流固热耦合条件下的径向变形量。随着叶高的增加,叶片的径向变形量随着增加,到达叶尖处变形量达到最大。涡轮叶片的最大径向变形量出现在叶高最大的最后一级动叶R3的叶尖处,约为1. 295 mm。此变形量对叶顶间隙的影响是较大的。
如图9所示,由于高温高压燃气对叶片的冲击, 涡轮每级动叶的轴向变形量由前缘至尾缘逐渐增加,也据此可知涡轮叶片的轴向变形主要是沿流道方向。由于涡轮流道内的压力及温度逐级降低,因此最小轴向变形量出现在R3叶片的前缘处。同时又由于R3叶片相较前两排动叶具有更大的叶高, 因此具有更大的离心力,故而最大轴向变形量出现在R3叶片叶顶尾缘处。由于机匣在沿流道方向具有正斜率,因此叶片的轴向变形将造成叶顶间隙增大。而相较于图8,涡轮的轴向变形量小于径向变形量。
图10为没有热载荷时动叶叶片的径向变形分布图。当没有热载荷时,动叶叶片径向变形量最大值为0. 254 mm,相对于有热载荷时变形量减少了80. 4% 。由此可见温度是影响叶片径向变形的最重要因素。
2.3涡轮轮盘结构特性分析
涡轮轮盘的材料参数为: 密度8 240 kg /m3,叶片杨氏模量1. 999 × 1011,泊松比为0. 3,热膨胀系数16 × 10- 6℃- 1,导热率为22 W/( m·℃) 。在流固耦合计算时,对内表面即最小半径处施加固定约束, 即限制其位移方向为外法线方向。对外表面施加温度、压力载荷,此外在表面即三个动叶的轮毂加载的压力、温度量来源于流场计算结果。
图11、图12分别为动叶轮盘在有热载荷和没有热载荷耦合条件下的径向变相量。没有温度载荷时径向变形量相对有温度载荷时减少86. 78% ,因此温度依然是影响涡轮轮盘变形的最主要原因。
如图11、图12所示,涡轮轮盘的最大径向变形出现在R1外表面,即靠近气动流道轮毂处,最大变形量约为3. 77 mm。涡轮轮盘的材料相对于叶片具有更大的密度,此外轮盘的体积远大于叶片,因而轮盘具有更大的质量及受到更高的离心力,故而对比图8,轮盘的径向变形量远大于叶片的径向变形量。 涡轮的变形是叶顶间隙变化的主要原因。
2.4涡轮机匣结构特性分析
涡轮的机匣的材料与叶片相同。对6排叶片对应的机匣分别施加气动计算得到的压力及温度载荷,图13为使用流固耦合方法下的机匣径向变形量,本文只对动叶排的机匣变形进行研究。由于在机匣两端施加固定法兰约束,因此最大变形量出现在中间R2机匣处,但由于温度及压力沿流向方向逐渐降低,因此最小变形出现在R3机匣处。
涡轮机匣的变形会使叶顶间隙变大,而叶片及轮盘的径向变形会使间隙减小。在实际中,叶片、轮盘、机匣共同的变形量决定了叶顶间隙高度的改变。 结合图8、图11、图13,得出各动叶排叶顶间隙高度变化值,如图14。对于3排动叶,在叶片、轮盘、机匣的共同作用下,叶顶间隙的净变化量皆为负值,间隙减小。叶顶间隙高度净变化量由R1至R3递减。
在实际运转环境中,涡轮叶片及轮盘在高温高压高转速的情况下会产生比较大的径向变形,径向变形量会导致叶尖间隙变小,而同时叶片也会产生轴向变形,从而造成了叶尖间隙控制的复杂性。在设计涡轮的时候一定要对叶顶间隙进行详尽的考虑及分析。在实际中为了确保涡轮的性能以及防止安全问题的出现,也应采取有效的叶顶间隙控制技术将叶顶间隙高度控制在合适的范围。
在涡轮设计时,要考虑叶顶间隙对性能的影响, 叶顶间隙不应过大。而叶片和轮盘的径向变形则制约着叶顶间隙的最小值,尤其是轮盘的径向变形量对间隙高度的影响更大。设计涡轮时需要考虑叶顶间隙的影响,尤其是叶片变形的影响,但同时需要兼顾气动性能的变化。
3结论
本文对叶顶间隙流动气动及涡轮叶片及轮盘的流固耦合分析,结论如下:
( 1) 涡轮叶顶间隙高度对涡轮的性能影响较大,随着叶顶间隙高度的增加,涡轮总压损失增加, 涡轮将压力转化为动能的能力降低,涡轮性能变差。
( 2) 随着叶顶间隙高度的增加,泄漏涡的大小及强度随之增加,尤其在大间隙情况下,泄漏涡甚至会扩散至整个叶片通道。
( 3) 涡轮动叶叶片的径向变形量要远大于轴向变形量,叶片的径向变形量对叶顶间隙高度有更显著的影响。由于轮盘的质量远大于叶片,故而在旋转中产生更大的离心力,所以轮盘的径向变形量要远大于叶片的径向变形量。
( 4) 温度是影响涡轮叶片和轮盘径向变形的最主要因素。相较有热载荷工况,无热载荷工况下涡轮动叶叶片及轮盘的径向变形量会减少80% 以上。
流固耦合基础知识 篇7
荆江河段航道整治工程是我国内河水运“十二五”期重点建设项目，受长江水位季节性涨落及三峡水库调节影响，水位变幅大，加之属于多雨地区，高滩岸坡渗流场、应力场产生周期型变化，极易产生较大变形甚至失稳，破坏水运航道，威胁过往船只的安全[7]。本文以荆江河段中洲子高滩岸坡为背景，结合该河段实际水位变化情况，考虑岸坡的渗流场与应力场完全耦合作用，对高滩岸坡进行了非线性有限元分析。应用Geo-Studio软件的Seep/W模块和Sigma/W模块，探究长江水位升降过程中岸坡位移场和应力场的变化规律，为长江荆江河段岸坡稳定性治理提供参考。
1 计算理论
流固耦合的本质是在实际渗流过程中，由于孔隙流体压力的变化引起土骨架有效应力的变化，直接导致土体的渗透率和孔隙度的变化;同时土体的渗透率和孔隙度的变化会反过来影响孔隙流体的流动和压力的分布。为准确模拟流固耦合，必须将渗流力学和土力学结合起来，全面考虑流体渗流与土体变形之间的耦合过程，建立合理的耦合数学模型。
1.1 土体平衡方程
1993年，著名学者Fredlund和Rahardjo提出对于非饱和土介质的增量应变-应力关系，对于二维问题，可以表示为下列应力-应变形式。
对于二维的平面应变问题，εz=0。其增量应力-应变关系可以写成下式:
式中:σ为法向应力，τ为剪切应力，ε为法向应变，γ为工程剪应变，υ为泊松比，H为基质吸力(uauw)有关的结构的非饱和土模量，E为弹性模量，D为排水条件下有效应力-应变本构矩阵。mH可按式(3)计算。
当任何时间空气压力为大气压力时，ua=0。
1.2 渗流方程
单元体积土体的二维渗流可以用达西定律来表示渗流方程。
式(5)中，kx、ky分别表示x和y方向的渗透系数;vw为渗流速度;γw为单位水的重量;θw为体积含水量;t为时间。
1.3 耦合分析的有限元公式
在完全耦合分析中，平衡方程和渗流方程同时被求解。有限元方程可以写成以下矩阵形式:
式(6)中:
式中:B为梯度矩阵;K为刚度矩阵;LD为耦合矩阵;N为形函数行矢量;Δδ为增量位移矢量;Δuw为增量孔隙水压力矢量。
对于使用虚功原理的有限元方程，渗流方程可以表示为孔隙水压力和体积应变的形式。将虚的孔隙水压力作用于渗流方程，并在整个体积上积分，使用有限元近似法，从时间t积分到时间t+Δt，通过时间差分技术，得到一个只包含增量孔隙水压力的方程:
2 有限元计算模型
2.1 工程概况
中洲子高滩岸坡位于长江中游荆江河段莱家铺水道，岸坡平均走向为140°。岸坡以上为中州子边滩，滩顶较平坦，滩顶高程33.0～35.8 m。水面线以上坡高11～12 m，后缘因崩岸呈多陡坎状，坎高8～10 m，整体坡度一般50°～70°，岸坡上部黏性土，中下部以下为粉细砂，抗冲能力差，岸坡近年来崩坍现象极严重。工程区的地层主要为第四系全新统冲积层，可细分为(2)1灰黄色粉质粘土、(3)1松散～稍密粉细砂、(3)2中密粉细砂、(3)密实粉细砂四层。
2.2 计算模型
本文采用加拿大Geo-Slope公司开发的岩土数值分析软件Geo-Studio，运用Seep/W模块和Sigma/W模块对水位升降下非饱和土坡渗流场和应力场进行完全耦合分析。首先在Seep/W模块中对边坡的渗流场进行稳态分析，得到初始条件下的孔隙水压力，并将其计算结果导入Sigma/W模块，进行初始应力场的模拟;再把初始孔隙水压力和初始应力导入下一步的耦合分析中，不断重复上述过程，直至计算到最终时段末的渗流场与应力场。
2.2.1 模型网格确定
计算断面位于长江荆江河段中洲子高滩守护工程上游200 m处，岸坡顶面高程34.5 m，河底高程约13.0 m，岸坡高度约21.0 m。岸坡断面共有三个土层:第一层土为粉质黏土，厚3.5 m;第二层土为松散～稍密粉细砂，厚17 m;第三层土为中密粉细砂。计算模型的几何长度取100 m，高54.5 m。采用四边形单元类型，二维有限元网格共剖分为2 034个单元，2 271个结点，计算模型网格如图1所示。
2.2.2 边界条件
利用Seep/W模块模拟渗流场时，模型两侧地下水位以上边界按零流量边界处理，地下水位以下为给定水头边界，模型左侧边界根据实测地下水位资料施加水头边界模拟地下水位，模型右侧坡面水位值根据长江水位的变化而定，模型底边取零流量边界。利用Sigma/W模块模拟应力场时，边坡底部采用水平和竖直方向的完全约束，左右采用水平约束，边坡涉水坡面为坡前水位静压力随时间的变化函数，其大小与坡前水位及高程相关。
2.3 计算参数
根据室内模型试验提供的土体力学参数，并参照工程地质手册等相关资料，本工程实例土体物理力学计算参数指标如表1。
2.4 计算工况
本文基于实测水位数据，研究水位升降下岸坡的应力应变规律。2014年长江中洲子河段水位变化过程如图2所示。计算工况考虑2014年4月1日至2014年12月28日，以1天为1个时步，共分271个时步进行。
3 计算结果与分析
为充分探讨水位升降下荆江河段高滩岸坡的变化规律，本文对水位上升过程、水位下降过程分别进行了岸坡的位移场、应力场分析。水位上升过程中取水位上升初期水位24.0 m和水位上升末期水位33.64 m的计算结果进行分析;水位下降过程中取水位下降中水位29.91 m和水位下降末期水位23.68 m的计算结果进行分析。
3.1 位移场分析
3.1.1 水位上升过程岸坡水平位移
图3为水位上升过程中不同水位下岸坡的水平位移等值线图。由图可知，当河水位为24.0 m时，位于高程25 m的坡面前缘渗流溢出点所在区域的位移最大，水平位移最大值为100 mm;坡顶位移在20～40 mm之间;当河水位上升至33.64 m时，坡内最大水平位移达到160 mm，位于坡面约25 m高程处，此时坡顶位移约为90 mm。
与初期低水位相比，末期高水位情况下，坡内最大水平位移以及坡顶位移均有大幅的增大，增大约60 mm，坡内大部分区域土的水平位移较低水位有明显的减小。随着浸润线上升，位移0 mm等值线也随之上升，坡前水压力挤压土体，土体颗粒明显向坡内方向移动。
3.1.2 水位下降过程岸坡水平位移场分析
图4为水位下降过程中不同水位下岸坡的水平位移等值线。由图可知，当河水位29.91 m时，位于高程25 m的坡面前缘产生最大水平位移180 mm，比高水位33.64 m时的位移增大了20 mm;坡顶位移约为90 mm;当河水位下降至低水位23.68 m时，坡内最大水平位移有所增大，增大为220 mm，而坡顶位移几乎不变。由此可见，水位下降过程中，浸润线不断下降，产生坡内指向坡外的渗流，位移0 mm等值线逐渐消失，20 mm等值线范围不断扩大，土体颗粒明显向坡外方向移动。
3.1.3 岸坡水平位移场随时间的变化
为揭示岸坡水平位移随水位变化而变化的规律，取特征点坡顶A点和坡脚B点的水平位移进行分析(A点水平位置为40 m，高程约33.5 m;B点水平位置为90 m，高程约14 m，略高于河底高程)，如图5、图6所示。
由于实时监测水位的波动性，水平位移计算值呈现出上下波动的结果。从图中可以看出，坡顶的水平位移随着水位的上升而增大，随着水位的下降略有减小。高水位33.64 m时，坡顶的水平位移达到最大，为100.54 mm;随着水位不断下降，坡顶位移也有所减小，水位降至23.68 m，此时水平位移为91.97 mm，减少了5.69 mm。
坡脚B点水平位移的变化随水位升降变化明显。河水位上升时，坡脚水平位移减小;河水位下降时，坡脚水平位移增大。随着水位的缓慢上升，地下水位与坡外水位基本同步上升，浸润线比较平缓。当水位达到最大值33.64 m时，坡脚水平位移为7.50 mm;当水位下降至23.68 m时，水平位移为32.34 mm，比高水位时增大了24.84 mm。
综上所述，随着水位上升，岸坡坡顶处的水平位移增大;而坡脚处的水平位移随着水位上升而减小。
3.2 应力场分析
3.2.1 水位上升过程岸坡应力场分析
图7为水位上升过程中不同水位下坡内的最大剪应力的分布情况。从图中可以看出，土体的最大剪应力的方向与坡面近似平行;随着水位的上升，坡顶及坡体前缘的剪应力等值线向下移动，坡体内土体所受剪应力逐渐减小。低水位时坡底最大剪应力值为160 k Pa，而高水位下坡底最大剪应力值为140 k Pa。
3.2.2 水位下降过程岸坡应力场分析
图8为水位下降过程中不同水位下坡内的剪应力的分布情况。从图中可以看出，土体的最大剪应力的方向与坡面近似平行;随着水位的下降，坡顶及坡体前缘的剪应力等值线向上移动，坡体内土体所受剪应力逐渐增加。高水位时坡底最大剪应力值为140 kPa，而低水位下坡底最大剪应力值为160 kPa。图8(b)中剪应力20 k Pa等值线区域明显比图8(a)减小，同时40 kPa等值线范围随着水位的下降向坡面扩展，坡内剪应力值增大，当坡内剪应力大于土体抗剪强度时，则可能产生岸坡失稳。上述结果表明，低水时岸坡稳定处于较危险状态。
3.2.3 岸坡特征点应力随时间的变化分析
图9分别给出了高程25 m、高程15 m下水平位置为10 m、30 m、50 m处岸坡最大剪应力随时间的变化过程。图9(a)中靠近坡内的两点(x=10 m、x=30 m)最大剪应力变化规律相似。水位上升时，剪应力值减小，高水位33.64 m时，最大剪应力分别为30.45 kPa和27.76 k Pa，比低水位时减小了约20kPa;水位下降过程，剪应力值增大，水位降至低水位23.68 m时，最大剪应力分别为50.67 k Pa和49.52kPa;靠近坡体前缘的点(x=50 m)剪应力值变化幅度较小，主要在15～20 kPa之间波动。
图9(b)中当水位上升时，剪应力减小，最大水位33.64 m时，最大剪应力分别为55.7 kPa、47.9kPa和37.7 k Pa;水位下降时，剪应力上升，水位23.68 m时，最大剪应力分别为80.0 kPa、74.2 kPa和54.8 kPa，比高水位时增加了约25 kPa;从图中可以看出坡内土体所受剪应力比坡体前缘剪应力大，且越靠近坡体前缘，剪应力越小，受水位变化影响越小。
4 结论
本文以非饱和土的流固耦合理论研究荆江河段中洲子高滩岸坡在水位变动下的位移场和应力场的变化规律。根据研究结果得到以下几点结论。
(1)在全年水位变化过程中，坡顶位移随着水位的上升逐渐增大，随着水位下降，坡顶位移变化不大;坡脚位移随水位升降的变化较为明显，随着水位的上升坡脚水平位移减小，随着水位的下降坡脚水平位移逐渐增大。
(2)岸坡内最大水平位移发生在水位下降末期，岸坡最大水平位移为220 mm，位于高程约25 m坡面处。
(3)岸坡内剪应力的大小随着水位的升降而变化，随着坡前水位上升，坡体剪应力减小;随着坡前水位降落，坡内最大剪应力逐步增加，水位下降末期，坡内剪应力最大，此时最易发生边坡失稳。
参考文献
[1]陈晓平,茜平一,梁志松,等.非均质土坝稳定性的渗流场和应力场耦合分析.岩土力学,2004;25(06):860—864Chen Xiaoping,Qian Pingyi,Liang Zhisong,et al.Coupling analysis of heterogeneous earth dam stability of seepage field and stress field.Rock and Soil Mechanics,2004;25(06):860—864
[2] 杨相如.坡前水位骤变对土坡特性的影响.杭州:浙江大学,2008Yang Xiangru.Study on the characteristics of earth slope with sudden change of water level in front of slope.Hangzhou:Zhejiang University,2008
[3] 张文杰,詹良通,凌道盛,等.水位升降对库区非饱和土质岸坡稳定性的影响.浙江大学学报(工学版),2006;40(08):1365—1370,1428Zhang Wenjie,Zhan Liangtong,Ling Daosheng,et al.Influence of reservoir water level fluctuations on stability of unsaturated soil banks.Journal of Zhejiang University(Engineering Science),2006;40(08):1365—1370,1428
[4] 汪斌,唐辉明,朱杰兵,等.考虑流固耦合作用的库岸滑坡变形失稳机制.岩石力学与工程学报,2007;26(S2):4484—4489Wang Bin,Tang Huiming,Zhu Jiebing,et al.Deformation and failure mechanisms of reservoir landside considering fluid-solid coupling effect.Chinese Journal of Rock Mechanics and Engineering,2007;26(S2):4484—4489
[5] 万文.水位循环升降下库岸边坡渗流与应力耦合分析及稳定性研究.南昌:南昌大学,2012Wan Wen.The analysis about seepage and stress coupling of the bank slope under the cycle of water lifting and the study of stability.Nanchang:Nanchang University,2012
[6] 李智慧.坡前水位降落对边坡稳定性影响模型试验与数值分析研究.北京:中国地质大学,2012Li Zhihui.Study on the model test and numerical analysis of slope stability caused by water drawdown.Beijing:China University of Geosciences,2012