一种蒸汽发生器管束两相横向流弹失稳分析方法
一种蒸汽发生器管束两相横向流弹失稳分析方法
【技术领域】
[0001] 本发明涉及一种蒸汽发生器管束两相横向流弹失稳分析方法，提出的两相横向流 弹失稳的解析模型，分析核反应堆中蒸汽发生器管束两相横向流弹失稳是否发生的一种方 法。
【背景技术】
[0002] 流弹失稳机理与航天工业中广泛研宄的颤振现象相同，不稳定性的主要参数是系 统的阻尼和流体弹性力。横向流作用下的管束，当流速增加到一定值时，系统中流体弹性力 所做的功将大于系统阻尼散耗的功，因此管子振动的振幅将急剧增大，这时就称出现了流 弹失稳。
[0003] 从理论上看，流体弹性力应该通过求解一组联立的流体力学和结构动力学方程得 到。但由于我们涉及的问题通常十分复杂，如大量的换热管，复杂的边界条件和高雷诺数的 流动等，使得从理论上进行求解十分困难，甚至是不可能的。通常的做法是根据大量的实验 结果，提出简化的数学物理模型。在这些模型中，流体力一般不再以独立方程的形式出现而 是表示成与结构的位移、速度和加速度等有关的函数，这种关系体现了流体作用使结构产 生位移和变形，而结构的位移和变形又改变了流体作用这样一个流固耦合振动的事实。
[0004]目前，已有分析管阵流弹失稳的模型，主要为拟静力（如Connors模型）、拟定常 (如Paidoussis模型）、非定常流模型（如Chen模型）及解析模型（Yetisir流管模型）。 分析管阵临界流速时，拟静力、拟定常和非定常流模型的流体弹性力都不同程度地依赖于 实验数据，且这些模型中的一些参数是在一定条件下、针对某个管阵形式得到的，有一定的 保守性和不确定性，并不能对所有实际问题给出可靠的临界流速。而解析模型则是根据大 量的实验结果，对流体弹性系统作一些简化，用比较成熟的数学、力学方法，得到流体弹性 力。尽管这类模型与实际情况还有一定差距，有待于进一步完善，但无疑是一类具有普遍意 义的数学模型。
[0005] 现有的解析方法计算过程比较繁琐，甚至需用到数值计算，且都是针对单相流弹 失稳的分析，如果用其分析两相流弹失稳问题，还需将气-液两相流等效为单相流，这一方 法的有效性也有待考证。
【发明内容】
[0006] 本发明的目的是：本发明是为了克服核反应堆蒸汽发生器换热管流弹失稳而导 致的管束破坏问题，而提出一种蒸汽发生器管束两相横向流弹失稳分析方法。包括以下步 骤：
[0007] 1两相横向流弹失稳参数确定
[0008] 热交换器通常有四种典型的管阵几何排列形式，对于所有的几何排列形式，解析 模型的基本概念是相同的，排列形式只是影响模型中使用的试验参数。具体以平行三角形 管束为例，而对于旋转三角形和旋转正方形管束（P，l. 7)，流体力Fy具有同样的形式，除了 它需要除以2c〇Sa，这是由于所定义的最小间距扰动面积函数a(sm，t)具有不同的形式。 [0009] 1.1基本定义
[0010] 引入"流管"概念，平行三角形管阵流管模型。曲线坐标S用以表示流动路径，流 管的平均横截面积在整个流管上假设为常数，且由-Sl位置的最小间隙面积A^确定。
[0011] A〇=min(Pcosa-D/2,P~D)
[0012] 瞬时流管面积A(s，t)会随着管子的振动而改变，A(s，t)包含两项：平均项怂，和 波动项a(s,t)。
[0013]A(s,t) =A〇+a(s,t)
[0014] 同样的，速度和压力定义为：
[0015] U(s,t) =U〇+u(s,t),P(s,t) =P〇+p(s,t)
[0016] 然而，在足够远的上游位置Si，扰动可以忽略，速度和压力可以分别视为常数UQ和 匕，流管假设在^位置开始附着于柔性管，并在ss位置分离。
[0017] 管子的振动不能立刻影响到流管的其它位置，因此，采用一个相函数供(.〇来考虑 扰动延迟效应。上游微扰函数a(s，t)可以表示为
[0019] f(s)为人工衰减函数，用以表示从柔性管的位置到上游点的微扰衰减。a(sm，t)是 最小间隙位置的面积微扰，是管子几何的函数。
[0020] 1. 2两相横向流弹失稳试验参数确定
[0021] 唯一依赖于试验输入的参数为开始产生附着和分离点的位置，这可以通过Weaver 和AbdRabbo以及Scott的可视化研宄进行估计。四种典型的管阵排列的可以归为以下两 组：
[0023] 注为流体附着表面的法向矢量角度
[0024] 1. 3两相横向流弹失稳计算参数确定
[0025] 1. 3. 1流场计算参数
[0026] (1)空泡份额
[0027] 根据均相流模型（HEM)假设，两相流的两相充分混合，并且密度和温度分布均匀， 不存在气相和液相的相对滑移，可将空泡份额定义为：
[0029] 上式中，\为真实含汽率，P8和P1分别为气相和液相的密度。
[0030] (2)流动参数
[0031] 节距流速％^，U"PAP-D)
[0032] 均相流密度P:P=Pid-ad+PgaH
[0033] 自由来流速度U":
[0034] 自由来流质量流速
[0035] 节距质量流速:
[0036] 式中，U"是自由来流速度（即假设不存在管子时的速度），P是管中心距，D是管 子直径（对于翅片管，用等效水力直径Dh代替），A为自由来流的流道截面积，Pg为气相 密度、Pi为气相密度，Vg为气相速度、Vi为气相速度。
[0037] 1. 3. 2动水参数
[0038] (1)水力直径
[0039] 对于非翅片管，水力直径就简单的为管子外径D。对于翅片管，水力直径等效为 Dh:
[0040] Dh=Dr+RF (D0-Dr)
[0041] 式中，D,为齿根圆直径，D。为翅片外径，比值RF为翅片占据的体积与管子外表面到 翅片外径表面的有效体积的比。对于翅片管，Dh要用于振动激励参数的计算，等效节径比P/ D也要变为P/Dh。
[0042] (2)水动力附加质量
[0043] 管束中一根管子单位长度的水动力附加质量可以表不为：
[0045] 式中，队是周围管子的等效直径，比值De/D是管束对管子约束程度的度量。
[0046] 对于三角形管阵：
[0047] De/D= (0? 96+0. 5P/D)P/D
[0048] 对于正方形管阵：
[0049] De/D= (1. 07+0. 56P/D)P/D
[0050] 单位长度管子总的动力质量m包括水动力附加质量％、管子本身的质量mt和管内 流体的质量mi:
[0051] m=mh+m^nii
[0052] (3)阻尼
[0053] 两相流中多跨热交换器管的阻尼包括支撑阻尼Gs、粘性阻尼^和两相阻尼 GTP:
[0054] ^T= ^ s+^y+^Tp
[0055] 通量很高、空泡份额很大的情况下，支撑就可能是干的，此时只有摩擦阻尼存在， 这种情况下，支撑阻尼与气体中的热交换器管阻尼基本相同：
[0056] ^ T= ^ s+ ^ v+ ^ TP
[0057] 管和支撑之间存在液体时，支撑阻尼Q既包括压膜阻尼GSF又包括摩擦阻尼
[0059]粘性阻尼：
[0061]两相阻尼比：
[0063] 式中，Pi是管和支撑间液体的密度，De是周围管子的等效直径，D是管子的外径， m是单位长度管子的总质量，包括用两相流均相密度计算的水动力附加质量，lm是管子的特 征长度，L是支承板宽度，f?是管子的频率，PTP是气液两相等效密度，uTP是等效运动粘度， uTP=ui/El+ajjOi/Ug-l)]，〇8和u身别为液相和气相的运动黏度，而f(aH)的确定
[0064] f ( a H) = a h/40 a H<40%
[0065] f ( a H) = 1 40%彡aH彡 70%。
[0066] f (a H) = l-(ah-70)/30 aH>70%
[0067] 2在Liming 提出的单相横向流弹失稳显式分析模型-"流管"模型的基础上，提出 针对两相横向流弹失稳的解析模型。热交换器通常有四种典型的管阵几何排列形式，对于 所有的几何排列形式，解析模型的基本概念是相同的，排列形式只是影响模型中使用的试 验参数。对于平行三角形管束模型：
[0068] 2. 1不稳定流速
[0070] 式中Y，很多研宄已经证明流体弹性不稳定通常首先在升力方向（y方向）产生。
[0071] 假设柔性管在两自由度方向以相同的固有频率做简谐运动
[0074]其中，aj^t) *a2(sm，t)分别是流管1和流管2在最小间隙处的面积微扰（面 积变小为正，变大为负。），它们由管子几何参数定义。
[0076]方程的解为y(t) =Yeiut
[0077] 式中，w是振动复频率。Ur=UQ/?nlQ是缩减速度，1 Q= 2s…频率比?/ ?n，是管子在静水中的固有频率，Sl是从振动管到压力扰动可以忽略的位置的距离。
[0078] 由于扰动衰减和时间衰减仅发生在"流管"路径上，以下给出两组针对平行三角形 管束的面积扰动方程
[0082] f(s) = 1,sa^s^ss
[0083] 而管子运动和流体力之间的时间延迟通过相位函数来表示
[0085] 为了便于使用，在接下来的公式中将炉〇7,.)写为识
[0091]
[0095] 2. 2非稳态压力的计算
[0097] 各分项计算：
[0100] B2= (s1-sa)A4/2
[0101] B3= (s「sa)A5
[0102]
[0103] 将以上代入非稳态压力的积分计算表达式得：
[0105]其中：
[0106] 〇!= -torl〇UriB2
[0109] 2. 3作用在管子上的流体力
[0111] 计算各分项：
[0116] 将以上各项代入流体力的计算公式可得：
[0118] Fy=Fu+FJ是复数形式表示的流体力，将上式代入公式
[0121] 中就可得到折合流速
?和质量阻尼比
之间的关系，从而得到稳定区图。对 于一个确定的流管模型，质量阻尼比是确定，只需要计算出具体折合流速，然后判断两者确 定的点，若在稳定区图内，则说明蒸汽发生器管束不会发生流弹失稳，若两者确定的点在稳 定区图外，则说明蒸汽发生器管束将发生流弹失稳。值得注意的是以上推导过程是根据平 行三角形管阵推导而来，对于旋转三角形管阵和旋转方阵，流体力具有相同的形式，只需将 上式流体力除以2c〇s(a)即可。
[0122] 本发明的有益效果是：具有显式表达、输入试验参数少、不依赖数值积分的优点。
【附图说明】
[0123] 图1为四种典型的管阵排列形式；
[0124] 图2为流管模型；
[0125] 图3为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出正方形管阵（p/d= 1. 25)的稳定区图；
[0126] 图4为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出旋转三角形管阵（p/d= 1. 25)的稳定区图；
[0127]图5为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出旋转方正方形管阵（p/d= 1. 25)的稳定区图；
[0128]图6为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出三角形管阵（p/d= 1. 25)的稳定区图。
【具体实施方式】
[0129] 为了使本发明的目的，技术方案及优点更加清晰，下面以聂清德等人风洞实验为 实施例结合附图对本发明作进一步的描述。
[0130] 本发明提供的一种蒸汽发生器管束两相横向流弹失稳分析方法，其中主要包括以 下步骤：1确定管子的基本参数，见表1、表2。
[0131] 表1各种管束的试验数据
[0133]注：m= 1. 24kg/m管径d= 0? 022m
[0134] 表2实验管束的几何参数
[0135]
[0136] 2确定两相横向流弹失稳参数确定
[0137] 2. 1两相横向流弹失稳试验参数确定
[0138] 唯一依赖于试验输入的参数为开始产生附着和分离点的位置，这可以通过Weaver 和AbdRabbo以及Scott的可视化研宄进行估计。四种典型的管阵排列的可以归为以下两 组：
[0140]注：0为流体附着表面的法向矢量角度
[0141] 2. 2两相横向流弹失稳计算参数确定
[0142] 2.2. 1流场计算参数
[0143] (1)空泡份额
[0144] 根据均相流模型（HEM)假设，两相流的两相充分混合，并且密度和温度分布均匀， 不存在气相和液相的相对滑移，可将空泡份额定义为：
[0146] 式中，\为真实含汽率，P8和Pi分别为气相和液相的密度。
[0147] (2)流动参数
[0148] 节距流速Up:Up=U"PAP-D)
[0149] 均相流密度P : P=PJl-aH) +pgah
[0150] 自由来流速度U"
[0151] 自由来流质量流速士
[0152] 节距质量流速
[0153] 式中，U"是自由来流速度（即假设不存在管子时的速度），P是管中心距，D是管 子直径（对于翅片管，用等效水力直径Dh代替）。A为自由来流的流道截面积；Pg为气相 密度、Pi为气相密度；Vg为气相速度、Vi为气相速度。
[0154] 2. 2. 2动水参数
[0155] (1)水力直径
[0156] 对于非翅片管，水力直径就简单的为管子外径D。对于翅片管，水力直径等效为 Dh:
[0157] Dh=Dr+RF (D0-Dr)
[0158] 式中，D,为齿根圆直径，D。为翅片外径，比值RF为翅片占据的体积与管子外表面到 翅片外径表面的有效体积的比。对于翅片管，Dh要用于振动激励参数的计算，等效节径比P/ D也要变为P/Dh。
[0159] (2)水动力附加质量
[0160] 管束中一根管子单位长度的水动力附加质量可以表不为：
[0162] 式中，队是周围管子的等效直径，比值De/D是管束对管子约束程度的度量。
[0163] 对于三角形管阵：
[0164] De/D= (0. 96+0. 5P/D)P/D
[0165] 对于正方形管阵：
[0166] De/D= (1. 07+0. 56P/D)P/D
[0167] 单位长度管子总的动力质量m包括水动力附加质量％、管子本身的质量mt和管内 流体的质量mi:
[0168] m=mh+mt+nii
[0169] (3)阻尼
[0170] 两相流中多跨热交换器管的阻尼包括支撑阻尼Gs、粘性阻尼^和两相阻尼 GTP:
[0171] ^T= ^s+^y+^Tp
[0172] 通量很高、空泡份额很大的情况下，支撑就可能是干的，此时只有摩擦阻尼存在， 这种情况下，支撑阻尼与气体中的热交换器管阻尼基本相同：
[0173] ^T= ^s+^y+^Tp
[0174] 管和支撑之间存在液体时，支撑阻尼Q既包括压膜阻尼GSF又包括摩擦阻尼
[0176] 粘性阻尼：
[0178] 两相阻尼比：
[0180] 式中，Pi是管和支撑间液体的密度，De是周围管子的等效直径，D是管子的外径， m是单位长度管子的总质量，N是折流板的数目，包括用两相流均相密度计算的水动力附加 质量，lm是管子的特征长度，L是支承板宽度，f是管子的频率，PTP是气液两相等效密度， UTP是等效运动粘度，UTP=ui/lil+ajjOi/Ug-l)]，〇8和ui分别为液相和气相的运动 黏度，而f(aH)的确定
[0181]f(aH) = ah/40 aH<40%
[0182]f(aH) = 1 40%彡aH彡 70 %。
[0183]f(aH) = l-(ah-70)/30 aH>70%
[0184] 3建立并求解两相横向流弹失稳的解析模型
[0185] 聂清德等人风洞实验：风洞宽0. 25m，高0. 35m。试验段中的自由流速可达40m/s， 其上游的紊流度接近1 %。放在试验段中的管束分别按正方形、转角三角形、转角正方形、三 角形四种典型的方式布置，节径比P/d为1. 25与1. 41。管子材料为T2紫铜，管长0. 39m， 直径0. 022m。上端以螺纹连接的方法与上管板固定后呈自由悬挂状，使管子的一阶固有频 率为90Hz。在距离管子的固定端20mm处，沿圆周且相应于流向以及与流动垂直的方向各粘 贴两片箔式应变片，采用半桥双臂互补法联接到动态应变仪的电桥盒上。管子振动时，应变 片电阻变化反映了管子自由端位移的大小，由此得到流速与管子振幅的关系曲线。相应于 振幅陡然急剧增大时的流速，即为发生流体弹性不稳定时的临界流速V。，实验数据参见表 1、表 2〇
[0186] 在Liming提出的单相横向流弹失稳显式分析模型-"流管"模型的基础上，提出 针对两相横向流弹失稳的解析模型。热交换器通常有四种典型的管阵几何排列形式，对于 所有的几何排列形式，解析模型的基本概念是相同的，排列形式只是影响模型中使用的试 验参数。对于平行三角形管束模型：
[0187] 3.1不稳定流速
[0189] 式中Y，很多研宄已经证明流体弹性不稳定通常首先在升力方向（y方向）产生。
[0190] 假设柔性管在两自由度方向以相同的固有频率做简谐运动
[0193]其中，aj^t) *a2(sm，t)分别是流管1和流管2在最小间隙处的面积微扰（面 积变小为正，变大为负。），它们由管子几何参数定义。
[0195]方程的解为y(t) =Yeiut
[0196] 其中w是振动复频率。14=UQ/?nlQ是缩减速度，1 Q= 2s^频率比'=?/? n， ?n是管子在静水中的固有频率，si是从振动管到压力扰动可以忽略的位置的距离。
[0197] 由于扰动衰减和时间衰减仅发生在"流管"路径上，以下给出两组针对平行三角形 管束的面积扰动方程
[0201]f(s) = 1,sa^s^ss
[0202] 而管子运动和流体力之间的时间延迟通过相位函数来表示
[0204] 为了便于使用，在接下来的公式中将写为供
[0207]
[0211]其中：
[0214] 3. 2非稳态压力的计算
[0216] 各分项计算：
[0218]
[0219] B2= (s「sa)A4/2
[0220] B3=(Sl-sa)A5
[0222] 将以上代入非稳态压力的积分计算表达式得：
[0224] 其中：
[0225] 〇!= -?rl〇UriB2
[0228] 3. 3作用在管子上的流体力
[0230] 计算各分项：
[0231] 为便于计算令
[0232]
[0235] 将以上各项代入流体力的计算公式可得：
[0237]Fy=Fh+FJ是复数形式表示的流体力，将上式代入公式
[0240] 中就可得到折合流速
和质量阻尼比
之间的关系，及可得到的折合流速 与质量阻尼比的关系曲线得到管束的稳定区图。对于一个确定的流管模型，质量阻尼比是 确定，只需要计算出具体折合流速，然后判断两者确定的点，若在稳定区图内，则说明蒸汽 发生器管束不会发生流弹失稳，若两者确定的点在稳定区图外，则说明蒸汽发生器管束将 发生流弹失稳。值得注意的是以上推导过程是根据平行三角形管阵推导而来，对于旋转三 角形管阵和旋转方阵，流体力具有相同的形式，只需将上式流体力除以2c〇S(a)即可。 [0241] 根据以上横向流弹失稳的原理编制计算程序，并根据聂清德等人的试验参数绘制 出各种管束的稳定区图。由图可见，对常见的四种管束，显式解析模型方法的结果与试验结 果基本一致。当mS/pD2>5,基于当前方法的结果相比试验结果更保守。当mS/pD2〈5,计 算结果相对试验结果偏大，但是两者之间的误差并不是很显著。结果表明，该显式解析模型 继承了Y-W"流管"模型的优点，即具有较好的计算精度及需要较少的试验数据。
【主权项】
1. 一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在于，蒸汽发生器管束两 相横向流弹失稳分析方法包括：建立两相横向流弹失稳的显式解析模型；确定求解该解析 模型所需要的两相流计算参数；最后求解解析模型得到管束的稳定区图。2. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，建立两相横向流弹失稳的显式解析模型中，为了考虑从中心管向外的面积微扰，采用以 下线性衰减函数f (s) - I, Sa^ s ^ s s 式中，曲线坐标S用以表示流动路径；足够远的上游位置为S1;流管在S a位置开始附着 于柔性管，并在Ss位置分离。3. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，建立两相横向流弹失稳的显式解析模型中，扰动衰减和时间衰减仅发生在"流管"路径 上，以下给出两组针对平行三角形管束的面积扰动方程式中，f (s)为线性衰减函数，-管子运动和流体力之间的时间延迟通 过相位函数，:为扰动衰减和时间衰减的相位函数，为扰动衰减的相位函数。4. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，确定求解该解析模型所需要的两相流计算参数，根据均相流模型（HEM)假设，两相流的 两相充分混合，并且密度和温度分布均匀，不存在气相和液相的相对滑移，可将空泡份额定 义为：式中，Xa为真实含气率，p 8和p i分别为气相和液相的密度。5. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，确定求解该解析模型所需要的两相流计算参数，两相流中多跨热交换器管的阻尼包括 支撑阻尼Gs、粘性阻尼ζν和两相阻尼ζ TP: ζ X= ζ S+ ζ V+ ζ TP 通量很高、空泡份额很大的情况下，支撑就可能是干的，此时只有摩擦阻尼存在，这种 情况下，支撑阻尼与气体中的热交换器管阻尼基本相同： ζ X= ζ S+ ζ V+ ζ TP 管和支撑之间存在液体时，支撑阻尼包括压膜阻尼ζ SF又包括摩擦阻尼ζ F: 粘性阻尼：两相阻尼比： 式中，P i是管和支撑间液体的密度，De是周围管子的等效直径，D是管子的外径，m是单 位长度管子的总质量，N是折流板的数目，包括用两相流均相密度计算的水动力附加质量， Im是管子的特征长度，L是支承板宽度，f是管子的频率，P TP是气液两相等效密度，υ τρ是 等效运动粘度，uTP= υ ^[l+ajjO^Ug-l)]，〇8和υ i分别为液相和气相的运动黏度， 而f ( a H)的确定 f(aH) = aH/40 aH<40% f(a H) = I 40% 彡 a H彡 70%。 f (a H) = l-(a h-70)/30 aH>70%。
【专利摘要】本发明涉及一种蒸汽发生器管束两相横向流弹失稳分析方法。本发明是为了克服核反应堆蒸汽发生器换热管流弹失稳而导致的管束破坏问题，而提出一种蒸汽发生器管束两相横向流弹失稳分析方法。该方法是在Liming提出的单相横向流弹失稳显示分析模型-“流管”模型的基础上，提出针对两相横向流弹失稳的显式解析模型，然后确定求解该解析模型所需要的两相流计算参数，最后求解解析模型得到蒸汽发生器管束临界折合流速和质量阻尼比，从而得到稳定区图。一个确定的流管模型，质量阻尼比是确定，只需要计算出具体折合流速，然后判断两者确定的点，若在稳定区图内，则说明蒸汽发生器管束不会发生流弹失稳，反之，则说明蒸汽发生器管束将发生流弹失稳。
【IPC分类】G06F17/50
【公开号】CN104899387
【申请号】CN201510335090
【发明人】柳军, 王伟, 桑果, 田晨, 郭小强, 何星, 黄陈
【申请人】西南石油大学
【公开日】2015年9月9日
【申请日】2015年6月17日
【技术领域】
[0001] 本发明涉及一种蒸汽发生器管束两相横向流弹失稳分析方法，提出的两相横向流 弹失稳的解析模型，分析核反应堆中蒸汽发生器管束两相横向流弹失稳是否发生的一种方 法。
【背景技术】
[0002] 流弹失稳机理与航天工业中广泛研宄的颤振现象相同，不稳定性的主要参数是系 统的阻尼和流体弹性力。横向流作用下的管束，当流速增加到一定值时，系统中流体弹性力 所做的功将大于系统阻尼散耗的功，因此管子振动的振幅将急剧增大，这时就称出现了流 弹失稳。
[0003] 从理论上看，流体弹性力应该通过求解一组联立的流体力学和结构动力学方程得 到。但由于我们涉及的问题通常十分复杂，如大量的换热管，复杂的边界条件和高雷诺数的 流动等，使得从理论上进行求解十分困难，甚至是不可能的。通常的做法是根据大量的实验 结果，提出简化的数学物理模型。在这些模型中，流体力一般不再以独立方程的形式出现而 是表示成与结构的位移、速度和加速度等有关的函数，这种关系体现了流体作用使结构产 生位移和变形，而结构的位移和变形又改变了流体作用这样一个流固耦合振动的事实。
[0004]目前，已有分析管阵流弹失稳的模型，主要为拟静力（如Connors模型）、拟定常 (如Paidoussis模型）、非定常流模型（如Chen模型）及解析模型（Yetisir流管模型）。 分析管阵临界流速时，拟静力、拟定常和非定常流模型的流体弹性力都不同程度地依赖于 实验数据，且这些模型中的一些参数是在一定条件下、针对某个管阵形式得到的，有一定的 保守性和不确定性，并不能对所有实际问题给出可靠的临界流速。而解析模型则是根据大 量的实验结果，对流体弹性系统作一些简化，用比较成熟的数学、力学方法，得到流体弹性 力。尽管这类模型与实际情况还有一定差距，有待于进一步完善，但无疑是一类具有普遍意 义的数学模型。
[0005] 现有的解析方法计算过程比较繁琐，甚至需用到数值计算，且都是针对单相流弹 失稳的分析，如果用其分析两相流弹失稳问题，还需将气-液两相流等效为单相流，这一方 法的有效性也有待考证。
【发明内容】
[0006] 本发明的目的是：本发明是为了克服核反应堆蒸汽发生器换热管流弹失稳而导 致的管束破坏问题，而提出一种蒸汽发生器管束两相横向流弹失稳分析方法。包括以下步 骤：
[0007] 1两相横向流弹失稳参数确定
[0008] 热交换器通常有四种典型的管阵几何排列形式，对于所有的几何排列形式，解析 模型的基本概念是相同的，排列形式只是影响模型中使用的试验参数。具体以平行三角形 管束为例，而对于旋转三角形和旋转正方形管束（P，l. 7)，流体力Fy具有同样的形式，除了 它需要除以2c〇Sa，这是由于所定义的最小间距扰动面积函数a(sm，t)具有不同的形式。 [0009] 1.1基本定义
[0010] 引入"流管"概念，平行三角形管阵流管模型。曲线坐标S用以表示流动路径，流 管的平均横截面积在整个流管上假设为常数，且由-Sl位置的最小间隙面积A^确定。
[0011] A〇=min(Pcosa-D/2,P~D)
[0012] 瞬时流管面积A(s，t)会随着管子的振动而改变，A(s，t)包含两项：平均项怂，和 波动项a(s,t)。
[0013]A(s,t) =A〇+a(s,t)
[0014] 同样的，速度和压力定义为：
[0015] U(s,t) =U〇+u(s,t),P(s,t) =P〇+p(s,t)
[0016] 然而，在足够远的上游位置Si，扰动可以忽略，速度和压力可以分别视为常数UQ和 匕，流管假设在^位置开始附着于柔性管，并在ss位置分离。
[0017] 管子的振动不能立刻影响到流管的其它位置，因此，采用一个相函数供(.〇来考虑 扰动延迟效应。上游微扰函数a(s，t)可以表示为
[0019] f(s)为人工衰减函数，用以表示从柔性管的位置到上游点的微扰衰减。a(sm，t)是 最小间隙位置的面积微扰，是管子几何的函数。
[0020] 1. 2两相横向流弹失稳试验参数确定
[0021] 唯一依赖于试验输入的参数为开始产生附着和分离点的位置，这可以通过Weaver 和AbdRabbo以及Scott的可视化研宄进行估计。四种典型的管阵排列的可以归为以下两 组：
[0023] 注为流体附着表面的法向矢量角度
[0024] 1. 3两相横向流弹失稳计算参数确定
[0025] 1. 3. 1流场计算参数
[0026] (1)空泡份额
[0027] 根据均相流模型（HEM)假设，两相流的两相充分混合，并且密度和温度分布均匀， 不存在气相和液相的相对滑移，可将空泡份额定义为：
[0029] 上式中，\为真实含汽率，P8和P1分别为气相和液相的密度。
[0030] (2)流动参数
[0031] 节距流速％^，U"PAP-D)
[0032] 均相流密度P:P=Pid-ad+PgaH
[0033] 自由来流速度U":
[0034] 自由来流质量流速
[0035] 节距质量流速:
[0036] 式中，U"是自由来流速度（即假设不存在管子时的速度），P是管中心距，D是管 子直径（对于翅片管，用等效水力直径Dh代替），A为自由来流的流道截面积，Pg为气相 密度、Pi为气相密度，Vg为气相速度、Vi为气相速度。
[0037] 1. 3. 2动水参数
[0038] (1)水力直径
[0039] 对于非翅片管，水力直径就简单的为管子外径D。对于翅片管，水力直径等效为 Dh:
[0040] Dh=Dr+RF (D0-Dr)
[0041] 式中，D,为齿根圆直径，D。为翅片外径，比值RF为翅片占据的体积与管子外表面到 翅片外径表面的有效体积的比。对于翅片管，Dh要用于振动激励参数的计算，等效节径比P/ D也要变为P/Dh。
[0042] (2)水动力附加质量
[0043] 管束中一根管子单位长度的水动力附加质量可以表不为：
[0045] 式中，队是周围管子的等效直径，比值De/D是管束对管子约束程度的度量。
[0046] 对于三角形管阵：
[0047] De/D= (0? 96+0. 5P/D)P/D
[0048] 对于正方形管阵：
[0049] De/D= (1. 07+0. 56P/D)P/D
[0050] 单位长度管子总的动力质量m包括水动力附加质量％、管子本身的质量mt和管内 流体的质量mi:
[0051] m=mh+m^nii
[0052] (3)阻尼
[0053] 两相流中多跨热交换器管的阻尼包括支撑阻尼Gs、粘性阻尼^和两相阻尼 GTP:
[0054] ^T= ^ s+^y+^Tp
[0055] 通量很高、空泡份额很大的情况下，支撑就可能是干的，此时只有摩擦阻尼存在， 这种情况下，支撑阻尼与气体中的热交换器管阻尼基本相同：
[0056] ^ T= ^ s+ ^ v+ ^ TP
[0057] 管和支撑之间存在液体时，支撑阻尼Q既包括压膜阻尼GSF又包括摩擦阻尼
[0059]粘性阻尼：
[0061]两相阻尼比：
[0063] 式中，Pi是管和支撑间液体的密度，De是周围管子的等效直径，D是管子的外径， m是单位长度管子的总质量，包括用两相流均相密度计算的水动力附加质量，lm是管子的特 征长度，L是支承板宽度，f?是管子的频率，PTP是气液两相等效密度，uTP是等效运动粘度， uTP=ui/El+ajjOi/Ug-l)]，〇8和u身别为液相和气相的运动黏度，而f(aH)的确定
[0064] f ( a H) = a h/40 a H<40%
[0065] f ( a H) = 1 40%彡aH彡 70%。
[0066] f (a H) = l-(ah-70)/30 aH>70%
[0067] 2在Liming 提出的单相横向流弹失稳显式分析模型-"流管"模型的基础上，提出 针对两相横向流弹失稳的解析模型。热交换器通常有四种典型的管阵几何排列形式，对于 所有的几何排列形式，解析模型的基本概念是相同的，排列形式只是影响模型中使用的试 验参数。对于平行三角形管束模型：
[0068] 2. 1不稳定流速
[0070] 式中Y，很多研宄已经证明流体弹性不稳定通常首先在升力方向（y方向）产生。
[0071] 假设柔性管在两自由度方向以相同的固有频率做简谐运动
[0074]其中，aj^t) *a2(sm，t)分别是流管1和流管2在最小间隙处的面积微扰（面 积变小为正，变大为负。），它们由管子几何参数定义。
[0076]方程的解为y(t) =Yeiut
[0077] 式中，w是振动复频率。Ur=UQ/?nlQ是缩减速度，1 Q= 2s…频率比?/ ?n，是管子在静水中的固有频率，Sl是从振动管到压力扰动可以忽略的位置的距离。
[0078] 由于扰动衰减和时间衰减仅发生在"流管"路径上，以下给出两组针对平行三角形 管束的面积扰动方程
[0082] f(s) = 1,sa^s^ss
[0083] 而管子运动和流体力之间的时间延迟通过相位函数来表示
[0085] 为了便于使用，在接下来的公式中将炉〇7,.)写为识
[0091]
[0095] 2. 2非稳态压力的计算
[0097] 各分项计算：
[0100] B2= (s1-sa)A4/2
[0101] B3= (s「sa)A5
[0102]
[0103] 将以上代入非稳态压力的积分计算表达式得：
[0105]其中：
[0106] 〇!= -torl〇UriB2
[0109] 2. 3作用在管子上的流体力
[0111] 计算各分项：
[0116] 将以上各项代入流体力的计算公式可得：
[0118] Fy=Fu+FJ是复数形式表示的流体力，将上式代入公式
[0121] 中就可得到折合流速
?和质量阻尼比
之间的关系，从而得到稳定区图。对 于一个确定的流管模型，质量阻尼比是确定，只需要计算出具体折合流速，然后判断两者确 定的点，若在稳定区图内，则说明蒸汽发生器管束不会发生流弹失稳，若两者确定的点在稳 定区图外，则说明蒸汽发生器管束将发生流弹失稳。值得注意的是以上推导过程是根据平 行三角形管阵推导而来，对于旋转三角形管阵和旋转方阵，流体力具有相同的形式，只需将 上式流体力除以2c〇s(a)即可。
[0122] 本发明的有益效果是：具有显式表达、输入试验参数少、不依赖数值积分的优点。
【附图说明】
[0123] 图1为四种典型的管阵排列形式；
[0124] 图2为流管模型；
[0125] 图3为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出正方形管阵（p/d= 1. 25)的稳定区图；
[0126] 图4为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出旋转三角形管阵（p/d= 1. 25)的稳定区图；
[0127]图5为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出旋转方正方形管阵（p/d= 1. 25)的稳定区图；
[0128]图6为根据以上横向流弹失稳的原理编制计算程序和聂清德等人的试验参数绘 制出三角形管阵（p/d= 1. 25)的稳定区图。
【具体实施方式】
[0129] 为了使本发明的目的，技术方案及优点更加清晰，下面以聂清德等人风洞实验为 实施例结合附图对本发明作进一步的描述。
[0130] 本发明提供的一种蒸汽发生器管束两相横向流弹失稳分析方法，其中主要包括以 下步骤：1确定管子的基本参数，见表1、表2。
[0131] 表1各种管束的试验数据
[0133]注：m= 1. 24kg/m管径d= 0? 022m
[0134] 表2实验管束的几何参数
[0135]
[0136] 2确定两相横向流弹失稳参数确定
[0137] 2. 1两相横向流弹失稳试验参数确定
[0138] 唯一依赖于试验输入的参数为开始产生附着和分离点的位置，这可以通过Weaver 和AbdRabbo以及Scott的可视化研宄进行估计。四种典型的管阵排列的可以归为以下两 组：
[0140]注：0为流体附着表面的法向矢量角度
[0141] 2. 2两相横向流弹失稳计算参数确定
[0142] 2.2. 1流场计算参数
[0143] (1)空泡份额
[0144] 根据均相流模型（HEM)假设，两相流的两相充分混合，并且密度和温度分布均匀， 不存在气相和液相的相对滑移，可将空泡份额定义为：
[0146] 式中，\为真实含汽率，P8和Pi分别为气相和液相的密度。
[0147] (2)流动参数
[0148] 节距流速Up:Up=U"PAP-D)
[0149] 均相流密度P : P=PJl-aH) +pgah
[0150] 自由来流速度U"
[0151] 自由来流质量流速士
[0152] 节距质量流速
[0153] 式中，U"是自由来流速度（即假设不存在管子时的速度），P是管中心距，D是管 子直径（对于翅片管，用等效水力直径Dh代替）。A为自由来流的流道截面积；Pg为气相 密度、Pi为气相密度；Vg为气相速度、Vi为气相速度。
[0154] 2. 2. 2动水参数
[0155] (1)水力直径
[0156] 对于非翅片管，水力直径就简单的为管子外径D。对于翅片管，水力直径等效为 Dh:
[0157] Dh=Dr+RF (D0-Dr)
[0158] 式中，D,为齿根圆直径，D。为翅片外径，比值RF为翅片占据的体积与管子外表面到 翅片外径表面的有效体积的比。对于翅片管，Dh要用于振动激励参数的计算，等效节径比P/ D也要变为P/Dh。
[0159] (2)水动力附加质量
[0160] 管束中一根管子单位长度的水动力附加质量可以表不为：
[0162] 式中，队是周围管子的等效直径，比值De/D是管束对管子约束程度的度量。
[0163] 对于三角形管阵：
[0164] De/D= (0. 96+0. 5P/D)P/D
[0165] 对于正方形管阵：
[0166] De/D= (1. 07+0. 56P/D)P/D
[0167] 单位长度管子总的动力质量m包括水动力附加质量％、管子本身的质量mt和管内 流体的质量mi:
[0168] m=mh+mt+nii
[0169] (3)阻尼
[0170] 两相流中多跨热交换器管的阻尼包括支撑阻尼Gs、粘性阻尼^和两相阻尼 GTP:
[0171] ^T= ^s+^y+^Tp
[0172] 通量很高、空泡份额很大的情况下，支撑就可能是干的，此时只有摩擦阻尼存在， 这种情况下，支撑阻尼与气体中的热交换器管阻尼基本相同：
[0173] ^T= ^s+^y+^Tp
[0174] 管和支撑之间存在液体时，支撑阻尼Q既包括压膜阻尼GSF又包括摩擦阻尼
[0176] 粘性阻尼：
[0178] 两相阻尼比：
[0180] 式中，Pi是管和支撑间液体的密度，De是周围管子的等效直径，D是管子的外径， m是单位长度管子的总质量，N是折流板的数目，包括用两相流均相密度计算的水动力附加 质量，lm是管子的特征长度，L是支承板宽度，f是管子的频率，PTP是气液两相等效密度， UTP是等效运动粘度，UTP=ui/lil+ajjOi/Ug-l)]，〇8和ui分别为液相和气相的运动 黏度，而f(aH)的确定
[0181]f(aH) = ah/40 aH<40%
[0182]f(aH) = 1 40%彡aH彡 70 %。
[0183]f(aH) = l-(ah-70)/30 aH>70%
[0184] 3建立并求解两相横向流弹失稳的解析模型
[0185] 聂清德等人风洞实验：风洞宽0. 25m，高0. 35m。试验段中的自由流速可达40m/s， 其上游的紊流度接近1 %。放在试验段中的管束分别按正方形、转角三角形、转角正方形、三 角形四种典型的方式布置，节径比P/d为1. 25与1. 41。管子材料为T2紫铜，管长0. 39m， 直径0. 022m。上端以螺纹连接的方法与上管板固定后呈自由悬挂状，使管子的一阶固有频 率为90Hz。在距离管子的固定端20mm处，沿圆周且相应于流向以及与流动垂直的方向各粘 贴两片箔式应变片，采用半桥双臂互补法联接到动态应变仪的电桥盒上。管子振动时，应变 片电阻变化反映了管子自由端位移的大小，由此得到流速与管子振幅的关系曲线。相应于 振幅陡然急剧增大时的流速，即为发生流体弹性不稳定时的临界流速V。，实验数据参见表 1、表 2〇
[0186] 在Liming提出的单相横向流弹失稳显式分析模型-"流管"模型的基础上，提出 针对两相横向流弹失稳的解析模型。热交换器通常有四种典型的管阵几何排列形式，对于 所有的几何排列形式，解析模型的基本概念是相同的，排列形式只是影响模型中使用的试 验参数。对于平行三角形管束模型：
[0187] 3.1不稳定流速
[0189] 式中Y，很多研宄已经证明流体弹性不稳定通常首先在升力方向（y方向）产生。
[0190] 假设柔性管在两自由度方向以相同的固有频率做简谐运动
[0193]其中，aj^t) *a2(sm，t)分别是流管1和流管2在最小间隙处的面积微扰（面 积变小为正，变大为负。），它们由管子几何参数定义。
[0195]方程的解为y(t) =Yeiut
[0196] 其中w是振动复频率。14=UQ/?nlQ是缩减速度，1 Q= 2s^频率比'=?/? n， ?n是管子在静水中的固有频率，si是从振动管到压力扰动可以忽略的位置的距离。
[0197] 由于扰动衰减和时间衰减仅发生在"流管"路径上，以下给出两组针对平行三角形 管束的面积扰动方程
[0201]f(s) = 1,sa^s^ss
[0202] 而管子运动和流体力之间的时间延迟通过相位函数来表示
[0204] 为了便于使用，在接下来的公式中将写为供
[0207]
[0211]其中：
[0214] 3. 2非稳态压力的计算
[0216] 各分项计算：
[0218]
[0219] B2= (s「sa)A4/2
[0220] B3=(Sl-sa)A5
[0222] 将以上代入非稳态压力的积分计算表达式得：
[0224] 其中：
[0225] 〇!= -?rl〇UriB2
[0228] 3. 3作用在管子上的流体力
[0230] 计算各分项：
[0231] 为便于计算令
[0232]
[0235] 将以上各项代入流体力的计算公式可得：
[0237]Fy=Fh+FJ是复数形式表示的流体力，将上式代入公式
[0240] 中就可得到折合流速
和质量阻尼比
之间的关系，及可得到的折合流速 与质量阻尼比的关系曲线得到管束的稳定区图。对于一个确定的流管模型，质量阻尼比是 确定，只需要计算出具体折合流速，然后判断两者确定的点，若在稳定区图内，则说明蒸汽 发生器管束不会发生流弹失稳，若两者确定的点在稳定区图外，则说明蒸汽发生器管束将 发生流弹失稳。值得注意的是以上推导过程是根据平行三角形管阵推导而来，对于旋转三 角形管阵和旋转方阵，流体力具有相同的形式，只需将上式流体力除以2c〇S(a)即可。 [0241] 根据以上横向流弹失稳的原理编制计算程序，并根据聂清德等人的试验参数绘制 出各种管束的稳定区图。由图可见，对常见的四种管束，显式解析模型方法的结果与试验结 果基本一致。当mS/pD2>5,基于当前方法的结果相比试验结果更保守。当mS/pD2〈5,计 算结果相对试验结果偏大，但是两者之间的误差并不是很显著。结果表明，该显式解析模型 继承了Y-W"流管"模型的优点，即具有较好的计算精度及需要较少的试验数据。
【主权项】
1. 一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在于，蒸汽发生器管束两 相横向流弹失稳分析方法包括：建立两相横向流弹失稳的显式解析模型；确定求解该解析 模型所需要的两相流计算参数；最后求解解析模型得到管束的稳定区图。2. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，建立两相横向流弹失稳的显式解析模型中，为了考虑从中心管向外的面积微扰，采用以 下线性衰减函数f (s) - I, Sa^ s ^ s s 式中，曲线坐标S用以表示流动路径；足够远的上游位置为S1;流管在S a位置开始附着 于柔性管，并在Ss位置分离。3. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，建立两相横向流弹失稳的显式解析模型中，扰动衰减和时间衰减仅发生在"流管"路径 上，以下给出两组针对平行三角形管束的面积扰动方程式中，f (s)为线性衰减函数，-管子运动和流体力之间的时间延迟通 过相位函数，:为扰动衰减和时间衰减的相位函数，为扰动衰减的相位函数。4. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，确定求解该解析模型所需要的两相流计算参数，根据均相流模型（HEM)假设，两相流的 两相充分混合，并且密度和温度分布均匀，不存在气相和液相的相对滑移，可将空泡份额定 义为：式中，Xa为真实含气率，p 8和p i分别为气相和液相的密度。5. 如权利要求1所述的一种蒸汽发生器管束两相横向流弹失稳分析方法，其特征在 于，确定求解该解析模型所需要的两相流计算参数，两相流中多跨热交换器管的阻尼包括 支撑阻尼Gs、粘性阻尼ζν和两相阻尼ζ TP: ζ X= ζ S+ ζ V+ ζ TP 通量很高、空泡份额很大的情况下，支撑就可能是干的，此时只有摩擦阻尼存在，这种 情况下，支撑阻尼与气体中的热交换器管阻尼基本相同： ζ X= ζ S+ ζ V+ ζ TP 管和支撑之间存在液体时，支撑阻尼包括压膜阻尼ζ SF又包括摩擦阻尼ζ F: 粘性阻尼：两相阻尼比： 式中，P i是管和支撑间液体的密度，De是周围管子的等效直径，D是管子的外径，m是单 位长度管子的总质量，N是折流板的数目，包括用两相流均相密度计算的水动力附加质量， Im是管子的特征长度，L是支承板宽度，f是管子的频率，P TP是气液两相等效密度，υ τρ是 等效运动粘度，uTP= υ ^[l+ajjO^Ug-l)]，〇8和υ i分别为液相和气相的运动黏度， 而f ( a H)的确定 f(aH) = aH/40 aH<40% f(a H) = I 40% 彡 a H彡 70%。 f (a H) = l-(a h-70)/30 aH>70%。
【专利摘要】本发明涉及一种蒸汽发生器管束两相横向流弹失稳分析方法。本发明是为了克服核反应堆蒸汽发生器换热管流弹失稳而导致的管束破坏问题，而提出一种蒸汽发生器管束两相横向流弹失稳分析方法。该方法是在Liming提出的单相横向流弹失稳显示分析模型-“流管”模型的基础上，提出针对两相横向流弹失稳的显式解析模型，然后确定求解该解析模型所需要的两相流计算参数，最后求解解析模型得到蒸汽发生器管束临界折合流速和质量阻尼比，从而得到稳定区图。一个确定的流管模型，质量阻尼比是确定，只需要计算出具体折合流速，然后判断两者确定的点，若在稳定区图内，则说明蒸汽发生器管束不会发生流弹失稳，反之，则说明蒸汽发生器管束将发生流弹失稳。
【IPC分类】G06F17/50
【公开号】CN104899387
【申请号】CN201510335090
【发明人】柳军, 王伟, 桑果, 田晨, 郭小强, 何星, 黄陈
【申请人】西南石油大学
【公开日】2015年9月9日
【申请日】2015年6月17日
最新回复(0)