输流管道线性流固耦合计算方法在核电工程中的应用

孙金雄 周瑞

摘 要：管道流致振动现象广泛存在于核电厂内，使管道产生疲劳和噪声，对管系安全造成威胁。文章对管道流致振动的原理进行探讨，并通过对两端简支输流管振动总响应计算方法的简化，得到输流管线挠度与流速的关系式，并与流固耦合数值仿真计算进行对比，两者结果非常吻合，这表明用线性流固耦合仿真技术为核电厂管道设计与验证提供参考是可行的。

关键词：流固耦合；流致振动；管道内流

1 概述

核电厂内分布大量管道系统，运行时发生流致振动现象十分普遍，振动会使管道产生疲劳和噪声，对管系安全造成威胁。常用工业规范（ASME、RCCM）对于管道流致振动的控制与验证未给出明确方法，因此需从原理上探索流致振动的原因，并找到可应用于工程的实用方法。虽然工业规范不完善，但流致振动问题却是国内外学者的研究热点，产生了很多重要研究成果，比如Blevins早在1990年就发表了经典专著[3]，对流致振动原因进行总结和分析；Paidoussis于2014年发表的管道内流流致振动专著[1]，对管道线性与非线性流固耦合方程进行全面介绍。以上研究都基于理论公式推导，并未与当前商用数值仿真软件融合。文章通过对输流管道线性运动方程的推导和简化，得到管内流速与管道最大挠度的关系，并与用仿真软件计算的结果对比，两者结果非常吻合，表明数值仿真方法在输流管道线性流固耦合分析中精度较高，可用于工程设计。

2 简支输流管线性流固耦合运动方程

2.1 简支输流管线计算模型

如图1所示，两端简支的输流管线偏离平衡位置，形成横向挠度Y，跨距为L的跨距。我们从管道上切下管道和流体两个微元，为每个微元进行受力分析。流体密度为？籽；压力为p；流速为v。管子横截面积A；长度为L；弹性模量为E；截面惯性矩为I；管道单位长度质量为M。

2.2 管道与流体微元平衡方程

由于管道振动，流体产生加速度、同时竖直方向有流体压力分量、管壁作用在流体上的压力F，上述力在竖直方向平衡可得：

流体沿管子长度方向压力梯度由管壁摩擦的剪切应力平衡，其方程为：

其中S为管内壁周界，q为管内表面剪切应力。

定义M为无水管道单位长度质量，管子横向剪切应力为Q，管子纵向拉应力为T，管线竖直方向加速度为：

定义T为管子纵向拉应力，Q为管子横向剪切应力，沿管轴线方向力平衡方程为：

联立（1）、（2）、（3）消去变量F和Q，剪切应力q从（2）（4）中消去，最终得到输流管道的自由振动的运动方程：

3 方程的解及管道变形的理论近似解

3.1 方程的边界条件和方程的解

对于两端简支输流管的边界条件为如下：

当x=0和x=L时，管道挠度与加速度为0。

由参考文献[4]知该方程的解的表达式如下：

其中j=1，2，3，…，wj是j振型的固有频率，表示为具有正弦和余弦的对称和反对称的各个空间振型的总和。

因为系统具有无限个自然振型，而实用解只需考虑头几阶振型。根据参考文献[3]，可以得到管内流体流动时管线一阶固有频率w1和管内流体静止时管道固有频率wN比值近似为：

其中vc为管道静态失稳时的临界流速：

3.2 管道最大变形与流速的关系推导

对于前面计算的两端简支管道，单位长度载荷为q。则该管道中心处挠度？啄n计算公式为：

管道固有频率wn计算公式为：

通常只考虑前两阶振型的响应，由参考文献[3]知第1阶振型a1 与第2阶振型a2比值为：

因此当管道以基波固有频率振动时，正弦波弯曲的管道第1阶振型在管道的响应中起主导作用。如果将第1阶的响应值近似等于管道总响应值，则可以得到最大挠度？啄n与第一阶固有频率wn1的关系如下：

