宫亚飞， 甄亚欣

(华北电力大学 数理学院，北京 102206)

在实际的工程问题中，流体流动引发输流管系统的振动现象非常普遍。1987年，Païdoussis[1]对输流管系统的线性振动问题进行了综述。考虑泊松耦合和摩擦耦合，曾国华[2]基于中线可伸长定理,应用哈密尔顿变分原理推导了输流管道轴向以及横向运动方程。基于高阶梁模型，朱晨光等[3]推导了管内外表面剪应力为零的边界条件下输流管系统的控制方程。

日本科学家井平敏雄于1984年提出功能梯度材料(functionally graded material)的概念。功能梯度材料是一种新型复合材料，其体积分数连续平稳变化，因此没有明显的材料界面[4]。功能梯度材料具有抗高温且保持结构完整性等优良性质。现如今，功能梯度材料已广泛应用到很多领域[5]。邹俭鹏等[6-9]也针对功能梯度输流管的非线性行为进行了很多有意义的研究。

均匀直流管常用于对输流管系统的研究，但在现实生活中，由于制造、安装等问题，大多数管道都是含有初始弯曲的管道[10]。因此对含有初始弯曲的输流管道的研究是非常必要的。Li等[11]对含有几何缺陷的输送脉动流体管道的非线性振动特性进行了研究。在油气开采过程中，深海油气输送结构得到了广泛的应用，张泽超等[12]分析了初始弯曲和管道长度等因素对系统屈曲过程中动态效应的产生及影响的变化规律。

输流管系统流速较高时，通常会产生屈曲失稳现象[13-14]，如大型水压机中的压力管道以及核电厂中的热交换管道都存在此现象。高国华等[15]通过对微元体变形进行受力分析，推导出弯曲井眼中受压管柱的屈曲方程。高德利等[16]基于管柱屈曲微分方程，给出在弯扭载荷作用下受井眼约束管柱的正弦屈曲和螺旋屈曲的临界载荷。刘祥康等[17]的研究为延长油管柱寿命和优化现场生产工艺提供了理论依据。

目前，针对有初始弯曲的功能梯度输液管的屈曲分岔行为分析尚不多见。本文基于哈密尔顿原理和欧拉-伯努利梁理论，研究了两端固支边界下含有初始弯曲的功能梯度输流管屈曲行为，给出了超临界平衡位形和临界流速的解析表达式，探究了功能梯度材料及初始弯曲幅值等参数对系统的影响。

1 数学建模

对含有初始弯曲的功能梯度输流管进行建模：其中管道内、外半径分别表示为rl、ro，平均半径为r。

假设功能梯度管道的材料性能服从幂律分布[18]

E(r)=ViEi+VoEo，

ρp(r)=Viρi+Voρo

(1)

式中：E是功能梯度管道杨氏模量;ρp是功能梯度材料的组合密度;下标o和i分别表示外层和内层。

材料的体积分数表示为

(2)

式中,n为幂律指数。

图1给出了含有初始弯曲的功能梯度输液管的示意图，其中x和z分别表示轴坐标和径向坐标。L表示管道的长度，t是时间坐标。假定流体为无黏不可压缩流体，Γ为流速。v0表示管道的初始弯曲。

图1 两端固定支撑的功能梯度输流管示意图Fig.1 Schematic diagram of a functionally graded transport pipe with fixed ends

管道动能与流体动能可以表示为

(3)

(4)

式中：u(x,t)和v(x,t)分别表示管道的纵、横向位移；ρf表示流体密度；Ap为管道截面积;Af表示流体的横截面积;下标x和t分别表示关于x和t的偏导数。

工程应变ε0可以表示为[19]

(5)

式中，ds0和ds分别表示变形前后梁中心线的长度微元。

(6)

根据梁的一阶剪切变形理论，弯曲应变εb表示为

εb=-zv,xx

(7)

总应变由式(5)和式(7)中工程应变和弯曲应变之和可得

(8)

管道的弹性本构关系为

σx=Eεx

(9)

变形功的变分为

δW=-∭VσxδεxdV

(10)

由初始张力P0引起的管道的势能和流体压力引起的势能分别为

(11)

(12)

(13)

式中，p为单位面积流体压力。

基于哈密顿原理

