复合材料螺旋桨流固耦合数值模拟

复合材料螺旋桨流固耦合数值模拟

大多数水面舰艘的螺母由锰-ni-铝-铜或铜-锌-银铜制成。尽管黄金材料具有承受强度高、可靠性好等优点，但很难将它们加工成带有旋转噪音的螺母。这是因为振动会产生噪音。纤维增强材料具有比较强的强度、抗腐蚀性好、强度好、耐腐蚀性好、材料设计性强等优点。因此，这种材料可以应用于船舶和海洋推进器系统中的螺钉的完整性。与金属螺旋管不同，金属螺旋管在水动力作用下会发生很大变形，其水动力性能将不可避免地改变。因此，为了获得复合螺旋的水动力性能，调整了螺钉流-固耦合算法是基础。传统的螺旋桨理论的设计和计算基于势流理论，没有充分考虑粘性的影响，也没有考虑旋转程度。因此，在实际情况下，无法正确预测桨叶边界层、舰尾流场结构以及叶尖断裂的形成等水流特征。对于基于rans方程的粘性流场，l-hj等人使用升力法和九节点退化地壳单元的联合算法来解决不同工作状态下螺钉的水动力性能。对于研究螺钉流-固耦合算法，lihj等人采用了提高性能法和九节点退化地壳单元的联合算法，以求解复合螺钉的水动力性能。对于年轻yl来说，面元法和软件abaqus的联合螺旋桨流-固耦合计算方法，但这些方法是基于势流理论的螺钉水动力计算。本文首先通过不同湍流模型和试验,验证求解螺旋桨的敞水性能,再将粘性流场计算软件和有限元耦合起来,推导变形后螺旋桨的几何参数,实现螺旋桨流-固耦合的稳态求解算法.本文计算采用实验测得的复合材料板拉伸模量和泊松比,将复合材料考虑为各向同性,没有考虑复合材料的铺层结构.1数学模型1.1静态压力密度法考虑螺旋桨在粘性湍流中旋转,其连续性方程和动量方程可表示为:连续方程:动量方程:式中:p是静态压力;μ是湍流粘度;ρ是液体密度;为Reynolds应力.RANS方程虽然不用求解流场中的瞬时量,但是方程却引入了新的未知量-雷诺应力,这时的方程解不封闭.要封闭求解方程,就必须引入新的方程,这些方程通过各种湍流模型来定义.1.2空车湍流水动力特性模型为使RANS方程封闭,将雷诺应力用低阶或时均量表达,即湍流模型.封闭RANS方程主要有REYNOLDS应力方程模型(RSM)与涡粘模型2种.涡粘模型以湍流各向同性为基础,认为雷诺应力和时均速度呈线性关系,该类模型求解简单,计算容易收敛.涡粘模型应用较广泛的有κ-ε,RANGκ-ε模型和SSTκ-ω模型.而RSM模型以湍流各向异性为基础,考虑雷诺应力的对流和扩散作用,直接寻找雷诺应力的输运方程.本文通过比较分析几种常用湍流控模型,通过求解时间和求解精度方面的对比以及考虑到流固耦合本身对计算机性能的要求等方面进行综合考虑.螺旋桨水动力性能计算模型选择SSTκ-ω湍流模型.该模型是利用混合函数将κ-ε和κ-ω2方程相结合而构建的湍流模型,在近壁区采用κ-ω方程,其他区域则采用κ-ε方程以获得湍流粘性作用,考虑了κ-ω方程近壁区模拟时的有效性及远场区无法准确模拟的不足.最终选择SSTκ-ω模型求解螺旋桨的流固耦合特性.1.3数值计算及求解计算螺旋桨流场的控制方程是一系列非线性偏微分数学物理方程,需借助数值方法对其进行求解.本文利用Ansys/cfx软件完成该数值计算.收敛判据设定为0.0001.计算区域分为内、外2个流域,在内流场建立一个固定于螺旋桨的旋转坐标系,采用MRF坐标系模型对螺旋桨周围的旋转流场进行计算,外围流场则在绝对静止坐标系下进行求解.将进口边界设置为速度入口条件,给定均匀来流的速度分量;出口边界给定静压分布,外边界设为开放界面;考虑粘性影响,螺旋桨表面定义为不可滑移壁面.2实体螺旋模型的建立文中采用DTMB4119螺旋桨,它是一种无侧斜无后倾分布的三叶螺旋桨,被ITTC选为验证数值方法预报精度的标准模型.桨叶直径为0.3048m,盘面比为0.6,桨叶剖面为NACA-66(mod)型,毂径比为0.2,螺距比(0.7R)为1.084.文中采用Excel编制计算过程文件,计算出螺旋桨叶面和叶背各个半径处的型值点,将计算出的型值点转换为文本文件,然后导入SolidWorks软件中,建立三维实体螺旋桨模型.螺旋桨几何模型和坐标如图1所示.坐标轴的定义为:x轴与螺旋桨的旋转轴一致,指向下游;y轴与桨叶参考线一致;z轴符合右手法则.计算采用全尺寸模型,转速为n=10r/s和20r/s,通过改变来流速度来实现不同的进速系数.网格质量直接决定计算的收敛性、效率和精度,因此,应根据流场中各物理量的分布特点对计算域进行合理的网格划分.将求解域分成旋转区和静止区2个区,采用非结构网格划分方法,2个区采用CFX的GGI方式连接.即首先在螺旋桨表面生成三角形网格,再通过值要控制在30~300,通过多次试算确定表面网格大小和桨叶边界层的过渡方式和层数.螺旋桨表面网格和壁面棱柱过渡层网格如图2所示.总网格数为120×104.3有限元网格划分通过流体计算软件求出螺旋桨的水动力载荷以后通过Ansys-cfx软件指定流固耦合界面将流体压力通过表面效应单元的方式传递给有限元单元.因流体计算和结构计算采用的是非同种单元类型,为保证求解精度和数据传递的准确性,在划分有限元网格时保证导边和随边以及叶梢附近的网格要密一些.经多次试算后的有限元网格如图3所示.考虑到螺旋桨自身质量和运转工况,在ANSYS中设置旋转轴和施加旋转速度即可施加离心力和重力作用,桨叶根部边界为固支端.计算模型采用的材料参数如表1所示,其中S玻璃纤维为实验所测板拉结果,计算出的铝青铜材料桨叶的变形比S玻纤桨叶变形低一个数量级,因此主要考虑高强S玻纤的桨叶结构.计算用复合材料桨几何和金属桨几何相同.4动态载荷分布对于流体动力载荷,由于桨叶工作于复杂的流场中,叶面和叶背受到分布载荷,这种载荷既不均匀,也不满足一些简单规律,因此如何尽可能真实地模拟桨叶的载荷分布是流固耦合分析的关键所在.4.1u3000螺旋回射线流固耦合流固耦合算法的确定为了实现结构载荷计算的准确性,首先采用几种常用的湍流模型(κ-ε,SSTκ-ω和RSM)对比求解螺旋桨的水动力敞水性能并和JESSUPSD实验结果进行了对比.对比结果如图4、5所示,螺旋桨敞水动力参数的计算值和实验结果吻合较好,计算所得的KT和10KQ与实验结果的最大误差除了κ-ε模型为12.5%以外,SSTκ-ω为9.6%和RSM为8%,而对于螺旋桨效率的计算值仅有κ-ε模型超过了11.5%,而其他情况误差均在5%.需要指出的是,最大误差值均出现于最大进速情况.原因在于相同转速和相同直径下,进速越大,相应的推力和扭矩越小,任何一个干扰就会导致预测值和实验值的误差增加很大.整体来看SSTκ-ω模型和RSM模型均能得到精度较高的计算结果.图6为在同台计算机上不同湍流模型在不同进速下的计算时间对比图.由图6可知,相同进速下RSM模型的求解时间最长,在低进速时尤为明显,κ-ε求解时间最短,而SSTκ-ω模型求解时间接近κ-ε模型.从以上对比可知,SSTκ-ω模型可以在保证求解精度的情况下,求解时间减少很多.另考虑到流固耦合计算本身对计算机求解性能和求解时间的严格要求,综合考虑,螺旋桨流体计算湍流模型采用SSTκ-ω模型.4.2叶截面弦向压力分布的预报值与试验值比较图7为螺旋桨在设计工况J=0.833时,r/R=0.3、0.7和0.9半径上JESSUPSD实验换算得到的弦向压力分布和采用SSTκ-ω计算值的比较.从压力系数比较结果可以看出,在外半径r/R=0.7和0.9处,叶切面弦向压力分布的预报值与试验值吻合较好;而在内半径r/R=0.3处,预报值与试验值出现一定的偏差.这主要是由于JESSUPSD的结果是采用势流的伯努利方程换算得到的,没有考虑旋涡和粘性的影响所致.在半径r/R=0.3处,桨轴边界层和桨叶边界层同时存在,旋涡和粘性的影响较大,造成粘性计算结果和势流换算结果的偏差.总体来看,预报压力值和试验压力结果趋势一致,吻合较好,再次说明采用SSTκ-ω湍流模型可以准确预报螺旋桨的压力分布.5桨叶锚点的确定桨叶参考线即叶剖面鼻尾线中点的坐标:设导边的坐标为(xl,r,θl)或(xl,yl,zl),随边的坐标为(xt,r,θt)或(xt,yt,zt),叶宽(即弦长)分布为C(r),则桨叶的轮廓线可表示:式中:C(r)为叶剖面弦长;φ(r)为叶剖面螺距角;下标l、t分别表示导边(取-)和随边(取+).选择无量纲弦长s,导边表示为0,随边表示为1,桨叶弦向中点为1/2,桨叶剖面的拱度和厚度分布分别表示为f(s)和t(s),桨叶拱弧面的表达式为:式中:下标c表示为拱弧面;δk=2π(k-1)/K,k=1,2…,K,为桨叶数.设变形后的导边和随边节点坐标则变形后螺旋桨桨叶参考线的坐标为变形后各叶切面半径为rref:将式(6)和式(7)代入式(3)得变形后的侧斜角再代入式(4),求出变形后的叶切面螺距角在推导式(9)的过程中,由于螺旋桨的变形量很小,所以认为C(r)不变.将代入式(4)就可得到变形后的纵倾最后将所有变形后的节点坐标代入到式(5)桨叶拱弧面的方程中,就可得到变形后的拱弧6求解流程及求解过程复合材料螺旋桨桨叶受到水动力载荷和离心力的作用,会产生较大的变形,那么桨叶结构的离散方程为式中:K为桨叶刚度矩阵,u桨叶节点位移矢量矩阵,F为桨叶所受外载荷矩阵.而式(11)中位移矢量矩阵u需采用有限元软件计算,而F求解需要流体计算软件.两者求解的平衡实现需借助流固耦合算法.因此设计了一种求解螺旋桨性能的流固耦合算法,求解流程如图8所示.计算采用稳态求解,不考虑瞬时效应.详细求解过程为:1)通过第2节中的方法建立螺旋桨的流体计算模型,计算出螺旋桨的敞水性能和压力分布.2)将1)计算的压力载荷结果通过流固耦合界面的方法传递给第3节中建立的有限元模型,计算出螺旋桨的变形,将变形后的桨叶节点坐标输出到文本文件.3)采用第5节中的方法确定变形后的螺旋桨几何参数,建立变形后的螺旋桨几何模型.4)将变形后的几何再次输入到CFD软件中求解螺旋桨的敞水性能和压力分布.5)判断是否满足平衡方程和收敛条件,如果满足则输出计算结果,包括结构变形量、应力场、推力、扭矩以及推进效率等.如不满足收敛条件,重复迭代2)、3)、4)步,直至结果收敛.文中指定的收敛准则为2个迭代步内的推力系数和扭矩系数小于2%.计算过程一般迭代2次就可达到收敛.7计算和变形后的压力分布、变形后的压力分布7.1最大mise应力分布图9是不同转速、不同进速下的MISE应力分布和变形对照图(变形均放大20倍).设计工况,n=10r/s和20r/s时,最大等效应力分别为2.3MPa和8.62MPa.由图可以看出,相同转速不同进速情况下,低进速情况下的MISE应力比高进速下的要大.不同转速,相同进速情况,转速高的应力要比转速低的要大.最大MISE应力分布均在叶根弦向中心位置.根据螺旋桨理论可知,相同转速下,在低进速时为螺旋桨的“重载”状态,推力和扭矩均较大;高进速时为“轻载”状态,推力和扭矩较小,所以在低进速时螺旋桨的变形要比高进速的大,应力分布也是一致的.相同进速,不同转速下,根据螺旋桨的无因次理论,推力和扭矩增大数为2的转速倍数次方.弹性范围应力和变形也增加同样倍数.图10和图11是不同转速、进速下的桨叶变形分布.从图10可以看出,桨叶变形引起了桨叶侧斜的改变.4119螺旋桨为无侧斜桨,变形后侧斜角为正值,即向随边倾斜.相同转速时,低进速的侧斜改变比高进速的大,这与前面分析的受力规律一致.而图11则显示了在水动力载荷的作用下,桨叶朝着船体方向变形,且在叶稍附近达到最大值,变形规律类似悬臂梁.变形后的桨叶剖面鼻尾线在螺旋桨轴向(x轴)方向的投影有零变为正值,即导致产生了纵倾的分布.另桨叶朝着船体方向变形,促使各半径处的螺距角分布改变.即在随边附件变形大,导致螺旋桨桨叶螺距的变小.7.2叶压力分布的对比图12为设计进