基于非结构化网格和浸入边界的流固耦合数值模拟 图文.doc

基于非结构化网格和浸入边界的流固耦合数值模拟_图文万方数据万方数据明平剑,等:基于非结构化网格和浸入边界的流同耦含数值模拟(a原始浸入式边界方法图l不同网格方法步的值;甜;和嵋均为上一时间步的速度值;D为单元体积。浸入式边界对流体的影响采用vOs法计算。质量力按下式计算:Ⅳ一“‘,:2i1,一1,‘2i其中”’和矿为速度二次修正值,可通过流体和浸入刚体边界速度线性插值得到‘11】:甜,=(1一口甜7+劬?“矿=(1一倪/+叫+1(b非结构化网格浸入式边界方法图2与浸入边界相交单元的固体体积分数325删IⅨ=Pl珀P2IⅨ挖肿2●_—--—-__l-●-—・_●._——●Ic2一几●—・———●・———_●_____●■_-_・—●lCl钾2I鼢2I啪姗俐挖图3单元内部形成链表示意图(82・2・2计算固体体积分数其中口是单元体中的固体体积分数(单元中固体体积与单元总体积的比值,。将材’和’,’代入(5、(6得,~“j一甜六钳tF,..匕一’,‘钳‘F(10压力速度耦合采用sIMPLEc算法,㈣。。定义了一个数据结构,形成链表来表示边界,其方向是使得固体区域始终在边界的左手侧。为了计算每个单元内部固体体积分数,定义了指针型结构数组,在每个单元内形成一个多边形的链表,多边形的边界具有逆时针方向,如图2阴影部分多边形,图3给出了单元ICl和IC2内的链表,图中Q岬表示网格节点,IBP为浸入边界点,IBCP表示浸入边界与网格面交点。将多边形剖分为一系列三角形计算其面积和,从而确定每个单元内固体体积分数。下面较详细的介绍多边形的形成。从图2中完成IBPl操作后,单元IC2内的链表包含节点GNP6.->IBCPl。完成IBP2后,链表为GNP6一>IBCPl一>IBP2,。找到交点后,IC2内链表增加节点IBCP2,同时确定位于固体区域内的网格节点GNP2。此时单元IC2的链万方数据万方数据明平剑,等:基于非结构化网格和浸入边界的流固耦合数值模拟图7浸入边界前沿示意图RP:型(13弛朋=坠(14fd其中‰圆柱运动最大速度,d圆柱直径厂振荡频率。圆柱直线振荡可以用简谐振动描述:x(f=一以sin(2妒(15其中x(f为圆柱中心位置,以为振荡幅值。由位移对时间求导可以得到圆柱运动速度327计算得到的单个周期内阻力系数变化与文献【14】的结果对比见图8。图9给出了不同相位角时圆柱周围的涡量分布与文献【ll】中数值计算结果的对比。从图中可以看出本文计算结果与文献中cartesian网格下计算结果一致,表明本文方法可以很好地模拟振荡圆柱周围流场。图lO分别给出了180、210、330三个相位角情形,不同空问位置的水平和竖直方向速度分布与文献【14】中实验数据的对比。从图中可以看出,除局部位置有些偏差外,计算结果与实验值吻合较好,如图lO(b右图中,数值计算结果具有较好的对称性,(乒0,∥0侧与实验值吻合较好,而下半部分与实验值不一致,原因可能是实验测量引起的误差。%(,=—2碱,cos(2妒(:‰=2彤k(17计算区域为30出20d。文献【13】采用了同样的计算域研究了类似问题,坐标原点位于计算区域中心,d=,流体密度p=,流体动力黏度∥=・s。给定RP=loo和xC=5,由此可以得到‰=,,=2Hz,4=》5时,阻力系数时间历程阀门是管路系统中重要元件之一,对管路系统有较重要的影响,如文献[15】利用商业软件Fluent研究了阀门的流场特性。本文对二维管道内的阀门动态工作过程进行模拟,各管段长度如图ll所示,,阀门为椭圆形,,长轴方向与水平方向夹角为450。,|jI:。先对阀门静止时进行了模拟。图12给出了阀门后的瞬时涡量分布图,横纵坐标分别为舶和m。图13给出了阀门周围速度矢量分布和固体体积分数分布。从图中可以看到,利用本文方法可以计算出阀门的流动特征,尽管由于阀门绕流在管道内形成了非常复杂流动。图14给出了阀门由图11所示位置顺时针转过90度角的动态过程。从图中可以明显地看到阀门引起的漩涡脱落过程,表明本文方法可以较好的模拟复杂结构运动的流固耦合问题。—(左为文献【l】的结果,右为本文计算结果——o.∞5On0∞