650 MHz双cell超导腔冷氦气自动降温方案研究

常正则，李 梅，韩瑞雄，孙良瑞，张祥镇，姜永诚，桑民敬，叶 瑞，李少鹏，葛 锐

(中国科学院 高能物理研究所，北京 100049)

超导腔是现代超导加速器的核心先进部件，可在单位长度上为束流提供更高的加速电压和高频功率，同时可节省设备占用的空间。超导腔多采用纯铌材料制成，使用时需将其安装在由液氮和液氦两层冷屏逐级热隔断的低温恒温器中，超导腔本体通常采用冷氦气和液氦冲刷浸泡的方式进行降温冷却。超导腔的降温过程具有以下特点。1) 超导腔是薄壁空腔结构，较易产生形变，无法承受过大的热应力和压力。2) 在不同的温度区间，超导腔的降温速率要求不同。120 K以上时要求降温平缓，表面温差小，避免热应力过大；120～40 K时要求降温速率加快，避免“氢中毒”现象产生[1]；40 K以下则可直接积液氦进行冷却。3) 超导腔形状复杂，有部分区域传感器无法测量，且传感器反应时间长，难以实现实时反馈。

国内外对超导腔降温过程的研究较少，王国平等[2]对BEPCⅡ超导腔垂直测试杜瓦的漏热进行了分析与实验研究，得到准确的漏热量和漏热区域。刘以勇[3]对超导波荡器进行了冷却测试，得到了各部分温度与时间的对应曲线。Gorbounov等[4]发展了一种结合数值和分析方法的Fortran计算程序，可用于计算单相氦冷却超导电缆的降温过程。Chang等[5]对超导腔运行过程中的热负荷进行了计算，并与实验得到的液氦蒸发量进行了对比。O’Rourke等[6]通过在低热负荷真空管内预冷却低温恒温管，可实现AIST超导加速器(SCA)低温恒温器在从室温(300 K)到4 K降温区间内液氦损失的蒸发量降为0。Wenskat等[7]对欧洲XFEL的AMTF垂直测试站的降温过程进行了研究，通过对贴在恒温器外壁面的3种不同TVO传感器和两个安装在顶部和底部的Cernox TM传感器进行监测，降温过程基于PLC可实现自动控制。Dhakal等[8]对Nb/Cu复合超导射频腔从室温冷却到1.6 K期间，超导腔的剩余磁场所受的影响进行了研究。Huang等[9]对大晶粒Nb超导射频腔的表面电阻在冷却过程中受温度梯度的影响进行了研究。

先进光源技术研发与测试平台(PAPS)中的650 MHz双cell超导腔采用波导型高阶模耦合器设计，可应用于高能环形正负电子对撞机(CEPC)、国际直线对撞机(ILC)和上海硬X射线自由电子激光等大型超导腔加速器项目中，因此对650 MHz双cell超导腔的自动降温进行研究具有重要意义。以上对超导腔的降温研究中，缺乏对其传热机理的分析研究，难以实现降温过程的自动化。因此，本文以PAPS中的650 MHz双cell超导腔为例，开展传热机理分析和数值仿真研究，建立完善的计算方法和计算程序，实现超导腔降温过程的自动化。

1 超导腔自动降温过程机理

1.1 超导腔结构模型

650 MHz双cell超导腔结构示意图如图1所示，其冷却结构主要由葫芦状的超导腔本身及外部包裹超导腔的外圆筒氦槽构成，超导腔的材质为纯铌(Nb)，氦槽的材质为工业纯钛(Ti-TA2)。

图1 650 MHz双cell超导腔结构示意图Fig.1 Scheme of 650 MHz double-cell superconducting cavity

650 MHz双cell超导腔降温过程的研究对象主要集中在超导腔、氦槽及其上游的输运管道上，忽略其他附属部件与它们之间的导热，这是由于超导腔内和氦槽通常处于真空度1×10-5Pa量级，基本可忽略空气的导热和对流热流。氦槽外部有液氮和液氦两层冷屏隔绝，其他部件与超导腔之间的接触面积较小，因此也可忽略接触面带来的导热热流。

1.2 自动降温的实现步骤

超导腔降温过程的基本思路为分段降温，即随降温时间的推移，逐步调整入口参数，避免腔体温差过大及降温速度过慢，最终实现整体平稳、快速、均匀的降温曲线。设每段的降温时间为1 h，即在该降温段内，入口参数保持不变，超导腔允许的最大温差为ΔTmax(K)，降温速率为Trate(K/h)，初始超导腔和氦槽壁面温度为Tinitial(K)，降温目标温度为Ttarget(K)，则总降温段数n为：

n=ceiling((Tinitial-Ttarget-ΔTmax)/Trate)+1

(1)

式中：ceiling为向上取整函数。

第1段冷氦气的入口温度为Tinitial-ΔTmax，可确保壁面上不会产生超过ΔTmax的温差。入口流量根据传热学的相关公式进行迭代计算得到恰当的来流流量，确保在第1个降温阶段(1 h)后，目标壁面温度恰好降至Tinitial-Trate。此时调整入口温度和入口流量，进入下一个降温阶段。第n段的冷氦气入口温度为Tinitial-ΔTmax-nTrate，同样迭代计算本段所需的流量，确保1 h后目标壁面温度降至Tinitial-nTrate，如此循环往复直至最后1个降温阶段。由于要将温度稳至目标温度，因此最后1段的入口温度即为Ttarget，以目标壁面温度降至Ttarget+0.1 K以内为整体降温过程结束的标志。降温时间不再是1 h，而是取决于上一个降温段结束时的目标壁面温度，因此最后1段降温时间τn=(Tn-1-Ttarget)/Trate，得到降温时间后，同样使用迭代计算方法得到最后1段降温过程所需的流量。

得到各降温阶段的腔前入口参数后，再根据上游流程管道的压降、掺混等公式进行反推，得到固定时刻兑温阀门的对应开度，最终通过得到的阀门开度时刻表实现降温过程的自动化。根据以上需求，本文开发了对应的Fortran程序Autocold V1.0，其中氦气的热物性通过实时调用NIST[10]的Refprop v9.0代码实现，对于金属热物性则通过实时调用Cryocomp V3.0代码实现。Autocold V1.0分为腔前参数计算和来流流程计算两部分。腔前参数计算通过输入初始温度、目标温度、最大允许温差、降温速率等参数对超导腔进行传热计算，得到合理的腔前流动参数。来流流程计算通过流程计算反推各降温阶段的兑温阀门开度参数。

1.3 腔前参数计算

1) 流场分析

采用CFD方法对超导腔流场进行稳态分析计算，对流场进行建模，两侧壁面设为绝热无滑移壁面，湍流模型采用k-ε模型[11]。采用ICEM对流场进行非结构化网格划分，经网格无关解验证后确定网格数量为160万个，求解器为Ansys CFX 14.0[12]。

入口边界设为质量流量入口，冷氦气流量为4 g/s，出口压力为101 325 Pa，计算不同入口温度下的截面流速分布，结果如图2所示。

由图2可知，随入口温度的降低，冷氦气流速不断减小。超导腔表面流速最快的区域在靠近入口处，表面流速最慢的区域在超导腔两侧的凹腔中。流速在入口下半部分分布不均匀，上半部分则较为均匀。氦槽表面流速的分布情况与之类似。

入口温度：a——260 K；b——220 K；c——180 K；d——140 K；e——100 K；f——60 K；g——40 K；h——20 K图2 不同入口温度下的截面流速分布Fig.2 Velocity distribution under condition of different inlet temperatures

2) 非稳态导热热流分析

考虑超导腔和氦槽在厚度方向(径向)的非稳态导热热流，厚度方向的毕渥数Bi为：

Bi=hδ/λ

(2)

式中：h为表面对流换热系数，W/(m2·K)；δ为厚度，m；λ为固体导热系数，W/(m·K)。

不同温度条件下两种材料(Nb和Ti-TA2)的导热系数如图3所示。可看出，在40 K温度以上时，两者的导热系数在20 W/(m·K)以上。因此Bimax=100×0.004/20=0.02<0.1，满足集总参数法的使用条件[13]，所以可忽略厚度方向上的温度变化。

图3 不同温度条件下Nb和Ti-TA2的导热系数Fig.3 Thermal conductivity factor of Nb and Ti-TA2 under condition of different temperatures

然后考虑轴向和周向的导热热流，在最极端的温差情况下，两侧凹腔表面流速为0，仅靠导热降温，同时冲击点处的对流换热系数无限大，冲击点的超导腔壁温等于来流冷氦气的温度，为定壁温边界。可将沿轴向(周向类似)的导热问题简化为第一类边界条件下的半无限大平板非稳态导热问题，其解析解为：

(3)

式中：T(x,τ)为某一指定τ时刻x处的温度；Tw为壁面温度；T0为初始壁温；α为热扩散率；η为函数自变量；erf为误差函数，为单调升函数。

以超导腔允许最大温差20 K、初始壁温为300 K、第1轮入口冷气温度为280 K计算，1 h后的腔体最高温度为287.4 K，下降12.6 K，因此设定的降温速率Trate不应大于该值。

3) 对流热流分析

超导腔降温过程可视为管内强迫对流过程，对于管内湍流，可使用Gnielinski公式计算平均努赛尔数Nuf[14]：

(1+(d/l)2/3)ct

f=(1.82lgRe-1.64)-2

(4)

式中：下标f表示流体；Re为雷诺数；Pr为普朗特数；l为管长，m；d为管内径，m；ct为热流方向上的温度校正项；f为管内湍流流动的Darcy阻力系数，按照Filonenko公式[14]计算。

对气体有：

ct=(Tf/Tw)0.45Tf/Tw=0.5～1.5

(5)

式中，Tf为流体温度，K。

对于管内层流，采用Sieder-Tate公式计算平均Nuf：

(6)

确定管内流动的平均Nuf后，则平均表面对流换热系数h为：

h=Nu×λ/d

(7)

分别计算氦槽侧壁面和超导腔侧壁面的平均对流换热系数，然后迭代求解得到冷氦气的出口温度。

(Tout-Tin)cpm=hqSq(Tq-Tf)+

hwSw(Tw-Tf)

(8)

式中：Tout和Tin分别为出、入口氦气温度；cp为比定压热容；m为质量流量；hq为超导腔侧的对流换热系数；Sq为超导腔侧的换热面积；Tq为超导腔温度；hw为氦槽侧壁面温度；Sw为氦槽侧换热面积。

腔前参数计算流程如图4所示。根据计算流程，可计算得到每段降温过程所需的腔前参数，为下一步来流流程计算提供了边界条件。

1.4 来流流程计算

图4 腔前参数计算流程Fig.4 Process diagram of cavity parameter

650 MHz双cell超导腔上游降温管路的2 500 W@4.5 K氦制冷机中，有300 K和40 K两个兑温阀门V1和V2，兑温后经一系列管路和弯头引至主分配阀箱，到达超导腔及氦槽位置。采用文献[15]的方法计算水平圆形直管流动摩擦阻力Δp，其中Δp=Δpf+Δpj，Δpf为直管摩擦阻力损失，Δpj为局部摩擦阻力损失。

直管摩擦阻力损失为：

(9)

式中：ρ为流体的密度；u为流体的速度。

局部摩擦阻力损失的计算采用阻力系数法，将局部能量损失表示为流体动能因子u2/2的1个倍数，即：

(10)

式中，ξ为局部阻力系数。管件与阀门的阻力系数需由实验确定，或查阅相关水力学手册[15]。

管线静态漏热及压力损失均会转化为内能，引起温度变化，若出口温度T2已知，则可求出入口温度T1，即：

(11)

(12)

式中：Q为热量，W；m为流体质量流量，kg/s。

对于竖直圆形直管阻力，计算方法与水平管类似，不同之处在于需考虑流动方向，相应地增加重力项，本工况中气流密度较小，且竖直管长度较短，因此可忽略重力项的影响，按水平管计算即可。阀门开度与流量对应关系取决于阀门自身的特性曲线，可查阅相关阀门的产品手册[16]。根据计算流程，可反推得到来流掺混所需的兑温路流量的配比关系及兑温阀门开度等参数。

2 超导腔自动降温的计算结果

2.1 腔前参数计算结果

首先计算650 MHz双cell超导腔冷却过程中最常见的一种工况(即工况1)。设超导腔允许的最大温差ΔTmax=20 K，设定的降温速率Trate=10 K/h，初始超导腔/氦槽壁面温度Tinitial=300 K，降温目标温度Ttarget=40 K，目标壁面为氦槽壁面。总降温时间为26 h，入口冷氦气的参数随降温段数的变化如图5所示。可看出，第1段降温过程流量略多，这是由于初始超导腔和氦槽的平均温度相同，之后的降温过程中，超导腔壁面的平均温度始终低于氦槽壁面的平均温度。入口流量在降温前14 h时变化不大，之后则大幅下降。这是由于金属热容在160 K以下会大幅下降，而氦气的热容却变化不大，因此此时降低同等温度所需的冷氦气流量会大幅下降。最后1段降温过程出现流量增加的现象，这是由于前n-1个降温段均是采用过余温差进行降温，而最后1段需以Ttarget的入口温度将目标壁面降温至Ttarget，其温差较小，导致的对流换热量也较小，因此需更大的流量。

图5 工况1下入口冷氦气参数随降温段数的变化Fig.5 Inlet cold helium parameter vs cooling section number at case 1

按照不同温区所需降温速率不同，计算另1组典型降温工况(即工况2)。在300～120 K的降温范围内设定降温速率Trate=10 K/h，超导腔允许的最大温差ΔTmax=20 K，为了使超导腔快速通过氢中毒温区，减少对超导腔品质的损害，在120 K处稳定片刻，随后在120～40 K降温范围内进行快速降温，设定降温速率Trate=30 K/h，超导腔允许的最大温差ΔTmax=40 K，目标壁面为氦槽壁面。入口冷氦气的参数随降温段数的变化如图6所示(总降温时间为20.67 h)。

由图6可知，前16段的入口流量基本不变，这是由于120 K以上时，超导腔和氦槽金属壁面的热物性变化不大(主要是热容)，氦气单位质量的比热容变化也不大，因此每段所需的入口流量较为稳定。第17段入口流量大幅增加的原因和上一个工况中末段流量增加的原因相同。第18段流量较大有两个原因：1) 降温速率加大，单位时间内需吸收更多热量导致的流量上升；2) 与上一个工况中第1段流量增加的原因相同。第19段的流量较第18段小，这是由于金属材料热容随温度的降低而减小导致的，但比前17段高，这是由于降温速率加大，单位时间内金属材料所需带走的热量较之前略有增加引起的。第20段流量大幅升高的原因则是末段温差较小。

图6 工况2下入口冷氦气参数随降温段数的变化Fig.6 Inlet cold helium parameter vs cooling section number at case 2

2.2 来流流程计算结果

图7 超导腔入口温度与流量的配比关系Fig.7 Relationship between superconducting cavity inlet temperature and flow rate

根据流动条件选择来流管道尺寸为φ26.9 mm×1.65 mm。根据流程计算反推得到超导腔入口温度与流量的配比关系，如图7所示。

若V1和V2两个兑温阀门均采用等百分比的理想阀门，其开度与流量的关系如图8所示。

图8 理想等百分比阀门开度与流量的关系Fig.8 Opening percent of ideal equal percentage valve vs flow rate

将流程计算得到的流量配比结果，根据阀门开度与流量的对应关系，即可反推出阀门开度与腔前参数一一对应的关系，最终建立阀门开度的时刻表，从而可实现降温过程的自动化。V1和V2阀门开度与降温时间的关系如图9所示。

图9 V1和V2阀门开度与降温时间的关系Fig.9 Valve opening percent of V1 and V2 vs cooling time

3 结论

本文以PAPS中的650 MHz双cell超导腔的降温过程为主要研究对象，采用准则关系式法及计算流体力学方法构建了整体换热模型并编程计算，得到了降温过程中超导腔所需的冷氦气入口腔前参数。根据腔前参数反推流程，考虑整体布局和局部阻力的损失，计算了兑温阀门开度与腔前参数的对应关系。研究结果表明，通过正向的流动传热计算可实现降温过程的自动化，因此可用于难以通过反馈控制的场合。在300 K初始温度、40 K目标温度、降温速率10 K/h条件下，总降温时间为26 h。第1段所需流量最大，为11.8 kg/h，倒数第2段所需流量最小，为3.7 kg/h。整体降温速率保持不变，由于金属热容随温度降低而减小，降温过程中所需的总流量基本呈下降趋势。