由上式可知，挠度仅与固有频率的平方成反比关系。当管内水流动时，其挠度为？啄i，一阶固有频率为wi1，则：

由公式（10），得出输流管线挠度？啄i与？啄n的关系为：

（14）

4 软件流固耦合计算

4.1 双向迭代耦合与直接耦合区别

双向迭代耦合是指流体和结构用各自求解器在时域积分，通过耦合面在流体和结构之间进行位移与压力的传递。

直接耦合求解是通过数学一致化的处理方法，将结构方程组、流体方程组、FSI 界面方程组等不同性质的方程组、借助于位移、压力在不同方程组中的等价性，形成一个方程组耦合系统方程组。

4.2 计算模型与相关参数

4.2.1 几何模型

计算几何模型与图1相同，结构为静态分析，流体为稳态分析。由于ADINA支持直接耦合求解，所以该算例在ADINA-FSI中进行。

4.2.2 材料参数

计算采用国际标准单位，材料参数见表1。

4.2.3 边界条件

流体为实体单元，管道为壳单元；管道与流体耦合面为有滑移耦合面，选择滑移面是为匹配理论近似解的假设，同时滑移边界对网格粗糙不敏感，無需边界层。重力加速度为9.81m/s2，在入口处施加线性递增流速，最大值为临界失稳速度；由公式（8）知vc=149.6931m/s。本次为对比迭代与直接耦合计算差异，进行了两次计算，耦合方式分别为迭代耦合和直接耦合，其他边界条件均未发生变化。

4.3 计算模型与相关参数

4.3.1 ADINA计算结果

（1）迭代耦合与直接耦合。迭代耦合计算出管道最大竖直位移值为-4.668mm；直接耦合计算出最大值为-4.713mm，与迭代耦合几乎一致。

（2）理论近似解由公式（9）和算例参数，可以得到该输流管道静止时，其中心处的挠度？啄n=2.1112×10-4m。根据公式（14）可知管内流体流动后，管道中心处的挠度表达式为：

其中v为入口处速度。将上式作为ADINA后处理自定义函数表达式，把管道入口流速作为该表达式自变量，可得管道最大挠度随入口速度的关系曲线。

4.3.2 结果对比

将上述3种分析方法得到的速度-处挠度曲线绘到同一图内进行对比，如图2所示。

4.3.3 结果差异分析

（1）计算终止时，入口速度未达到临界速度149.6931m/s，因为接近临界速度时，管道刚度为0，产生失稳，计算由于不收敛终止。（2）近似精确解与ADINA计算结果相比，在接近临界速度时，有一定的差别，但差别较小。造成差别的主要原因是在推导公式时，为简化计算，只考虑管道的第1阶振型对流体激励的响应，后续的n阶响应都被忽略，造成估算算出的位移小于实际的响应。（3）线性流固耦合问题分析时，迭代耦合与直接耦合结果差别较小，直接耦合精度更高，但实际计算时，迭代耦合计算效率比直接耦合高很多。

5 结束语

综合以上结论可认为，基于数值模拟的线性流固耦合仿真技术，计算精度比简化公式更高，且当结构复杂，公式无法适用时，该方法优势更加明显，因此可将其推广到核电厂的管道设计应用中。

参考文献

[1]Michael P.Paidoussis.Fluid-structure Interactions Slender structures and Axial Flow：Second Edition VOLUME 1[M].UK：The Academic Press of Elsevier，2014：63-332.

[2]Shigehiko Kaneko，Tomomichi Nakamura. Flow-induced Vibrations Classifications and lessons from Practical Experiences： Second Edition[M]. USA： The Academic Press of Elsevier，2014：157-269.

[3]Robert D.Blevins.Flow-Induced Vibration：Second Edition[M]. Malabar： Krieger Publishing Company，1990.

[4]Housner，G.W.Bending Vibration of a Pipe Line Containing Flowing Fluid：[M].J.Appl.Mech，1954：205-208.

[5]張阿漫.流固耦合动力学[M].北京：国防工业出版社，2010.