CN112100892B - 一种不同边界条件下带内流柔性管道固有频率的预测方法 - Google Patents

一种不同边界条件下带内流柔性管道固有频率的预测方法 Download PDF

Info

Publication number
CN112100892B
CN112100892B CN202010987423.2A CN202010987423A CN112100892B CN 112100892 B CN112100892 B CN 112100892B CN 202010987423 A CN202010987423 A CN 202010987423A CN 112100892 B CN112100892 B CN 112100892B
Authority
CN
China
Prior art keywords
formula
flexible
expression
flexible pipeline
virtual
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202010987423.2A
Other languages
English (en)
Other versions
CN112100892A (zh
Inventor
高云
张壮壮
姜泽成
彭庚
潘港辉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Harbin Institute of Technology Weihai
Original Assignee
Harbin Institute of Technology Weihai
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Harbin Institute of Technology Weihai filed Critical Harbin Institute of Technology Weihai
Priority to CN202010987423.2A priority Critical patent/CN112100892B/zh
Publication of CN112100892A publication Critical patent/CN112100892A/zh
Application granted granted Critical
Publication of CN112100892B publication Critical patent/CN112100892B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/08—Fluids
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/14—Pipes
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Computing Systems (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Operations Research (AREA)
  • Fluid Mechanics (AREA)
  • Measuring Volume Flow (AREA)

Abstract

一种不同边界条件下带内流柔性管道固有频率预测方法。目前在考虑内部运动的流体的情况下，难以准确可靠地估算带内流运动立管的固有频率。本发明首先基于弹性体虚功原理建立内外流耦合激励下的柔性管道振动偏微分方程，然后基于伽辽金法将柔性管道的振动偏微分方程转化为常微分方程，再将基于特征值法得到不同边界条件下带内流运动的柔性管道振动固有频率，最后基于分析数据，对实例进行计算分析，从而预测固有频率随无量纲内流速度的变化关系的过程。本发明用于海洋工程领域中。

Description

一种不同边界条件下带内流柔性管道固有频率的预测方法
技术领域
本发明涉及一种预测方法，属于电数字数据处理技术领域。
背景技术
海洋蕴含着极其丰富的油气资源，这些丰富的油气资源不仅是未来世界经济与科技发展的重点，同时也是影响经济全球化的重要因素之一。为了加快我国对海洋油气资源的开发进度，国家已明确提出了“提高海洋资源开发能力、促进海洋经济发展、维护国家海洋权益、推进海洋强国战略”的发展目标。立管是深水海洋油气开发中必不可少的设备，同时也是海洋油气开发系统中最为薄弱的结构之一。因此，有必要对立管进行安全可靠的设计。在早期设计过程中，为了避免流体激励与结构响应发生共振现象，需要满足：立管的固有频率应避开流体激励(波浪以及海流等)的主频率区。
因此，对立管的固有频率进行准确估算显得尤为重要。目前绝大多数针对柔性圆柱体固有频率的研究，均是基于仅考虑外流激励展开，并没有考虑内部运动的流体。此时，外流激励对立管产生的载荷主要包括三部分：由尾部流场漩涡发放引起的升力和阻尼力、以及由结构外部流体引起的附加质量力；结构固有频率只与阻尼力和附加质量力有关，与升力无关。而对于真实的海洋油气运输立管，除了承受外部海流激励外，还得承受内部流体激励。运动的内部流体对立管产生的载荷同样包括三部分：内流惯性力、内流科氏力以及内流离心力，且这三部分载荷均会改变立管固有频率属性以及大小。因此为了更为准确可靠地估算带内流运动立管的固有频率，必须得考虑内部运动的流体，但目前在考虑内部运动的流体的情况下，缺少准确可靠地估算带内流运动立管的固有频率的相关方法。
发明内容：
针对上述问题，本发明公开了一种不同边界条件下带内流柔性管道固有频率的预测方法。
本发明所采用的技术方案为：
一种不同边界条件下带内流柔性管道固有频率的预测方法，所述预测方法首先基于弹性体虚功原理建立内外流耦合激励下的柔性管道振动偏微分方程，然后基于伽辽金法将柔性管道的振动偏微分方程转化为常微分方程，再将基于特征值法得到不同边界条件下带内流运动的柔性管道振动固有频率，最后基于分析数据，对实例进行计算分析，从而预测固有频率随无量纲内流速度v的变化关系的过程。
作为优选方案：所述预测方法包括以下步骤：
步骤一：建立内外流耦合激励下柔性管道结构振动偏微分方程：
取一长度为L、直径为D的柔性管道，柔性管道在外部均匀来流Ue以及内部均匀来流Ui耦合作用下产生的频率为固有频率，柔性管道的单位长度管道质量为mr，柔性管道单位长度内部流体质量为mi，柔性管道单位长度外部流体附加质量为mf；
建立坐标系：柔性管道的底面中心点为坐标原点，x轴方向为外部来流方向，z轴方向为内流流动方向，y轴方向为横流振动方向，依据柔性管道弹性体虚功原理建立表达式为：
δUa+δUb＝δWL+δWc+δWin (1)
上式(1)中，δUa为柔性管道轴向变形引起的虚应变能；δUb为柔性管道弯曲变形引起的虚应变能；δWL为升力在虚位移上所作的虚功；δWc为阻尼力在虚位移上所作的虚功；δWin是作用在柔性管道上的惯性力在虚位移上所作的虚功，通过对δUa、δUb、δWc、δWL、δWin的求解得到预测固有频率的相关指标。
作为优选方案：柔性管道的应变为ε，δε为应变ε的变分形式，柔性管道的曲率为κ；
δε和δκ的求解过程依次为：
取柔性管道上长度为dz的微元段，则变形后的微元长度ds表示式为：
上式(2)中，dy为微元长度ds在y方向上的投影，y′为y对空间z求一阶偏导数，ε0为初始应变，柔性管道上的动应变ε的表达式为：
上式(3)中，为因受力而产生的应变，得到应变ε的变分形式δε的表达式为：
δε＝y′·δy′ (4)
在柔性管道上任意一点的曲率κ的表达式为：
y″为y对空间z求二阶偏导数，当|y″|≤1时,κ≈y″,即得到曲率κ的变分形式δκ的表达式为：
δκ＝δy″ (6)
即得到了曲率κ的变分形式δκ。
作为优选方案：柔性管道的轴向变形引起的虚应变能δUa的求解过程为：
因柔性管道轴向变形引起的应变能的表达式为：
上式(7)中E为柔性管道的弹性模量、A为柔性管道的截面积，EA即为轴向抗拉刚度；
张力Te和应变ε之间的关系式为：
Te＝EAε
对上式(7)两边求变分，同时结合表达式(4)以及张力Te和应变ε之间的关系式，得到：
由于[δy]0＝[δy]L＝0，因此得到：
即上述过程为虚应变能δUa求解过程。
作为优选方案：柔性管道弯曲变形引起的虚应变能δUb的求解过程为：
因柔性管道弯曲变形引起的变形能的表达式为：
上式(10)中，E为柔性管道弹性模量，I为截面惯性矩，EI即为抗弯刚度；
对式(10)两边求变分，并结合表达式(6)，得到：
上式(11)中，y″′为y对空间z求三阶偏导数，y″″为y对空间z求四阶偏导数；由于[δy]0＝[δy]L＝[δy′]0＝[δy′]L＝0，因此:
作为优选方案：作用在柔性管道上的惯性力在虚位移上所作的虚功δWin的求解过程为：
作用在柔性管道上的惯性力在虚位移上所作的虚功δWin包括三个部分，分别为管道惯性力作的虚功、内部流体惯性力作的虚功以及外部附加流体惯性力作的虚功。
柔性管道的绝对加速度ar的表达式为：
上式(13)中，表示y对时间坐标t求二阶偏导数；柔性管道内部流体速度的矢量表达式为：
上式(14)中，和/>表示z以及y方向的单位矢量；/>为管道的振动速度；/>为管道内某一点单位切向量；由于不考虑轴向振动，因此/>且/>得到如下表达式：
为简便推导，将速度矢量形式改写为：
上式(16)中，uz＝Ui,由速度的全导数公式得到：
由于uz＝Ui为恒值，因此由uy的表达式可知uy与变量y无关，因此/>因此式(17)进一步表示为：
将uz＝Ui以及代入上式得到：
因此柔性管道内部流体的绝对加速度ai表达式为：
上式(20)中，表示/>因此作用在柔性管道(50)上的惯性力在虚位移上所作的虚功δWin表达式为：
即上述过程为虚功δWin的求解过程。
作为优选方案：主动力所作的虚功求解过程为升力在虚位移上所作的虚功δWL和阻尼力在虚位移上所作的虚功δWc的求解过程为：
阻尼力在虚位移上所作的虚功δWc的表达式为：
上式(22)中，cs为结构阻尼；cf为流体阻尼；
升力在虚位移上所作的虚功δWL的表达式为：
上式(23)中，FL为升力，δ为变分符号，将式(9)、(12)、(21)、(22)以及式(23)代入式(1)得到：
上式(24)中，为系统的惯性力，/>为系统的阻尼力，为内流运动所产生的科氏力，EIy(4)为弯曲变形产生的弹性力，Tey″为拉伸变形所产生的弹性力，/>为内流运动所产生的离心力；
结合流体阻尼cf,且cf的表达式为：cf＝γΩfρD2＝(2πStU/D)γρD2，其中Ωf为漩涡脱落频率；St为斯脱哈尔数；ρ为流体密度；γ为黏滞力系数，与流体阻力系数CD的关系为：γ＝CD/4πSt，代入式(24)进一步整理得到：
上式(25)中，CL0为柔性圆柱体处于静止状态下的升力系数，q(Z,T)为与柔性圆柱体上的升力系数有关的无量纲尾流变量。为系统的惯性力，为系统的阻尼力，/>为弯曲变形所产生的弹性力，为拉伸变形所产生的弹性力，/>为内流运动所产生的科氏力，为内流运动所产生的离心力，/>为升力。
作为优选方案：所述预测方法还包括以下步骤：
步骤二为基于伽辽金法将结构振动偏微分方程转化为常微分方程的过程；
首先将式(25)化为无量纲形式，转化过程中涉及的表达式为：
上式(26)中，η、ξ以及τ为无量纲振动位移、无量纲坐标位置以及无量纲时间，将式(26)代入式(25)得到带内流运动的柔性管道结构振动无量纲方程为：
上式(27)中，c为无量纲拖曳力系数；β为内流质量比；v为无量纲内流速度；α为无量纲静水拖曳力系数；μ为无量纲张力；分别表示如下：
在求柔性管道的固有频率时，略去式(27)中的阻尼力项以及升力项(αq)，并将剩下的5项同时乘以振型函数φi(ξ)，并在区间[0,1]上进行定积分得到：
基于伽辽金法取前四阶振型将η展开如下：
上式(30)中，φi(ξ)是第i阶横向位移振型函数，其具体表达式由边界条件确定；为第i阶广义坐标；将式(30)代入式(29)得到以下矩阵形式：
上式(31)中，为4×1矩阵列向量，可表示为：/>[M]、[C]以及[K]分别为4×4矩阵，且矩阵中各元素的表达式如下：
即上述过程为矩阵中各元素的表达式。
作为优选方案：所述预测方法还包括以下步骤：
步骤三为基于特征值法求不同边界条件下带内流运动的柔性管道结构振动固有频率；
为了便于求解，将上式(31)的二阶微分方程转换为一阶微分方程形式，具体为：
上式(33)中：
设{Z}的表达式为：
{Z}＝{A}eλt (35)上式(35)中，{A}为特征向量，λ为特征值，将式(35)代入式(33)得到：
(λ[I]-[Y])·{A}＝{0},[Y]＝-[B]-1[E] (36)
上式(36)中，[I]为单位矩阵，由式(36)看出：λ为[Y]的特征值，这里进一步将[Y]进行展开，得到：
求出特征值λ后，λ的虚部对应结构的固有频率ω；λ的实部对应结构的阻尼特性；
结合式(37)与式(32)得出：λ值只与β、v和μ有关。
作为优选方案：所述预测方法还包括以下步骤：
步骤四：基于分析数据，对实例进行计算分析：
取β＝0.3，研究不同边界条件以及不同无量纲张力下对固有频率随内流无量纲速度v的变化特性，分为以下两种情况：
第一种情况：当柔性管道两端处于铰接的连接状态时，此时结构振型写作：
第二种情况：当柔性管道两端处于固定连接的连接状态时，此时结构振型写作：
上式(39)中，β1、β2、β3以及β4写作：
β1＝4.73,β2＝7.853,β3＝10.996,β4＝4.5π (40)
从而得出柔性管道两端铰接以及两端固接这两种边界下柔性管道前4阶固有频率随无量纲内流速度v的变化关系，以实现不同边界条件下带内流柔性管道的固有频率的预测过程。
本发明的有益效果为：
一、本发明的计算过程科学合理，以柔性管道作为结构振型，依次基于弹性体虚功原理和伽辽金法建立并转化振动偏微分方程，便于简化计算步骤，同时还能够确保计算结果的准确性，为后续计算提供准确的数据基础，此后再利用特征值法得到不同边界条件下带内流运动的柔性管道振动固有频率，最后基于分析数据，对实例进行计算分析，从而预测固有频率随无量纲内流速度v的变化关系的过程。从而实现考虑内部运动的流体的情况下，准确可靠地估算带内流运动立管的固有频率，预测效果有效准确，具有重要的指导意义。
二、本发明具体针对两种不同边界条件以及三种不同张力下带内部流体运动的柔性管道固有频率展开了数值研究。建立了完整的带内流运动柔性管道固有频率的数值预测模型，用于分析不同边界条件以及不同张力下固有频率随内流速度的变化特性，该数值预报模型可以很好地模拟出不同边界条件下、不同内部流速以及不同张力比下的柔性结构固有频率特性。
三、本发明中两种不同边界条件为典型的边界条件，对于实际工作有重要指导意义，两种不同边界条件分别为柔性管道的两端处于固定连接状态以及柔性管道的两端处于铰接状态。两种不同边界条件覆盖边界条件的范围典型，使本发明具有通用性。
四、本发明的推导过程简便易算，计算结果全面且准确可靠，适于实际使用及推广。
附图说明：
为了易于说明，本发明由下述的具体实施及附图作以详细描述。
图1为内外流耦合作用下的柔性管道处于坐标系中的立体结构示意图；
图2a为柔性管道前4阶固有频率随v的变化关系示意图，图中柔性管道处于两端铰接的边界条件；
图2b为柔性管道前4阶固有频率随v的变化关系示意图，图中柔性管道处于两端固接的边界条件；
图3a为不同边界条件下柔性管道前1阶固有频率的第一对比图；
图3b为不同边界条件下柔性管道前2阶固有频率的第一对比图；
图4a为不同无量纲张力下柔性管道第1阶固有频率对比图，图中柔性管道处于两端铰接的边界条件；
图4b为不同无量纲张力下柔性管道第1阶固有频率对比图，图中柔性管道处于两端固接的边界条件。
具体实施方式：
为使本发明的目的、技术方案和优点更加清楚明了，以振动台混合试验原理为基础，说明应用本发明方法开展振动台混合试验的基本原理，但是应该理解，这些描述只是示例性的，而并非要限制本发明的范围。此外，在以下说明中，省略了对公知结构和技术的描述，以避免不必要地混淆本发明的概念。
具体实施方式一：结合图1、图2a、图2b、图3a、图3b、图4a和图4b说明本实施方式，本实施方式中所述预测方法首先基于弹性体虚功原理建立内外流耦合激励下建立柔性管道50的振动偏微分方程，然后基于伽辽金法将柔性管道50的振动偏微分方程转化为常微分方程，再将基于特征值法得到不同边界条件下带内流运动的柔性管道50振动固有频率，最后基于分析数据，对实例进行计算分析，从而预测固有频率随无量纲内流速度v的变化关系的过程。
具体实施方式二：本实施方式为具体实施方式一的进一步限定，所述预测方法包括以下步骤：
步骤一：建立内外流耦合激励下柔性管道结构振动偏微分方程：
取一长度为L、直径为D的柔性管道50，柔性管道50在外部均匀来流Ue以及内部均匀来流Ui耦合作用下产生的频率为固有频率，柔性管道50的单位长度管道质量为mr，柔性管道50单位长度内部流体质量为mi，柔性管道50单位长度外部流体附加质量为mf；
建立坐标系：柔性管道50的底面中心点为坐标原点，x轴方向为外部来流方向，z轴方向为内流流动方向，z轴方向也就是柔性管道50的轴向方向，y轴方向为横流振动方向，x轴方向和y轴方向均为柔性管道50的径向方向，二者之间的夹角为90度；
依据柔性管道50弹性体虚功原理建立表达式为：
δUa+δUb＝δWL+δWc+δWin (1)
上式(1)中，δUa为柔性管道50轴向变形引起的虚应变能；δUb为柔性管道50弯曲变形引起的虚应变能；δWL为升力在虚位移上所作的虚功；δWc为阻尼力在虚位移上所作的虚功；δWin是作用在柔性管道50上的惯性力在虚位移上所作的虚功，通过对δUa、δUb、δWc、δWL、δWin的求解得到预测固有频率的相关指标。
具体实施方式三：本实施方式为具体实施方式一或二的进一步限定，柔性管道50的应变为ε，δε为应变ε的变分形式，柔性管道50的曲率为κ；
δε和δκ的求解过程依次为：
取柔性管道50上长度为dz的微元段，则变形后的微元长度ds表示式为：
上式(2)中，dy为微元长度ds在y方向上的投影，y′为y对空间z求一阶偏导数，ε0为初始应变，柔性管道50上的动应变ε的表达式为：
上式(3)中，为因受力而产生的应变，得到应变ε的变分形式δε的表达式为：
δε＝y′·δy′ (4)
在柔性管道50上任意一点的曲率κ的表达式为：
y″为y对空间z求二阶偏导数，当|y″|≤1时,κ≈y″,即得到曲率κ的变分形式δκ的表达式为：
δκ＝δy″ (6)
具体实施方式四：本实施方式为具体实施方式一、二或三的进一步限定，柔性管道50的轴向变形引起的虚应变能δUa的求解过程为：
因柔性管道50轴向变形引起的应变能的表达式为：
上式(7)中E为柔性管道50的弹性模量、A为柔性管道50的截面积，EA即为轴向抗拉刚度；
张力Te和应变ε之间的关系式为：
Te＝EAε
对上式(7)两边求变分，同时结合表达式(4)以及张力Te和应变ε之间的关系式，得到：
由于[δy]0＝[δy]L＝0，因此得到：
具体实施方式五：本实施方式为具体实施方式一、二、三或四的进一步限定，柔性管道50弯曲变形引起的虚应变能δUb的求解过程为：
因柔性管道50弯曲变形引起的变形能的表达式为：
上式(10)中，E为柔性管道50弹性模量，I为截面惯性矩，EI即为抗弯刚度；对式(10)两边求变分，并结合表达式(6)，得到：
上式(11)中，y″′为y对空间z求三阶偏导数，y″″为y对空间z求四阶偏导数；由于[δy]0＝[δy]L＝[δy′]0＝[δy′]L＝0，因此:
具体实施方式六：本实施方式为具体实施方式一、二、三、四或五的进一步限定，作用在柔性管道50上的惯性力在虚位移上所作的虚功δWin的求解过程为：
作用在柔性管道50上的惯性力在虚位移上所作的虚功δWin包括三个部分，分别为管道惯性力作的虚功、内部流体惯性力作的虚功以及外部附加流体惯性力作的虚功。
柔性管道50的绝对加速度ar的表达式为：
上式(13)中，表示y对时间坐标t求二阶偏导数；柔性管道50内部流体速度的矢量表达式为：
上式(14)中，和/>表示z以及y方向的单位矢量；/>为管道的振动速度；/>为管道内某一点单位切向量；与横向振动相比，轴向振动可看作是小量，因此为了简便计算，这里不考虑轴向振动，此处的轴向振动运动为坐标系中z方向振动运动，因此/>且得到如下表达式：
为简便推导，将速度矢量形式改写为：
上式(16)中，uz＝Ui,由速度的全导数公式得到：
由于uz＝Ui为恒值，因此由uy的表达式可知uy与变量y无关，因此/>因此式(17)进一步表示为：/>
将uz＝Ui以及代入上式得到：
因此柔性管道50内部流体的绝对加速度ai表达式为：
上式(20)中，表示/>因此作用在柔性管道50上的惯性力在虚位移上所作的虚功δWin表达式为：
具体实施方式七：本实施方式为具体实施方式一、二、三、四或五的进一步限定，主动力所作的虚功求解过程为升力在虚位移上所作的虚功δWL和阻尼力在虚位移上所作的虚功δWc的求解过程为：
阻尼力在虚位移上所作的虚功δWc的表达式为：
上式(22)中，cs为结构阻尼；cf为流体阻尼；
升力在虚位移上所作的虚功δWL的表达式为：
上式(23)中，FL为升力，δ为变分符号，将式(9)、(12)、(21)、(22)以及式(23)代入式(1)得到：
上式(24)中，为系统的惯性力，/>为系统的阻尼力，为内流运动所产生的科氏力，EIy(4)为弯曲变形产生的弹性力，Tey″为拉伸变形所产生的弹性力，/>为内流运动所产生的离心力。
结合流体阻尼cf,且cf的表达式为：cf＝γΩfρD2＝(2πStU/D)γρD2，其中Ωf为漩涡脱落频率；St为斯脱哈尔数；ρ为流体密度；γ为黏滞力系数，与流体阻力系数CD的关系为：γ＝CD/4πSt，代入式(24)进一步整理得到：
/>
上式(25)中，CL0为柔性圆柱体处于静止状态下的升力系数，q(Z,T)为与柔性圆柱体上的升力系数有关的无量纲尾流变量。为系统的惯性力，为系统的阻尼力，/>为弯曲变形所产生的弹性力，为拉伸变形所产生的弹性力，/>为内流运动所产生的科氏力，为内流运动所产生的离心力，/>为升力。
具体实施方式八：本实施方式为具体实施方式一、二、三、四、五、六或七的进一步限定，所述预测方法还包括以下步骤：
步骤二为基于伽辽金法将结构振动偏微分方程转化为常微分方程的过程；
首先将式(25)化为无量纲形式，转化过程中涉及的表达式为：
上式(26)中，η、ξ以及τ为无量纲振动位移、无量纲坐标位置以及无量纲时间，将式(26)代入式(25)得到带内流运动的柔性管道结构振动无量纲方程为：
上式(27)中，c、β、v、α以及μ表达式分别为：
在求柔性管道50的固有频率时，略去式(27)中的阻尼力项以及升力项(αq)，并将剩下的5项同时乘以振型函数φi(ξ)，并在区间[0,1]上进行定积分得到：
基于伽辽金法取前四阶振型将η展开如下：
上式(30)中，φi(ξ)是第i阶横向位移振型函数，其具体表达式由边界条件确定；为第i阶广义坐标；将式(30)代入式(29)得到以下矩阵形式：
上式(31)中，为4×1矩阵列向量，可表示为：/>[M]、
[C]以及[K]分别为4×4矩阵，且矩阵中各元素的表达式如下：
即上述计算过程为矩阵中各元素的表达式。
具体实施方式九：本实施方式为具体实施方式一、二、三、四、五、六、七或八的进一步限定，所述预测方法还包括以下步骤：
步骤三为基于特征值法求不同边界条件下带内流运动的柔性管道结构振动固有频率；
为了便于求解，将上式(31)的二阶微分方程转换为一阶微分方程形式，具体为：
上式(33)中：
假设{Z}的表达式为：
{Z}＝{A}eλt (35)上式(35)中，{A}为特征向量，λ为特征值，将式(35)代入式(33)得到：
(λ[I]-[Y])·{A}＝{0},[Y]＝-[B]-1[E] (36)
上式(36)中，[I]为单位矩阵，由式(36)可看出：λ为[Y]的特征值，这里进一步将[Y]进行展开，得到：
求出特征值λ后，λ的虚部对应结构的固有频率ω；λ的实部对应结构的阻尼特性；
结合式(37)与式(32)得出：λ值只与β、v和μ有关；
具体实施方式十：本实施方式为具体实施方式一、二、三、四、五、六、七、八或九的进一步限定，所述预测方法还包括以下步骤：
所述预测方法还包括以下步骤：
步骤四：基于分析数据，对实例进行计算分析：
取β＝0.3，研究不同边界条件以及不同无量纲张力下对固有频率随内流无量纲速度v的变化特性，分为以下两种情况：
第一种情况：当柔性管道50两端处于铰接的连接状态时，此时结构振型写作：
第二种情况：当柔性管道50两端处于固定连接的连接状态时，此时结构振型写作：
上式(39)中，β1、β2、β3以及β4写作：
β1＝4.73,β2＝7.853,β3＝10.996,β4＝4.5π (40)
从而在柔性管道50两端铰接以及两端固接这两种边界下得出柔性管道50前4阶固有频率随无量纲内流速度v的变化关系，以实现不同边界条件下带内柔性管道50的固有频率的预测过程。
具体实施方式十一：本实施方式中根据本发明的预测方法得出两端铰接以及两端固接这两种边界下柔性管道50前4阶固有频率随无量纲内流速度v的变化关系数据并绘制成图，即图2a和图2b给出了两端铰接以及两端固接这两种边界下的柔性管道50前4阶固有频率，即第1阶固有频率为ω1，第2阶固有频率为ω2，第3阶固有频率为ω3，第4阶固有频率为ω41，在β＝0.3，μ＝0，且保持不变的前提下，随无量纲内流速度v的变化关系，通过结合图2a和图2b能够得出：无论是铰接边界还是固接边界，第1阶固有频率ω1均随着无量纲速度v的增加而逐渐减小；当无量纲速度v增加到一定值时，ω1＝0，此时对应的速度v称为临界速度，记为：vcr，当两端铰接时，vcr＝π；当两端固接时，vcr＝2π。对于两端铰接边界，当无量纲内流速度v超过2π时，则会出现第1阶固有频率ω1和第2阶固有频率ω2重合的现象。
图3a和图3b给出了不同边界条件下柔性管道50第1阶固有频率ω1和第2阶固有频率ω2的对比图，由图3a和图3b可知：无论边界调节是铰接还是固接，第2阶固有频率ω2均随着无量纲内流速度的增加呈下降趋势；当无量纲内流速度v保持不变时，相对于两端铰接边界条件，两端固接边界下柔性管道50的固有频率更大。
图4a和图4b给出了不同无量纲张力下柔性管道50的第1阶固有频率ω1对比图。由图4a和图4b得出：无论是铰接边界还是固接边界条件，临界速度vcr均随着无量纲张力μ的增加呈上升趋势；相对于固接边界，铰接边界下无量纲张力μ对临界速度vcr的影响更为明显。其他未提及的内容与具体实施方式一、二、三、四、五、六、七、八、九或十一相同。
具体实施方式十二：本实施方式为具体实施方式一、二、三、四、五、六、七、八、九、十或十一的进一步限定，这里仅以β＝0.3为代表研究固有属性随内流变化的特性。当内流为油时，根据β的公式进行相关计算，得到β＝0.3。
取β＝0.3的确定过程如下：
柔性管道50外径Do为0.6m，柔性管道50内径Di为0.56m，柔性管道50材料密度ρp为7850kg/m3，内部流体密度ρf为1000kg/m3，外部流体密度ρo为1025kg/m3。
对于柔性管道50的形状为圆柱体时，附加质量系数Ca取1.0。
基于上述推导过程得到的β＝0.3，β＝0.3在计算过程中保持不变能够确保本发明预测的准确性。
以上显示和描述了本发明的基本原理和主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。

Claims (6)

1.一种不同边界条件下带内流柔性管道固有频率的预测方法，其特征在于：所述预测方法首先基于弹性体虚功原理建立内外流耦合激励下的柔性管道(50)振动偏微分方程，然后基于伽辽金法将柔性管道(50)的振动偏微分方程转化为常微分方程，再将基于特征值法得到不同边界条件下带内流运动的柔性管道(50)振动固有频率，最后基于分析数据，对实例进行计算分析，从而预测固有频率随无量纲内流速度v的变化关系的过程；
所述预测方法包括以下步骤：
步骤一：建立内外流耦合激励下柔性管道结构振动偏微分方程：
取一长度为L、直径为D的柔性管道(50)，柔性管道(50)在外部均匀来流Ue以及内部均匀来流Ui耦合作用下产生的频率为固有频率，柔性管道(50)的单位长度管道质量为mr，柔性管道(50)单位长度内部流体质量为mi，柔性管道(50)
单位长度外部流体附加质量为mf；
建立坐标系：柔性管道(50)的底面中心点为坐标原点，x轴方向为外部来流方向，z轴方向为内流流动方向，y轴方向为横流振动方向，依据柔性管道(50)弹性体虚功原理建立表达式为：
δUa+δUb＝δWL+δWc+δWin(1)
上式(1)中，δUa为柔性管道(50)轴向变形引起的虚应变能；δUb为柔性管道(50)弯曲变形引起的虚应变能；δWL为升力在虚位移上所作的虚功；δWc为阻尼力在虚位移上所作的虚功；δWin是作用在柔性管道(50)上的惯性力在虚位移上所作的虚功，通过对δUa、δUb、δWc、δWL、δWin的求解得到预测固有频率的相关指标；
步骤二为基于伽辽金法将结构振动偏微分方程转化为常微分方程的过程；
首先将式(25)化为无量纲形式，转化过程中涉及的表达式为：
上式(26)中，η、ξ以及τ为无量纲振动位移、无量纲坐标位置以及无量纲时间，将式(26)代入式(25)得到带内流运动的柔性管道结构振动无量纲方程为：
上式(27)中，c为无量纲拖曳力系数；β为内流质量比；v为无量纲内流速度；α为无量纲静水拖曳力系数；μ为无量纲张力；分别表示如下：
在求柔性管道(50)的固有频率时，略去式(27)中的阻尼力项以及升力项(αq)，并将剩下的5项同时乘以振型函数φi(ξ)，并在区间[0,1]上进行定积分得到：
基于伽辽金法取前四阶振型将η展开如下：
上式(30)中，φi(ξ)是第i阶横向位移振型函数，其具体表达式由边界条件确定；为第i阶广义坐标；将式(30)代入式(29)得到以下矩阵形式：
上式(31)中，为4×1矩阵列向量，可表示为：/>[M]、[C]以及[K]分别为4×4矩阵，且矩阵中各元素的表达式如下：
即上述表达式为矩阵中各元素的表达式；
步骤三为基于特征值法求不同边界条件下带内流运动的柔性管道结构振动固有频率；
为了便于求解，将上式(31)的二阶微分方程转换为一阶微分方程形式，具体为：
上式(33)中：
设{Z}的表达式为：
{Z}＝{A}eλt(35)
上式(35)中，{A}为特征向量，λ为特征值，将式(35)代入式(33)得到：
(λ[I]-[Y])·{A}＝{0},[Y]＝-[B]-1[E](36)
上式(36)中，[I]为单位矩阵，由式(36)可看出：λ为[Y]的特征值，这里进一步将[Y]进行展开，得到：
求出特征值λ后，λ的虚部对应结构的固有频率ω；λ的实部对应结构的阻尼特性；
结合式(37)与式(32)得出：λ值只与β、v和μ有关；
步骤四：基于分析数据，对实例进行计算分析：
取β＝0.3，研究不同边界条件以及不同无量纲张力下对固有频率随内流无量纲速度v的变化特性，分为以下两种情况：
第一种情况：当柔性管道(50)两端处于铰接的连接状态时，此时结构振型写作：
第二种情况：当柔性管道(50)两端处于固定连接的连接状态时，此时结构振型写作：
上式(39)中，β1、β2、β3以及β4写作：
β1＝4.73,β2＝7.853,β3＝10.996,β4＝4.5π(40)
从而得出柔性管道(50)两端铰接以及两端固接这两种边界下柔性管道(50)前4阶固有频率随无量纲内流速度v的变化关系，以实现不同边界条件下带内流柔性管道(50)的固有频率的预测过程。
2.根据权利要求1所述的一种不同边界条件下带内流柔性管道固有频率的预测方法，其特征在于：柔性管道(50)的应变为ε，δε为应变ε的变分形式，柔性管道(50)
的曲率为κ；
δε和δκ的求解过程依次为：
取柔性管道(50)上长度为dz的微元段，则变形后的微元长度ds表示式为：
上式(2)中，dy为微元长度ds在y方向上的投影，y′为y对空间z求一阶偏导数，ε0为初始应变，柔性管道(50)上的动应变ε的表达式为：
上式(3)中，为因受力而产生的应变，得到应变ε的变分形式δε的表达式为：
δε＝y′·δy′(4)
在柔性管道(50)上任意一点的曲率κ的表达式为：
y″为y对空间z求二阶偏导数，当y″≤1时,κ≈y″,即得到曲率κ的变分形式δκ的表达式为：
δκ＝δy″(6)
即得到了曲率κ的变分形式δκ。
3.根据权利要求1或2所述的一种不同边界条件下带内流柔性管道固有频率的预测方法，其特征在于：柔性管道(50)的轴向变形引起的虚应变能δUa的求解过程为：
因柔性管道(50)轴向变形引起的应变能的表达式为：
上式(7)中E为柔性管道(50)的弹性模量、A为柔性管道(50)的截面积，
EA即为轴向抗拉刚度；
张力Te和应变ε之间的关系式为：
Te＝EAε
对上式(7)两边求变分，同时结合表达式(4)以及张力Te和应变ε之间的关系式，得到：
由于[δy]0＝[δy]L＝0，因此得到：
即上述过程为虚应变能δUa求解过程。
4.根据权利要求3所述的一种不同边界条件下带内流柔性管道固有频率的预测方法，其特征在于：柔性管道(50)弯曲变形引起的虚应变能δUb的求解过程为：因柔性管道(50)弯曲变形引起的变形能的表达式为：
上式(10)中，E为柔性管道(50)弹性模量，I为截面惯性矩，EI即为抗弯刚度；
对式(10)两边求变分，并结合表达式(6)，得到：
上式(11)中，y″′为y对空间z求三阶偏导数，y″″为y对空间z求四阶偏导数；由于[δy]0＝[δy]L＝[δy′]0＝[δy′]L＝0，因此得到：
即上述过程为虚应变能δUb的求解过程。
5.根据权利要求4所述的一种不同边界条件下带内流柔性管道固有频率的预测方法，其特征在于：作用在柔性管道(50)上的惯性力在虚位移上所作的虚功δWin的求解过程为：
作用在柔性管道(50)上的惯性力在虚位移上所作的虚功δWin包括三个部分，分别为管道惯性力作的虚功、内部流体惯性力作的虚功以及外部附加流体惯性力作的虚功；
柔性管道(50)的绝对加速度ar的表达式为：
上式(13)中，表示y对时间坐标t求二阶偏导数；柔性管道(50)内部流体速度的矢量表达式为：
上式(14)中，和/>表示z以及y方向的单位矢量；/>为管道的振动速度；/>为管道内某一点单位切向量；由于不考虑轴向振动，因此/>且/>得到如下表达式：
为简便推导，将速度矢量形式改写为：
上式(16)中，uz＝Ui,由速度的全导数公式得到：
由于uz＝Ui为恒值，因此由uy的表达式可知uy与变量y无关，因此/>因此式(17)进一步表示为：
将uz＝Ui以及代入上式得到：
因此柔性管道(50)内部流体的绝对加速度ai表达式为：
上式(20)中，表示/>因此作用在柔性管道(50)上的惯性力在虚位移上所作的虚功δWin表达式为：
即上述过程为虚功δWin的求解过程。
6.根据权利要求5所述的一种不同边界条件下带内流柔性管道固有频率的预测方法，其特征在于：主动力所作的虚功求解过程为升力在虚位移上所作的虚功δWL和阻尼力在虚位移上所作的虚功δWc的求解过程为：
阻尼力在虚位移上所作的虚功δWc的表达式为：
上式(22)中，cs为结构阻尼；cf为流体阻尼；
升力在虚位移上所作的虚功δWL的表达式为：
上式(23)中，FL为升力，δ为变分符号，将式(9)、(12)、(21)、(22)以及式(23)代入式(1)得到：
上式(24)中，为系统的惯性力，/>为系统的阻尼力，为内流运动所产生的科氏力，EIy(4)为弯曲变形产生的弹性力，Tey″为拉伸变形所产生的弹性力，/>为内流运动所产生的离心力；
结合流体阻尼cf,且cf的表达式为：cf＝γΩfρD2＝(2πStU/D)γρD2，其中Ωf为漩涡脱落频率；St为斯脱哈尔数；ρ为流体密度；γ为黏滞力系数，与流体阻力系数CD的关系为：γ＝CD/4πSt，代入式(24)进一步整理得到：
上式(25)中，CL0为柔性圆柱体处于静止状态下的升力系数，q(Z,T)为与柔性圆柱体上的升力系数有关的无量纲尾流变量；为系统的惯性力，/>为系统的阻尼力，/>为弯曲变形所产生的弹性力，/>为拉伸变形所产生的弹性力，/>为内流运动所产生的科氏力，/>为内流运动所产生的离心力，/>为升力。
CN202010987423.2A 2020-09-18 2020-09-18 一种不同边界条件下带内流柔性管道固有频率的预测方法 Active CN112100892B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010987423.2A CN112100892B (zh) 2020-09-18 2020-09-18 一种不同边界条件下带内流柔性管道固有频率的预测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010987423.2A CN112100892B (zh) 2020-09-18 2020-09-18 一种不同边界条件下带内流柔性管道固有频率的预测方法

Publications (2)

Publication Number Publication Date
CN112100892A CN112100892A (zh) 2020-12-18
CN112100892B true CN112100892B (zh) 2023-08-18

Family

ID=73760411

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010987423.2A Active CN112100892B (zh) 2020-09-18 2020-09-18 一种不同边界条件下带内流柔性管道固有频率的预测方法

Country Status (1)

Country Link
CN (1) CN112100892B (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113033122B (zh) * 2021-04-19 2022-04-22 哈尔滨工业大学（威海） 内流作用下柔性管道非线性响应预测方法、系统及装置
CN113111420B (zh) * 2021-04-19 2022-04-22 哈尔滨工业大学（威海） 一种边界激励细长张力梁不稳定区间的快速预测方法
CN114544125B (zh) * 2021-10-16 2024-04-26 西南石油大学 一种内外流与平台耦合作用的立管振动测试装置及方法
CN114819347B (zh) * 2022-04-26 2023-07-04 哈尔滨工业大学（威海） 多相内流激励下的输流管道振动响应预测方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH10259630A (ja) * 1997-03-18 1998-09-29 Kubota Corp 排水集合管
US7103520B1 (en) * 2000-08-31 2006-09-05 General Electric Company Methods and apparatus for model based shrouded bellows stiffness determinations
JP2012137116A (ja) * 2010-12-24 2012-07-19 Ihi Corp ポンプバレルの振動予測方法
CN109505585A (zh) * 2018-11-15 2019-03-22 西安建筑科技大学 段塞流引起的海洋立管耦合振动响应模拟装置及方法
CN109635500A (zh) * 2019-01-02 2019-04-16 西北工业大学 航空管道三维流固耦合参数共振响应特性预测方法和装置
CN109740211A (zh) * 2018-12-21 2019-05-10 西北工业大学 一种功能管道流固耦合固有特性的预测方法
CN109827734A (zh) * 2019-03-29 2019-05-31 西安建筑科技大学 一种评估内外流作用下深海立管涡激振动的方法
CN111382503A (zh) * 2020-02-27 2020-07-07 中南大学 弹性支承下旋转的柔性圆环的振动分析方法及系统

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH10259630A (ja) * 1997-03-18 1998-09-29 Kubota Corp 排水集合管
US7103520B1 (en) * 2000-08-31 2006-09-05 General Electric Company Methods and apparatus for model based shrouded bellows stiffness determinations
JP2012137116A (ja) * 2010-12-24 2012-07-19 Ihi Corp ポンプバレルの振動予測方法
CN109505585A (zh) * 2018-11-15 2019-03-22 西安建筑科技大学 段塞流引起的海洋立管耦合振动响应模拟装置及方法
CN109740211A (zh) * 2018-12-21 2019-05-10 西北工业大学 一种功能管道流固耦合固有特性的预测方法
CN109635500A (zh) * 2019-01-02 2019-04-16 西北工业大学 航空管道三维流固耦合参数共振响应特性预测方法和装置
CN109827734A (zh) * 2019-03-29 2019-05-31 西安建筑科技大学 一种评估内外流作用下深海立管涡激振动的方法
CN111382503A (zh) * 2020-02-27 2020-07-07 中南大学 弹性支承下旋转的柔性圆环的振动分析方法及系统

Also Published As

Publication number Publication date
CN112100892A (zh) 2020-12-18

Similar Documents

Publication Publication Date Title
CN112100892B (zh) 一种不同边界条件下带内流柔性管道固有频率的预测方法
Wu et al. Experimental and numerical investigation of hydroelastic response of a flexible hydrofoil in cavitating flow
Housner Bending vibrations of a pipe line containing flowing fluid
Fu et al. Features of vortex-induced vibration in oscillatory flow
Meng et al. Nonlinear free vibrations and vortex-induced vibrations of fluid-conveying steel catenary riser
Gu et al. Effect of aspect ratio on the dynamic response of a fluid-conveying pipe using the Timoshenko beam model
Gu et al. Prediction of vortex-induced vibration of long flexible cylinders modeled by a coupled nonlinear oscillator: integral transform solution
Li et al. Vector form intrinsic finite element method for nonlinear analysis of three-dimensional marine risers
CN106934138B (zh) 内输气液两相段塞流的海底管道悬跨段振动特性分析方法
CN109827734A (zh) 一种评估内外流作用下深海立管涡激振动的方法
Klaycham et al. Large amplitude vibrations of a deepwater riser conveying oscillatory internal fluid flow
CN114819347B (zh) 多相内流激励下的输流管道振动响应预测方法
Klaycham et al. Nonlinear vibration of marine riser with large displacement
Hu et al. Dynamic behaviors of a marine riser with variable length during the installation of a subsea production tree
CN115238975A (zh) 深海采矿水力提升管多场耦合非线性振动预测方法
He et al. Vortex-induced vibrations of a pipe subjected to unsynchronized support motions
Yuan et al. Free vibration analysis and numerical simulation of slightly curved pipe conveying fluid based on Timoshenko beam theory
Heidari-Soureshjani et al. Critical speed and frequency behavior of rotating joined FG-CNTRC conical-conical shells
Czerwiński et al. Experimental and numerical study on vibrations of a helical pipe with fluid flow
Huang et al. Numerical investigations of the transient cavitating vortical flow structures over a flexible NACA66 hydrofoil
Wu et al. Free and forced vibration of fluid-filled laminated cylindrical shell under hydrostatic pressure
Athisakul et al. Critical top tension for static equilibrium configuration of a steel catenary riser
Ma et al. Natural frequencies and stability for Z-shaped fluid-conveying pipes with constrained-layer damping–Semi-analytic modeling with experimental validation
Inoue et al. Experimental and Numerical Study of Motion of Rotating Drill Pipe Owing to Magnus Effect
Toorani et al. Dynamic analysis of anisotropic cylindrical shells containing flowing fluid

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant