003控制棒组件落棒时间与历程分析计算程序开发.doc

控制棒组件落棒时间与历程分析计算程序开发孙 磊 于建华2 魏永涛2 谷芳毓 李天勇（1：中国核动力研究设计院二所，2：四川大学力学所）摘要 控制棒落棒时间是保证核电站安全运行的重要参数之一，本文在已有研究成果的基础上对落棒的计算分析作了进一步的研究，其主要的方法是：采用直接约束法处理接触问题、用反力构造法确定碰撞力、用理论及数值分析完成流体阻力计算等。在此基础上，编制了控制棒落棒分析有限元程序CRFD，并针对秦山600MW反应堆控制棒驱动线给出了试算结果。关键词 控制棒组件落棒时间分析； 接触和碰撞； 有限元1 引言 反应堆控制棒的落棒时间是核电厂安全分析的重要参数之一，也是驱动线设计考核的重要指标之一。控制棒落棒系统结构复杂，落棒时间与其结构设计、制造和安装直接相关，并受到堆内流场、温度、压力和外载荷等的影响，控制棒下落时受到的主要作用力有：重力、浮力、流体阻力和摩擦力（由控制棒与导向组件间发生碰撞而产生）。其下落过程涉及到的力学机理有：控制棒的下落运动、控制棒与各导向部件的耦合振动和接触碰撞以及流体效应（流固耦合和流体阻力）等，使落棒时间分析计算有相当的难度。 由于控制棒落棒时间对核电站安全运行的重要性，故各核电设计国在其计算程序的开发上都做了大量工作12。从目前得到的资料看，国内外现有的落棒分析软件都存在某些不足，因此本文在现有国内外研究成果的基础上，对落棒行为作了进一步的研究，开发出控制棒落棒时间有限元计算程序CRFD。分析中假设：结构小变形；流体阻力计算和振动方程计算解耦，即流体阻力与结构的变形状态无关，而只与运动部件的下落位置和速度有关。 分析方法2.1 控制方程 首先将控制棒组件、导向筒和导向管等都视为弹性体建立相互耦合的横向振动方程3： （1）式中，部件的材料密度,Kg/m3；部件横截面积，m2;部件的梁抗弯刚度,N/m；粘性阻尼系数,Ns/m； 应变速率阻尼系数,Ns/m； 地震载荷,N； 流体附加质量, kg/m3； 在t时刻部件间碰撞力,N；部件横向变形,m； 部件沿轴向长度局部坐标,m。再建立控制棒组件的竖向运动方程： （2）式中，控制棒竖向运动位移,m；控制棒组件总质量，kg；重力加速度，m/s2； 摩擦系数； 流体阻力,N； 控制棒总长,m； 控制棒总数。 方程（1）和（2）通过碰撞力而相互耦合。2.2 流体固体耦合 导向管和导向筒内的流固耦合属于弱耦合机理（控制棒与导向管之间的缝隙和控制棒直径比大于1/10），可用附加质量来较准确地模拟。对于控制棒导向管形式的结构，其流体附加质量矩阵4为： （3）其中：； 流体密度, kg/m3； 环向模态数； 控制棒半径,m； 导向管内径,m。2.3 接触问题 发生接触的两个物体必须满足无穿透约束条件，本文采用直接约束法解决接触问题。即直接跟踪物体的运动轨迹，一旦探测出接触发生，就将接触所需的运动约束（法向无相对运动，切向可滑动）和节点力（法向压力和切向摩擦力）作为边界条件直接施加在发生的节点上。这种方法对接触的描述精度高，具有普遍的适应性，它不需要增加特殊的界面单元，也不涉及复杂的接触条件变化，且不增加系统的自由度数。2.4碰撞力学模型控制棒组件与导向组件之间的接触碰撞为协调的非完全弹性碰撞5、6，且碰撞的接触区域是未知的，本文采用力函数法6来确定碰撞力，其理论如下： 假定碰撞过程满足下述条件：1）碰撞过程分为变形发生阶段和变形恢复阶段，且在此二阶段的变形是对称的；2）碰撞过程满足能量守恒原理、动量守恒原理和冲量定理；3）力与变形的关系为线性，变形与时间的关系为非线性；4）采用恢复系数e（碰撞后与碰撞前的速度比）来模拟非完全弹性碰撞。 设结构的结点在时刻开始与其他物体发生碰撞，碰撞持续时间为（s），根据上述假设，碰撞力（N）的表达式为： （4） 式中，a0、a1和a2为待定系数，可由能量守恒、动量守恒和冲量定理推导得到。2.5 摩擦力 摩擦力的计算采用库仑摩擦模型。 当两物体在接触面上无滑动时：； （5) 当两物体在接触面上有滑动时：； (6)为摩擦力（N）， 为沿接触面法向的正压力（N），最大静摩擦系数，动摩擦系数。2.6 流体阻力 从流体力学理论可知，控制棒在流体中运动所受的作用力有三类：浮力、因物体加速运动而产生的流体阻力和绕流阻力，详细讨论见文献7。3 数值求解3.1 振动方程的有限元离散 将结构动力方程（1）式有限元离散化后有： （7） 式中、和分别为结构的质量矩阵、附加质量矩阵、阻尼矩阵和刚度矩阵； 为外力向量（包括地震载荷和流体脉动压力）；为时刻的碰撞力向量；、和分别为结构的加速度、速度和位移响应向量。和由单元质量和单元刚度矩阵组装得到。设与质量和刚度矩阵成正比，即，和为比例阻尼系数，可由振动理论计算得到。用Newmark-法直接积分法解方程（7）。3.2 接触判断采用直接约束法，实时跟踪运动物体，探测发生接触的位置。在实施的过程中，首先根据运动部件和导向部件各节点的位置建立可能发生接触的节点对，节点对分为两组：运动部件对导向部件（如节点B与节点D和C）和导向部件对运动部件（节点C与节点A和B），如图1示。再对每一节点对判断其是否接触或穿透，以节点对B与D和C为例，其判断方法如下：ACDAB导向部件运动部件图1 接触判断示意图 （8）式中，为节点B的横向位置，为节点D和C所在单元对应节点B水平位置点的横向位置，为接触容限值。如发生穿透则返回到上一时间步，减少，重新计算横向位移。4 程序流程 本文采用面向对象编程技术编制完成有限元计算程序CRFD。CRFD程序主要由两大部分组成：一是固体振动计算（包括横向振动、碰撞判别和碰撞力计算等）；二是流体阻力计算。在固体振动计算部分，采用二维梁单元建立有限元模型，能完成外载荷（力或位移）作用下带导向约束结构的振动和运动状态分析。在流体阻力计算模块中，可以分别采用8节点轴对称流体单元和理论及经验公式完成在圆形导管中运动的圆柱体所受的流体阻力的计算。4.1 CRFD程序功能简介 CRFD程序使用C+语言，采用面向对象编程方法，在C+Builder3平台上编制而成，程序由约4150行语句组成。 CRFD程序有以下主要特点： CRFD程序为二维有限元分析程序，基本单元为二维梁单元； 输入载荷的形式为：节点位移和节点力； 外载荷作用方式为任意多点的载荷时程输入； 程序输出的计算结果有：总下落时间、每一时间步各节点的横向振动位移、速度和加速度以及每一时间步控制棒组件的下落位置、下落速度、加速度、所受的流体总阻力、摩擦总阻力和各碰撞点及其碰撞力大小。 CRFD程序能完成反应堆控制棒组件在自重和有（或无）外载作用下从任意初始位置下落到导向管缓冲段底部的时间与历程分析。CRFD程序运行环境：PC电脑；Windows 98或 Windows2000操作系统；C+Builder3或以上的运行平台。4.2 CRFD程序流程CRFD程序的计算思路是：首先，将运动部件和导向部件都视为弹性体并进行有限元离散化，在外载荷的作用下计算出当前步的结构动态响应，在此基础上来判别是否发生接触和穿透，如发生接触则确定接触节点位置，如发生穿透则返回到上一时间步，并减少时间增量重新计算；在接触点上约束其法向位移，并计算碰撞力和摩擦力；碰撞力并入下一步的结构振动计算，摩擦力再和根据上一步的下落位置和速度计算出的流体阻力一起代入下落运动方程，计算出当前步的下落位置、速度和加速度。循环以上计算直至运动部件到达规定位置而完成整个计算，并得到运动部件到达规定位置的运动时间。程序流程见图2。5 算例在以上计算理论和方法的基础上，本文编制出控制棒落棒分析有限元程序CRFD，应用该程序对秦山二期600MWe核电站反应堆驱动线进行落棒时间的计算，驱动线的运动部件为：一束（24根）控制棒、星形架和驱动杆，导向部件为：燃料组件中的导向管、下部导向筒、上部导向筒和驱动杆外导管。计算的载荷模拟抗震试验的加载方式，即在整个驱动线结构上的四个不同点（燃料组件下管座、燃料组件上管座、下部导向筒和上部导向筒的交接处和上部导向筒的顶端）施加四个不同的位移地震载荷，计算控制棒从最高悬挂位置在自重作用下落到导向管缓冲段入口的历程和时间（）以及下落到导向管底部的时间（t5+t6），地震输入为SSE设计波形，落棒状况又分为静态落棒和SSE地震、1.414倍SSE地震及1.6倍SSE地震作用。在计算模型上，将控制棒组件和导向组件中的各个部件（控制棒束、驱动杆、导向管、上下导向筒和驱动杆外导管等）都简化为不同特性的梁，图2 程序流程图形成节点、单元、材料和载荷数据文件开始读取节点、单元和材料参数形成单元质量和刚度矩阵形成整体质量、刚度和阻尼矩阵t=t0+Dt t0=0读取边界条件解横向振动方程判断运动部件与导向部件是否接触判断是否穿透减少Dt重新计算确定碰撞力摩擦力计算下落运动判断运动部件是否到位输出时间结束流体阻力计算否前一时刻下落位置和速度否是否是是而建立起二维有限元分析模型，见图3。计算结果和试验结果比较见表1，计算结果图示见图48。从图表中可以看到：控制棒在下落的过程中受到的主要阻力来自于流体，但随着外载荷的增大，摩擦力对落棒时间的延迟作用增大。表1 CRFD计算的落棒时间和试验结果的比较 （冷却剂纵向流速为0） 单位：st5t5+t6计算值试验值计算值试验值静态落棒1.251.241.831.80SSE地震落棒1.341.341.901.89SSE地震落棒1.421.412.011.98图3 控制棒落棒计算模型示意图导向组件控制棒组件S1S2S3S4地震载荷输入zyo附加结构1.6SSE地震落棒1.541.502.102.07 图4 静态落棒位移时程数值计算和试验测试比较 图5 四种输入工况下落棒位移时程计算结果比较 图6 四种输入工况下落棒速度时程计算结果比较图7 SSE地震下落棒控制棒受摩擦力时程计算结果 图8 SSE地震下落棒导向管内压力变化时程计算结果6结论 通过本文的研究开发，只需确定几个基础参数（摩擦系数和碰撞恢复系数及接触时间由试验确定）就能完成落棒时间与历程分析，和已有分析软件相比CRFD程序有以下特点：a. 采用通过能量守恒、动量守恒和冲量定理确定的反力构造法来计算碰撞力，使碰撞力学模型能较好地体现动态碰撞特性；b. 采用直接约束法处理接触问题，很好地解决了动态接触判断问题；c. 采用理论和经验公式相结合的方法计算流体阻力，使流体阻力的计算初步脱离了试验测定；参考文献：1窦一康 姚伟达 杨仁安 姜南燕. 事故工况下控制棒落棒时间分析方法及其工程应用. 第十届全国反应堆结构力学会议议论文集，第十届全国反应堆结构力学会议，上海，1998.11， 3713762 J. Stabel H. P. Fuchs. Influence of seismic excitation on the drop time of control assemblies in PWRs. Transactions of the 9th International Conference on Structural Mechanics in Reactor Technology, Lausanne, 1987, P112111253 R.W.Clough， J.Penzien. 结构动力学. 王光远等译. 科学出版社. 19814 Rusan. Krajcinoric. Vibration of two coaxial cylindrical shell containing fluid. Nuclear Engineering and Design, 1974, 30：2422485金栋平 胡海岩 .碰撞振动及其典型现象. 力学进展，1999.5，第29卷第2期：1551646孙焕纯 宋亚新 刘魏 .结构碰撞的动力响应分析. 计算结构力学及其应用，1994.2，第11卷第1期：31427于建华 魏永涛等 .控制棒组件在流体环境中下落时所受阻力的计算. 核动力工程，2001，第22卷第3期：2362418胡海昌 .弹性力学的变分原理及其应用. 科学出版社，1981年9汪明德 .面向对象程序设计. 电子工业出版社，1993年Software Developing of Control Rod Assembly Drop-time and Course Analysis The Drop-time of control assembly was one of the most important parameters to ensure Nuclear Power Plant running safety. Upon researched works, this paper made a deep study on the issue. Main research methods in this paper were to adapt directly restraint method to deal with contact problem, impact reacting force function to modify impact force and numerical analysis to calculate fluid resistance. At last, the analysis result of control assembly drop-time test of Qingsan 600MWe PWR was presented.Keywords: Control Rod Assembly Drop-time, Contact and Impact, FEM作者简介：孙 磊（1968），男，副研究员，硕士，1996年在中国核动力研究设计院获硕士学位，现为反应堆工程与安全专业博士研究生，主要从事反应堆结构力学的计算和试验研究工作。于建华（1944），男，教授，力学所所长，1986年在法国贝桑松大学获博士学位，主要从事工程结构动力学和非线性计算的研究。魏泳涛（1971），男，讲师，2000年在四川大学获博士学位，主要从事计算力学的研究。