国产SA508-3钢拉伸力学性能的近场动力学模拟

刘贞谷， 郭居上， 李丽娟， 庾明达， 陈建国， 石凯凯， 江小州

(1. 中国核动力研究设计院核反应堆系统设计技术国家重点实验室， 四川成都 610213；2. 哈尔滨工业大学航天学院航天科学与力学系， 黑龙江哈尔滨 150001)

随着国产SA508-3钢研制的成功，中国学者们对其进行大量的研究，其力学性能成为关注的焦点，我国的压水堆核电站—回路压力容器也广泛采用这种材料[1]。利用试验以及有限元方法等手段，研究者们已经取得了大量的研究成果。余美芳等人收集了大量的数据，对国产SA508-3 钢的力学性能有了较为全面的了解，通过和西方RPV 钢的对比，在拉伸、断裂韧性等方面展开了性能评估。基于耦合损伤的本构模型，张丽屏等[2]模拟了SA508-3钢的循环塑性变形过程，在评价以此材料制造的核电设备的力学性能方面提供了相关经验。林赟等[3]研究了国产SA508-3钢辐照性能，在研究堆上进行了加速辐照试验，进而评估其拉伸和冲击等力学性能。在国产SA508-3钢的断裂韧性方面，彭啸等[4]使用了主曲线方法，取得了1/2PCVN试样以及PCVN试样的研究成果。

有限元、边界元、无网格法、分子动力学的发展，使得针对材料力学性能的表征和模拟也引入了各种新方法，在某些特定方面能取得较好的模拟效果。如针对材料及结构的断裂方面，传统的有限元面临极大的困难，Belytschko和Black[5]以及Moes[6]提出并发展了扩展有限元(XFEM)的思想和方法，而分子动力学模拟(MDS)也解决了这些困难，Cox等人[7]针对动态断裂的物理过程给出了更为直观的理解。研究了原子之间的相互作用，Eringen和Edelen[8]、Kroner[9]以及Kunin[10]进一步提出了非局部连续理论，把具有特征长度参数的长程影响考虑在内。为了从理论底层解决不连续的难题，同时考虑了以上理论的特点和不足，Silling于2000年[11]及2007年[12]提出新的非局部框架下的固体力学理论——近场动力学理论(Peridynamic Theory，简称PD)，黄丹等国内学者[13]在国内首先对其进行综述，译为“近场动力学理论”[14]。本文利用近场动力学理论，编写Fortran程序，模拟核电设备中SA508-3钢材料受拉的力学性能，并对模拟结果进行评价。

1 近场动力学基本理论

近场动力学把连续介质力学中的微分方程用空间积分方程来代替，由于积分在不连续处也是可积的，因此自然在断裂破坏等问题上也好处理。

如图1所示，近场动力学的基本思想为，考虑某个占据体积V(k)的物质点x(k)，与有限范围内物质点x(j)(j=1,2,...,)存在着相互作用，作用范围δ称作近场范围，反映了近场动力学的非局部特性，即，范围内存在作用，超出该范围意味着作用消失。

图1 物质点相互作用范围示意

在作用域内，PD运动方程可以写为：

∀x∈R,t≥0

(1)

式中：H代表在近场范围内进行积分，若使用以下记法：

ξ=x′-x,η=u(x′,t)-u(x,t)

(2)

(3)

则如图2所示，ξ、η分别代表两个作用点之间的相对位置以及相对位移，而η+ξ则表示结构在变形后两个点之间的相对位置。f为近场动力学中的对点力函数，可以形象地把两个点之间的作用看成“键”，满足牛顿第三定律

图2 PD理论中物质点间作用

(4)

以及角动量守恒定律

(ξ+η)×f(η,ξ)=0 ∀η,ξ

(5)

对微势能函数进行求导可以得到f，即：

(6)

微势能函数w的量纲为能量/体积2，反映了单个键储存的能量。对于任意一个物质点，可以把其应变能密度写作：

(7)

为了模拟材料及结构发生变形后的破坏，PD引入了伸长率s来反映键的状态，即：

(8)

当s在变形过程中不断增大直至达到临界伸长率s0时，键就判定为永久断裂，代表该作用力消失。例如，标准的微弹性脆性(PMB)的对点力函数可以写为：

f(η,ξ)=g(s(t,ξ))μ(t,ξ)

(9)

式中：g是一个线性的标量函数，μ的取值为0或1(图3)，即：

g(s)=cs∀s

(10)

(11)

图3 PMB材料的PD本构关系

则物质点x处的损伤可以由键数量的统计值给出，大小为已断开的键占所有键的数量比值，即：

(12)

因此，可以看出，对点力函数天然适用于不连续，在整个过程中，裂纹自发萌生并扩展，不需要其它的断裂判据或者扩展条件等。

2 SA508-3钢拉伸力学性能的近场动力学模拟

同时利用PD理论和经典理论计算物质点的应变能密度，通过对比可以建立两种理论基本物理量之间的关系，则推导出微模量函数为：

(13)

将SA508-3钢考虑为一块薄板，材料属性取自ASME规范[15]，其基本参数见表1～表3。板的两端承受均匀拉伸载荷，大小为p0=170MPa，编写程序模拟其弹性变形过程(图4)。

表1 SA-508 Grade 3 Class 2的材料特性

表2 薄板的几何及材料参数

表3 PD模型离散参数

图4 SA508-3钢薄板受拉伸载荷作用

时间步长取nt=2000，计算得到了薄板离散的所有物质点在x、y向的位移值大小(图5)。

(a) X方向位移

(b) Y方向位移图5 X、Y向位移云图

针对该算例，存在弹性力学上的解析解，以水平和竖直两个方向中线上各物质点为参考，将PD计算的位移解与其解析解作比较(图6)。

从图6中可以看出，程序的PD计算值与解析解基本一致，进行误差分析后，结果如表4所示。

(a) 位移X向解析解与PD解对比

(b) 位移Y向解析解与PD解对比图6 位移解析解与PD解对比

表4 中线误差分析

程序收敛性方面，考察任意一物质点，提取出其每一计算步下的位置解，得到位移解随时间的变化曲线(图7)，可以看出，在时间步到800后，计算结果趋于稳定。

图7 位移结果收敛性

3 结 论

从计算结果可以看出，基于近场动力学理论编写的计算程序在SA508-3钢受拉弹性阶段的计算结果与解析解误差很小，为后续分析较大的反应堆设备结构积累经验，同时也可以进一步模拟损伤断裂等不连续力学特征。