(14)

将式(3)、(4)、(10)、(12)、(13)代入到式(14)，导出了含有初始弯曲的输流管道的非线性控制方程

(15)

式中，Ip为惯性矩。EAp和EIp可以表示为

由于两端固定支撑，管道的纵向位移远小于横向位移，因此假定u=o(v2)。令

(16)

将式(16)积分两次并代入边界条件可以得到

(17)

将式(17)代入到式(15)中，可以得到系统横向振动积分-偏微分型动力学方程

(ρpAp+ρfAf)v,tt+2ΓρfAfv,xt+ρfAfΓ2v,xx-

(Afp+P0)(v0,xx+v,xx)+EIv,xxxx-

(18)

固支边界条件为

(19)

为简化计算，引入无量纲参数

(20)

将式(20)代入到式(18)中推导出输流管横向振动无量纲控制方程为

v,tt+2γMrv,xt+γ2v,xx-Λ(v,xx+v0,xx)+v,xxxx=

(21)

边界条件为

(22)

2 直流管道的非平凡静态平衡位形

由于系统的平衡位形仅依赖于空间坐标，忽略控制方程中时间t和v0，输流管平衡方程转化为

(23)

其边界条件为

(24)

将式(23)整理为

v,xxxx+λ2v,xx=0

(25)

其中：

(26)

由式(24)，设式(25)的解为

v(x)=C1[1-cos(2mπx)]

(27)

将式(27)代入式(26)中可得

(28)

将式(24)与式(25)代入到式(22)中解得

(29)

由此得到功能梯度输液管平衡位形的解析解

v(x)=

(30)

由式(30)可得功能梯度输液管的临界流速表达式

(31)

取表1中的参数值，绘制不同幂律指数下，管道中点随流速变化的平衡位形图。从图2可以得到，随着幂律指数的增大，系统屈曲失稳的临界流速变小。

表1 物理参数Tab.1 Physical parameter

图2 不同幂律指数下位移随流速变化的分岔图(直流管)Fig.2 Bifurcation solutions of mid-point displacement with fluid velocity for different power law exponents (straight pipe)

根据临界流速的表达式(30)，绘制临界流速随功能梯度材料的幂律指数的变化图。显然图3验证了图2的结论，即幂律指数越大时临界流速越小。且当幂律指数较小时，对临界流速的影响越明显。

图3 临界流速随幂律指数的变化图(直流管)Fig.3 Critical velocity variation with power law exponent (straight pipe)

取定n=5，图4和图5描述了管道长度和初始轴向力对管道平衡位形的影响。从图4、图5可以看出，随着管道长度和初始轴向张力的增大，系统的临界流速减小。

图4 不同管道长度下位移随流速变化的分岔图(直流管)Fig.4 Bifurcation solutions of mid-point displacement with fluid velocity for different power law exponents different pipe lengths (straight pipe)

图5 不同初始张力下中点位移随流速变化的分岔图(直流管)Fig.5 Bifurcation solutions of mid-point displacement with fluid velocity for different initial axial force (straight pipe)

图6和图7研究了不同幂律指数下，管道中心点的位移及系统临界流速随流体密度的变化。图6表明，幂律指数越大，系统分岔点越小；在图7中临界流速随流体密度的增大而减小，验证了图6中所得结论。

图6 不同幂律指数下位移随流体密度变化的分岔图(直流管)Fig.6 Bifurcation solutions of mid-point displacement with fluid density for different power law exponents (straight pipe)

图7 不同幂律指数下临界流速随流体密度的变化图(直流管)Fig.7 Critical velocity variation with fluid density for different power law exponent (straight pipe)

3 含有初始弯曲的输流管道的非平凡平衡位形

忽略时间t，式(18)简化为含有初始弯曲的欧拉-伯努利梁静力平衡方程

γ2v,xx+v,xxxx-Λ(v,xx+v0,xx)-

(32)

将式(32)整理为

v,xxxx+λ2v,xx=(H2+Λ)v0,xx

(33)

其中

(34)

λ2=γ2-Λ-H2

(35)

为了描述功能梯度管道的初始弯曲，我们选择无量纲函数

v0=A0[1-cos(2πx)]

(36)

式(33)的解可以表示为

(37)

