脉动流下两端简支细长圆柱动力学行为分析_舒亚锋.pdf

振 动 与 冲 击 第 39 卷第 13 期JOURNAL OF VIBRATION AND SHOCKVol. 39 No.13 2020 基金项目 中国科学院先导专项XDA21010202 收稿日期 2018 -11 -09 修改稿收到日期 2019 -03 -01 作者简介 舒亚锋 男,博士生,1990 年生 通信作者 武建军 男,博士,教授,博士生导师,1964 年生 E-mailwujjun lzu. edu. cn 脉动流下两端简支细长圆柱动力学行为分析 舒亚锋1,2, 武建军1, 杨永伟2, 李佳骏1 1. 兰州大学 西部灾害与环境力学教育部重点实验室,土木工程与力学学院,兰州 730000; 2. 中国科学院近代物理研究所,兰州 730000 摘 要在考虑由附加轴向力引起非线性项基础上,建立了轴向脉动流下两端简支圆柱的非线性动力学模型。 为 了研究该描述系统的动力学特性,进一步讨论关键参数脉动频率对系统动力学行为的影响。 从分岔图可以看出系统在一 定的频率范围内可能存在混沌,并且通过三个特定脉动频率下的相图、庞加莱映射、功率谱、以及最大 Lyapunov 指数等统 计特征来分析该系统的动力学特性。 该研究结果对反应堆中脉动流引起燃料棒流致振动的安全评估具有指导意义。 关键词 非线性项;脉动流;分岔图;动力学特性 中图分类号 O32;TL3 文献标志码 ADOI10. 13465/ j. cnki. jvs. 2020. 13. 008 Dynamic behavior analysis for a two-end simply-supported slender cylinder subjected to pulsating flow SHU Yafeng1,2, WU Jianjun1, YANG Yongwei2, LI Jiajun1 1. MOE Key Lab of Mechanics on Disaster and Environment in Western China, College of Civil Engineering and Mechanics, Lanzhou University, Lanzhou 730000, China; 2. Institute of Modern Physics, Chinese Academy of Sciences, Lanzhou 730000, China Abstract The nonlinear dynamic model of a two-end simply-supported cylinder subjected to axial pulsating flow was established considering a nonlinear term caused by additional axial force. In order to study dynamic characteristics of the system, effects of key parameter’s pulsating frequency on dynamic behaviors of the system were further discussed. From bifurcation diagram, it was observed that chaos may exist in the system within a certain frequency range. Phase diagrams, Poincare mapping and statistical features, such as, power spectrum and the maximum Lyapunov exponent under 3 special pulsating frequencies were used to analyze dynamic characteristics of the system. It was shown that the results provide a guide for safety assessment of flow induced vibration of a fuel rod caused by pulsating flow in a reactor system. Key words nonlinear term; pulsating flow; bifurcation diagram; dynamic characteristics 轴向流作用下细长圆柱结构的动力学行为研究在 工程中具有广泛应用,如对核能工程中核燃料元件,控 制棒,中子测量管以及海洋中柔性皮筏拖拽器,石油勘 探用的拖拽式阵列声纳等结构的设计和流致振动控制 都具有参考价值[1-2]。 Padoussis[3-5]分别建立了轴向流 中悬臂、简支、固支圆柱的动力学模型,发现在足够高 的流速下Ucr≈50 m/ s,圆柱才会发生屈曲失稳,并且 理论模型计算结果得到了实验验证,但是反应堆工程 中冷却剂流速低于临界流速,无量纲流速 u 3. 14 被 定义为次临界流速。 王琳等[6-7]主要研究了输流管稳 定性、分岔及混沌等非线性动力学特征,轴向流作用下 圆柱的动力学模型与之相似但动力行为却不一样,受 到的外力非常复杂。 Modarres-Sadeghi 等[8-9]建立了比 较精确的恒定外流下的非线性动力学模型,并且分析 了特定参数下的圆柱失稳以及首次通过理论计算发现 高流速下细长圆柱颤振失稳等非线性动力学特性。 实 际应用中,在次临界流速下脉动流是引起结构微动磨 损失效或破坏的主要原因[10-11]。 所以,研究脉动流下 细长圆柱的非线性动力学建模与数值方法研究对反应 堆安全运行具有重要意义。 本文首先建立脉动流下两端简支细长圆柱非线性 动力学方程,进行数值计算研究脉动频率参数对动力 行为的影响,通过频率参数的分岔图,以及三个特定频 率参 数 下 系 统 的 相 图, 庞 加 莱 图, 功 率 谱 和 最 大 Lyapunov 指数混沌振动数值识别方法描述圆柱结构在 脉动流中的动力学特征。 这些工作为反应堆工程中冷 却剂主泵的脉动流速设计提供参考。 ChaoXing 1 建立动力学微分方程 核反应堆中,燃料棒束可以简化为图 1a所示细 长圆柱簇模型。 为了研究单根孤立燃料棒动力学行 为,认为棒之间的流道足够宽而忽略相互影响。 图 1 b为单根棒简化的二维简支梁模型示意图,Dh为周 围流场域的虚拟边界直径,虚线为变形后形状。 冷却 剂脉动流速 U U01 μsin Ωt。 下面对流体和圆柱 体做合理假设① 流体不可压缩;② 圆柱为细长体, Euler-Bernoulli 梁理论适用;③ 圆柱横向变形很小;④ 圆柱轴线是可伸长的。 a 圆柱簇模型 b 单根圆柱模型 图 1 反应堆堆芯燃料棒简化模型示意图 Fig. 1 Fuel rod simplified model sketch for reactor core 直接应用纳维斯托克斯方程求得流体力很复杂, 因此将圆柱受力分为非黏性流体力、流体黏性力、静水 压力三部分分别计算。 对图 1b模型沿轴向取微元 受力分析如图 2 所示。 图 2 圆柱微元受力分析 Fig. 2 Cylindrical element force analysis 根据 Lighthill 细长体理论[12],非黏性流体力表达 式为 FA ∂ ∂t U ∂ ∂x M ∂y ∂t U ∂y ∂x [] 1 因此,脉动流引起的附加非黏性流体力 Fadd MU ∂y ∂x 2 黏性流体力表达式由 Taylor[13]通过实验给出 FN 1 2 ρDU2CNsin i CDpsin2i3 FL 1 2 ρDU2CTcos i4 根据小变形假设,i 是微元变形前后轴线的夹角, 如图 2 所示,非常小。 i arctan ∂y ∂x arctan ∂y ∂t U [] 5 当∂y ∂x和 ∂y ∂t U 很小时,sin i≈∂y ∂x ∂y ∂t U,cos i≈1, 式3、式4可写作 FN 1 2 ρDUCN ∂y ∂t U ∂y ∂x 1 2 ρDCDp ∂y ∂t 2 2U ∂y ∂t ∂y ∂x U2 ∂y ∂x 2 [] 6 FL 1 2 ρDU2CT7 将二次黏性力 1 2 ρDCDp ∂y ∂t 2 2U ∂y ∂t ∂y ∂x U2 ∂y ∂x 2 []在 流速为零线性平均化为[14]1/2ρDCD ∂y ∂t ,其中 CDp为 形状阻力系数,CD为零流速的形状阻力系数有量纲。 那么式6为 FN 1 2 ρDUCN ∂y ∂t U ∂y ∂x 1 2 ρDCD ∂y ∂t 8 FN、FL分别为某轴向位置处沿法向和轴向的黏性流体 力的表达式,单位N/ m。 对于足够大流道中单根圆柱,由达西公式 ΔP fρLV2/2/ De 知,流体流动引起的压力差是可以忽略的, 因此 ∂P ∂x ρg9 对于图 1a所示圆柱簇,由流体流动引起的压力 差是不可以忽略的,它与作用在圆柱体的摩擦力有关。 截取长度为 δx 的流体单元,近似有 - Ach ∂P ∂x - Ff ρgAch 010 Ach为流道中流体的横截面面积,Ff为总摩擦力。 假设 各个横截面上压力是均匀的,得 Ff∑ i FLi Fch∝ ∑ i S Sch 11 S 为单位长度单个圆柱的侧面面积,Sch为整个流道的 侧面积,从式11可得到 Ff FL Stot S 12 其中,Stot为单位长度总的侧表面积,FL由式7给出。 05振 动 与 冲 击 2020 年第 39 卷 ChaoXing 将式12代入式10得 - Ach ∂P ∂x - FL Stot S ρgAch 013 两边同乘以 A/ Ach得 - A ∂P ∂x - FL 4A/ S 4Ach/ Stot ρgA 0 14 因为 D 4A/ S,Dh 4Ach/ Stot,其中 D 为圆柱直径,Dh 水力学直径。 把式7代入式14可得 A ∂P ∂x ρgA - 1 2 ρDU2CT D Dh 15 对上式从 x 到 L 进行积分 APx APL 1 2 ρDU2CT D Dh - ρgA[]L - x 16 Modarres-Sadeghi 给出了圆柱端部截面压力表达式 ApL 1 - 2υp -A ρgAL 2 []δ 17 式中p - 为 x 1/2L 处的静压值; - 2υA 为径向收缩引 起的压缩载荷;δ 1 为圆柱两端轴向固定,δ 0 为圆 柱下端可以轴向滑动。 由微元所受的力和力矩平衡推导圆柱运动微分 方程。 x 方向的力平衡方程 ∂T ∂x mg FL FA FN ∂y ∂x - FPx 018 y 方向的力平衡方程 ∂Q ∂x - FA FN FPy FL ∂y ∂x ∂ ∂x T ∂y ∂x - m ∂2y ∂t2 019 弯矩平衡 Q - ∂M∗ ∂x 20 其中,弯矩 M∗ EI ∂2y ∂2x 。 假设圆柱为黏弹性材 料,且符合 Kelvin-Voigt 假定,即应力-应变关系满足关 系式 σ E E∗ ∂ ∂t ε,则式20可写作 Q -E E∗ ∂ ∂t I ∂3y ∂x3 21 Fpx和 Fpy是由作用在单元圆柱外表面的静水压力 引起的力。 通过积分可得,但用以下方法更为方便。 将图 2 所示圆柱微元“冻结”并将其表面浸没在流体 中,如图 3 所示微元受流体静压。 那么其上、下表面就 会产生压力 PA 和 ∂PA ∂x []δx。 这些力与 Fpxδx,Fpyδx 的合力等于浮力。 图 3 等效刚体微元所受外力 Fig. 3 The external resultant force of equivalent rigid body element 假设 Px是 x 的线性函数,有 - ∂PA ∂x - Fpx[]δx􀭴i Fpy- ∂ ∂x PA ∂y ∂x []δx􀭴 j -∯PndA -∭ vol ▽Pdvol - ∂P ∂xAδx􀭴 i22 由式22可得 - Fpx - ∂P ∂xA ∂PA ∂x P ∂A ∂x Fpy ∂ ∂x PA ∂y ∂x 23 若横截面积不变,即∂A ∂x 0,则 - Fpx 0 Fpy A ∂P ∂x ∂y ∂x PA ∂2y ∂x2 24 将式24代入式18和19得 x 方向力平衡方程 ∂T ∂x mg FL FA FN ∂y ∂x 025 y 方向力平衡方程 ∂Q ∂x - FA FN A ∂P ∂x ∂y ∂x AP ∂2y ∂x2 FL ∂y ∂x ∂ ∂x T ∂y ∂x - m ∂2y ∂t2 026 根据小变形理论,假设横向位移是一阶,即 y οε,那么转角、弯矩、剪力均为一阶小量,所以FA FN∂y ∂x οε2为二阶小量。 将式25中的二阶小量忽略后得 ∂T ∂x mg FL 027 对其两边从 x 到 L 积分得 Tx TL mgL - x ∫ L x FLdx28 TL是圆柱体下端所受拉力。 Modarres-Sadeghi 给出 表达式 TL Tδ 1 2 ρD2U2Cb1 - δ - 15第 13 期舒亚锋等 脉动流下两端简支细长圆柱动力学行为分析 ChaoXing L 2 δ 1 2 ρDU2CT1 D Dh mg 29 式中T 是外部施加的均匀拉力;Cb是基本拉力系数。 对于两端铰支或者固定的圆柱,必须考虑弯曲引 起的轴向伸长对轴向力的影响。 圆柱体轴向微段长度 为 dx 并考虑其弯曲变形,如图 4 所示。 图 4 圆柱微元变形 Fig. 4 Element deation 圆柱微段 dx 的伸长量为 ds 1 y′2- 1dx30 由泰勒级数展开并取一阶近似得 ds 1 2 y′2dx31 其中 y′ ∂y ∂x。 于是,圆柱的应变为 ε ds dx 1 2 y′232 由此可见每个位置处的应变是不同的,由应力-应 变关系得附加应力 σ E E∗ ∂ ∂t 1 2 y′233 附加轴向力项 T∗EA E∗A ∂ ∂t 1 2 y′234 由于离散化时,附加轴向力项计算比较复杂主要 涉及到化为状态方程。 因此,这里使用平均应变代替 一般应变,总伸长量为 s 1 2∫ L 0 y′2dx ,那么平均应变为 ε s L 1 2L∫ L 0 y′2dx35 平均附加应力为 σ E E∗ ∂ ∂t 1 2L∫ L 0 y′2dx36 附加轴向力为 T∗ EA 2L∫ L 0 y′2dx E∗A L ∫ L 0 y′y ′dx 37 所以考虑附加轴向力后式28变为 Tx TL mgL - x ∫ L x FLdx EA 2L∫ L 0 y′2dx E∗A L ∫ L 0 y′y ′dx 38 利用式21、25、2638可以推导得圆柱横 向运动微分方程 E E∗ ∂ ∂t I ∂4y ∂x4 MU2 ∂2y ∂x2 2MU ∂2y ∂x∂t MU ∂y ∂x M m ∂2y ∂t2 - EA 2L∫ L 0 y′2dx E∗A L ∫ L 0 y′y ′dx []∂ 2y ∂x2 - {[1 - 2υPA T]δ [ 1 2 ρDU2CT1 D Dh m - ρAg]L - δ 2 L - x 1 2 ρD2U2Cb1 - δ} ∂2y ∂x2 m - ρAg 1 2 ρDU2CN CT D Dh []∂y ∂x 1 2 ρDUCN ∂y ∂t 1 2 ρDCD ∂y ∂t 039 式中m 为单位长度圆柱质量;ρ 为流体密度;I 为截面 惯性矩;E 为弹性模量;M 为单位长度流体附加质量;A 为圆柱横截面积;δ 为支撑端可滑动取值 0,不可滑动 取 1;CN为法向摩擦阻力系数;CT为切向摩擦阻力系 数;g 为重力加速度;L 为圆柱长度。 为了给出无量纲运动方程,下面首先引入无量纲 变量和无量纲参数 η y L , ξ x L , u0 ρA EIU0L, β ρA m ρA, γ m - ρAgL3 EI , τ EI m ρA 1/2 t L2 , Π0 EAL2 EI , ω m ρA EI 1/2 ΩL2, cn 4 π CN, c 4 π CD ρA EI 1/2 L, ct 4 π CT, cb 4 π Cb, ε L D , Π PAL2 EI , Γ TL2 EI , α EI ρA m E∗ L2 , h D Dh 40 无量纲化脉动流速 u u01 μsin ωτ41 代入式39,得到无量纲化处理后的方程 αη 4 η4 χu2η″2χβ1/2uη ′ χβ1/2u 0μωcos ωτη′ [1 χ -1β]η - Π0 2 ∫ 1 0 η′2dξ αΠ0∫ 1 0 η′η ′dξ []η″ - {[1 - 2υΠ Γ]δ 1 2 εctu21 h γ[] 1 - 1 2 δ- ξ[] 1 2 εcbu21 - δ}η″ 1 2 εcnβ1/2uη 1 2 εcβ1/2η 1 2 εu2cn hct γη′ 042 2 动力学方程的 Galerkin 离散 为便于求解运动微分方程,使用 Galerkin 方法进行 离散,将无量纲化方程42进行离散化处理并降阶。 ηξ,τ ∑ ∞ i 1 ϕiξqiτ ΦTq43 方程矩阵形式如下 25振 动 与 冲 击 2020 年第 39 卷 ChaoXing αΦ4Tq Φ4Tq χu2Φ″Tq 2χβ1/2uΦ′Tq χβ1/2u0μωcos ωτΦ′Tq [1 χ - 1β]ΦTq - Π0 2 ∫ 1 0 qTΦ′Φ′Tqdξ αΠ0∫ 1 0 qTΦ′Φ′Tq dξ []Φ″ Tq - {[1 - 2υΠ Γ]δ 1 2 εctu21 h γ[] 1 - 1 2 δ- ξ[] 1 2 εcbu21 - δ}Φ″Tq 1 2 εcnβ1/2uΦTq 1 2 εcβ1/2ΦTq 1 2 εu2cn hct γΦ′ Tq 0 44 两边乘 Φ 得 αΦΦ4Tq ΦΦ4Tq χu2ΦΦ″Tq 2χβ1/2uΦΦ′Tq χβ1/2u0μωcos ωτΦΦ′Tq [1 χ - 1β]ΦΦTq - Π0 2 ∫ 1 0 qTΦ′Φ′Tqdξ αΠ0∫ 1 0 qTΦ′Φ′Tq dξ []ΦΦ″ Tq - {[1 - 2υΠ Γ]δ 1 2 εctu21 h γ[] 1 - 1 2 δ- ξ[] 1 2 εcbu21 - δ}ΦΦ″Tq 1 2 εcnβ1/2uΦΦTq 1 2 εcβ1/2ΦΦTq 1 2 εu2cn hct γΦΦ′ Tq 0 45 在区间[0, 1]上关于 ξ 积分,再利用模态函数的正交 性,作如下替代 ∫ 1 0 ΦΦTdξ I,∫ 1 0 ΦΦ′Tdξ B, ∫ 1 0 Φ′Φ′Tdξ - C,∫ 1 0 ΦΦ″Tdξ C, ∫ 1 0 ξΦΦ″Tdξ D,∫ 1 0 ΦΦ4Tdξ Λ46 其中,I 为单位矩阵,Λ 为对角阵,它的对角元素对应梁 频率方程的特征值。 αΛq Λq χu2Cq2χβ1/2uBq χβ1/2u0μωcos ωτBq [1 χ - 1β]Iq Π0 2 qTCqCq αΠ0qTCq Cq [] 1 2 εctu21 h γ[]Dq -{[1 - 2υΠ Γ]δ 1 2 εctu21 h γ[]1 - 1 2 δ 1 2 εcbu21 - δ}Cq 1 2 εcnβ1/2uIq 1 2 εcβ1/2Iq 1 2 εu2cn hct γBq 047 两端铰支梁模态函数为 ϕi2sinλiξ, λ1 π, λ2 2π, , λi iπ48 式46 中 B,C,D 矩阵各元素 bij,cij,dij具体表达式 如下 bij 2λiλj λ2 j - λ2 i [ - 1 ij - 1],i ≠ j 0,i j { , cij 0,i ≠ j - λ2 j, i j { , dij 4λiλ3 j λ2 j - λ2 i 2[1 - - 1 ij], i ≠ j cjj/2,i j , Λij λ4 iδij 49 根据式47得到降阶后离散化方程矩阵形式 Mq Cq Kq hq,q 0 50 作如下变换 q M-1Cq M-1Kq M-1hq,q 0 51 Mij 1 χ -1βδij,δij为 Kronecker 记号。 Cij αλ4 iδij 2χβ1/2u0bij 2χβ1/2u0μsin ωτbij 1 2 εcnβ1/2u0 1 2 εcβ1/2δij 1 2 εcnβ1/2u0μsin ωτδij Kij λ4 iδij χu2 0cij 2χu2 0μsin ωτcij χu2 0μ 2sin2 ωτcij χβ1/2u0μωcos ωτbij-{[1 - 2υΠ Γ]δ 1 2 εctu2 01 h γ[] 1 - 1 2 δ 1 2 εcbu2 01 - δ}cij 1 2 εctu2 01 h γ[]dij 1 2 εu2 0cn hct γbij- u2 0μ 2sin2 ωτ{ 1 2 εct1 h[]1 - 1 2 δ 1 2 εcb1 - δ}cij- u2 0μsin ωτ{[εct1 h] 1 - 1 2 δ εcb1 - δ}cij u2 0μ 2sin2 ωτ 1 2 εct1 hdij u2 0εctμsin ωτ1 hdij 1 2 εu2 0μ 2sin2 ωτcn hctbij εu2 0μsin ωτcn hctbij h Π0 2 qTCqCq αΠ0qTCq Cq 为方便计算,引入状态变量 z [q q ]T,将式 51化为状态方程 z Az μωcos ωτB1 sin ωτB2z μ2sin2ωτB3z Hz52 其中系数矩阵的具体表达式如下 A 0I - M-1K′ ij - M-1C′ ij , B1 00 - χβ1/2u0M-1bij0 , H 0N1 - M-1h 35第 13 期舒亚锋等 脉动流下两端简支细长圆柱动力学行为分析 ChaoXing K′ ij λ4 iδij χu2 0cij 1 2 εctu2 01 h γ[]dij 1 2 εu2 0cn hct γbij-{[1 - 2υΠ Γ]δ 1 2 εctu2 01 h γ[] 1 - 1 2 δ 1 2 εcbu2 01 - δ}cij, C′ij αλ4 iδij 2χβ1/2u0bij 1 2 εcnβ1/2u0 1 2 εcβ1/2δij B2 00 - M-1[2χu2 0cij - u2 0{[εct1 h] 1 - 1 2 δ εcb1 - δ}cij u2 0εct1 hdij εu2 0cn hctbij] - M-12χβ1/2u0bij 1 2 εcnβ1/2u0δij B3 00 M-1[- χu2 0cij u2 0{ 1 2 εct1 h[]1 - 1 2 δ 1 2 εcb1 - δ}cij- 1 2 εctu2 01 hdij - 1 2 εu2 0cn hctbij] 0 3 数值计算分析 首先验证本文建立模型的正确性,下面对恒定流 u u0下本文模型的特征值、屈曲失稳幅值计算结果 进行验证。 在 相 同 参 数 下, 图 5 给 出 本 文 简 化 模 型 与 Modarres-Sadeghi 模型特征值实部、虚部随流速变化的 关系对比,发现本文模型计算的系统稳定区域为 0 u0 π,与 Modarres-Sadeghi 模型计算稳定区域一致。 图 6 给出系统发生一阶模态发散失稳的屈曲幅值,其 值与 Modarres-Sadeghi 模型计算值基本相等,从而图 5、 图 6 验证了本文模型正确性。 a 实部 b 虚部 图 5 特征值实部、虚部随流速变化 Fig. 5 Real part and imaginary part of eigenvalues with flow velocity 图 6 简支圆柱的分岔图 Fig. 6 Bifurcation diagram of a simply-supported cylinder 为了研究脉动频率参数作用下系统的动力特性, 脉动 幅 值 参 数 μ 取 值 为 0. 47, 其 它 给 定 参 数 与 Modarres-Sadeghi 取值一致如下β 0. 47,ε 16,γ 0. 838,Π04 000,χ 1,α 0,υ 0. 47,Γ 0,h 0, cb0,u0 3,cn ct 0. 025。 当 ε ≥10 时, Euler- Bernoulli 梁与铁木辛柯梁第一、第二阶固有频率比较 相对误差在 5内[15]。 本文主要分析低阶固有频率的 影响,因此使用了 Euler-Bernoulli 梁理论,忽略剪切变 形和转动效应。 采用 Runge-Kutta 数值积分法就脉动 频率对系统动力学行为影响进行数值研究,由文献[7] 知,Galerkin 离散阶数 N 2 可得到较高的计算精度。 首先,给出梁中点位移关于脉动频率参数的全局分岔 图,如图 7 所示,横坐标为可变的激励频率值,纵坐标 为梁中点振动的位移幅值,当中点速度为零时记录此 时的瞬时位移。 从图 7 可以看出,ω 在频率[1,10]范围内,存在单 周期、二周期以及混沌运动,但是不能明显看出更高倍 图 7 频率参数 ω 全局分岔图 Fig. 7 Frequency parameter ω global bifurcation diagram 周期分岔,与普通标准分岔过程不同,周期倍化过程不 连续。 为了确定系统的动力学特征,下面对比较典型三 个频率参数 ω 6. 5,ω 6. 774,ω 1 分别作相轨迹 图,庞加莱映射图,功率谱图以及 Lyapunov 指数计算如 图 8、图 9、图 10 所示。 庞加莱映射图采用脉动频率 ω 作为制作截面的依据,并且速度方向取正;根据功率谱 45振 动 与 冲 击 2020 年第 39 卷 ChaoXing 定义 Pω lim T→∞ 1 2πT XTω 2, 表示某一物理量功 率强度在频率坐标上的分布,本文采用快速傅里叶 变换方法计算。 Lyapunov 指数是表征混沌运动的重要 方法之一[16-18],本文采用正交法进行计算。 a 相图 b 庞加莱映射图 c 功率谱密度图 d 最大 Lyapunov 指数图 图 8 系统单周期运动 Fig. 8 Single periodic motion a 相图 b 庞加莱映射图 c 功率谱密度图 d 最大 Lyapunov 指数图 图 9 系统二周期运动 Fig. 9 Two periodic motion a 相图 b 庞加莱映射图 c 功率谱密度图 d 最大 Lyapunov 指数图 图 10 系统混沌运动 Fig. 10 Chaotic motion 从图 8a、图 9a看出,相轨迹出现闭合轨线,对 应的庞加莱图 8b、9b 分别为单点和两点,图 8 c、9c 功率谱曲线仅有竖线形式的谱尖点,图 8 d、9d所示,最大 Lyapunov 指数均收敛于零。 这些 非线性动力学特征都表明 ω 6. 5,ω 6. 774 系统分别 为单周期和二周期运动。 而图 10b所示庞加莱映射 既不是有限点集也不是封闭曲线,某一局部放大与整 体具有自相似结构;图10c功率谱曲线具有连续宽频 性质; 图 10 d 所 示, 最 大 Lyapunov 指 数 收 敛 于 0. 056 3大于零,这些特征表明 ω 1 系统处于混沌 运动状态。 综上所述,通过频率参数对系统动力学特性的影 响分析,在一定脉动幅值 μ 和脉动频率 ω 下,系统将会 发生混沌运动。 例如,秦山一期 RER700 型主泵的转速 1 488 r/ min,对应无量纲频率 ω 7 此种情况下系统处 于混沌发生区域。 因此,根据上述分岔图可知,需要改 变其它设计参数来避开混沌发生区域。 4 结 论 考虑轴向力引起的非线性项和脉动流,本文建立 了轴向脉动流作用下圆柱的非线性运动微分方程,数 值算例研究了较低流速无量纲流速 u03. 14下脉动 流对系统非线性动力学特性的影响。 通过频率参数 ω 全局分岔图分析系统混沌发生可能的区域,以及系统 周期运动,倍周期直至混沌发生,并且发现系统的脉动 频率在低频下,很容易发生混沌。 最后,分别通过三个 低频脉动频率点 ω 6. 5,ω 6. 774,ω 1 的相图、庞 加莱映射、功率谱及最大 Lyapunov 指数计算,给出这种 流致振动系统混沌振动的数值识别方法。 因而本文建 立的动力学模型和研究频率参数对动力学系统影响为 反应堆结构动力学系统安全设计和分析奠定了可靠的 理论基础。 55第 13 期舒亚锋等 脉动流下两端简支细长圆柱动力学行为分析 ChaoXing 参 考 文 献 [ 1] HOFSTEDE E T, KOTTAPALLI S, SHAMS A. Numerical prediction of flow induced vibrations in nuclear reactor applications[ J].Nuclear Engineering and Design, 2017, 319 81-90. [ 2] MARCUM W R, WOODS B G.Predicting the onset of dynamic instability of a cylindrical plate under axial flow conditions[J]. Nuclear Engineering and Design, 2012, 250 81-100. [ 3] PADOUSSIS M P. Dynamics of flexible slender cylinders in axial flow, Part 1 Theory[J]. Journal of Fluid Mechanics, 1966, 264 717-736. [ 4] PADOUSSIS M P. Dynamics of flexible slender cylinders in axial flow, Part 2 Experiments [ J].Journal of Fluid Mechanics, 1966, 264 737-751. [ 5] PADOUSSIS M P. The dynamical behaviour of cylindrical structures in axial flow[J]. Annals of Nuclear Science and Engineering, 1974183-106. [ 6] 王琳. 输流管道的稳定性、分岔与混沌行为研究[D]. 武 汉华中科技大学,2006. [ 7] JIN J D. Stability and chaotic motions of a restrained pipe conveying fluid[J]. Journal of Sound and Vibration, 1997, 208 427-439. [ 8] MODARRES-SADEGHI Y. Nonlinear dynamics of a slender flexible cylinder subjected to axial flow [ D].Montreal McGill University, 2006. [ 9] MODARRES-SADEGHI