将式(37)与式(35)联立，得到

(38)

为了得到临界屈曲速度的解析解，将式(37)代入到式(33)中，可得关于C2的一元三次方程为

(39)

由一元三次方程根的判别式，给定参数A0等可导出确定的临界流速

(40)

其中

d=4A0Λπ2

(41)

选取表1中的数据，对含有初始弯曲的功能梯度输液管进行屈曲分析。

在图8中，选取初始弯曲振幅A0=0.001，根据式(38)绘制屈曲分析图，结果表明：随着流速的增大，当幂律指数n=1时比n=5时稳定。

图8 不同幂律指数下的分岔对比图(有初始弯曲)Fig.8 Bifurcation comparison for different power-law exponents (with initial curvature)

图9描绘了不同初始弯曲幅值下，幂律指数对临界流速的影响。从图9可以看出，随着幂律指数的增大，系统的临界流速减小，且随着幂律指数对临界流速的影响随幂律指数的增大而减小。同时，幂律指数不变时，初始弯曲振幅越大，临界流速越小。

图9 不同初始弯曲幅值下临界流速随幂律指数的变化图(有初始弯曲)Fig.9 Critical velocity variation with power-law exponent for different initial curvature(with initial curvature)

选取确定的流速Γ=100 m/s以及初始弯曲的振幅A0=0.001,图10描述了管道长度对分岔点的影响。结果表明：管道长度取L=11 m时比取L=13 m时稳定，即管道越短，输液管系统越稳定。

图10 不同管道长度下位移随流速变化的分岔图(有初始弯曲)Fig.10 Bifurcation solutions of mid-point displacement with fluid velocity for different pipe lengths (with initial curvature)

由于初始弯曲幅值较小时临界流速数值差距不明显，为了更清晰地体现出初始弯曲项对管道临界流速的影响，因此选用表2呈现初始弯曲对系统临界流速的影响。结果表示：直流管比具有初始弯曲的输流管道更稳定，与图9所得结论一致。

表2 不同初始弯曲幅值下临界流速随幂律指数的变化Tab.2 Variation of critical velocity with power law exponent under different initial curvature amplitude

图11绘制了两端固支和两端简支两种边界条件下含有初始缺陷的功能梯度输流管在不同幂律指数条件下临界流速随初始弯曲幅值的变化趋势。结果显示，两端简支条件下，临界流速随初始弯曲的增大而减小；同时幂律指数越大临界流速越小，管道越不稳定；两端固支条件下输流管模型临界流速趋势同两端简支条件下变化趋势一致，且两端简支条件下临界流速恒大于两端固支条件下，即两端固支情况下输流管系统更为稳定。

图11 临界流速随初始弯曲幅值的变化图(有初始弯曲)Fig.11 Critical velocity variation with initial curvature amplitude (with initial curvature)

4 结 论

通过对含有初始弯曲的功能梯度输流管平衡分岔问题进行解析研究，得到如下结论：

(1) 应用广义哈密顿原理以及欧拉-伯努利梁理论导出了固支边界条件下含有初始弯曲的功能梯度输流管的纵横耦合的非线性方程。

(2) 对不含初始缺陷的功能梯度输流管模型进行解析求解，给出了输流管道非平凡平衡位形及临界流速的解析表达式，并分析了各物理参数对系统分岔的影响。结果显示随着幂律指数、管道长度以及初始压力的增大模型的分岔点越小，且流体临界流速随着幂律指数以及流体密度的增大而减小。

(3) 推导出含有初始弯曲的功能梯度输流管模型的非平凡平衡位形及临界流速的解析表达式，并分析各物理参数对系统的影响。结果显示随着幂律指数、管道长度的增大模型的分岔点变大。流体临界流速随着幂律指数的增大逐渐趋于一个稳定值。

(4) 将直流管与含有初始弯曲的输流管系统进行对比。结果表明初始曲率越小，临界流速越大，管道越稳定。当初始弯曲幅值为0的时候，临界流速最大，管道分岔发生得越晚。

(5) 选取不同幂律指数，在两种边界条件下对输流管系统进行研究。结果显示，临界流速随着初始曲率的增大逐渐变为0；两端固支条件下临界流速恒大于两端简支，即两端固支比两端简支更加稳定。