基于热膨胀标定高瞬时性固态相变JMAK方程参数的方法

基于热膨胀标定高瞬时性固态相变jmak方程参数的方法
技术领域
1.本发明涉及一种基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，属于冶金学领域。


背景技术：

2.近年来，随着航空航天技术的发展，针对具有高瞬时性加热、冷却过程的金属相变过程的研究需求在不断增加，但由于高瞬时性加热、冷却过程中的相变过程与传统认知的等温相变过程不同，不符合经典形核理论，也无法通过常用的ttt曲线或cct曲线进行描述，故对这种过程的研究比较困难。直到目前为止，还没有任何一种可以完全准确的描述金属整体在受到高瞬时性加热、冷却作用后的金属相变体积分数的数学模型，因此针对现有的经典相变理论进行改进是一种比较好的解决办法。由于人类科技的进步，在外太空中以及在特种环境下作业的需求越来越大，因此针对金属的高瞬时性整体固态相变的研究是十分必要的。
3.目前在大多数研究与工程实践中，针对固态相变的体积分数的计算与预测都是通过经典jmak方程进行的，但这种理论仅在慢速加热条件下比较有效，在高瞬时性加热、冷却过程中这种方法的计算值与实际实验值之间存在着较大的差距。
4.另外，由于高瞬时性加热、冷却过程的实验条件比较苛刻，目前仅能通过热膨胀实验或热台实验得到相关的参数，而其中由于热膨胀仪的精确度较高，故在实验研究中多是通过热膨胀实验来进行分析。但在热膨胀实验中如何通过试样膨胀的长度变化数据准确推算出试样中奥氏体体积分数依然是一个难点。
5.针对于第一个问题，以volkerschulze为代表的一些学者提出了基于经典jmak方程的微分形式对具有高瞬时性相变特点的过程进行计算，但最终其拟合参数并不能达到理想效果。针对第二个问题，以liu为首的研究者提出了通过类似于铁-碳合金相图中的杠杆法的方法来推算试样膨胀过程中的奥氏体体积分数，但如何准确确定相变发生的长度变化拐点以及如何准确定位相变过程中杠杆的位置依然是一个难题。


技术实现要素：

6.针对现有冶金方法无法准确描述高瞬时加热、冷却过程中金属相变的体积分数的问题，本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，基于差分形式的jmak方程和杠杆法，通过对热膨胀实验数据进行向前、向后差分，准确的判断出热膨胀实验中试样长度变化的拐点，通过建立杠杆求得奥氏体体积分数，并通过获得的奥氏体体积分数以及相应的温度数据拟合、测算出jmak方程中的参数n、c、δh，即基于热膨胀实验拟合、测算、标定出高瞬时性固态相变jmak方程参数，提升高瞬时性固态相变jmak方程对金属在高瞬时性加热过程中的奥氏体相变的预测精度。本发明具有测算、标定、预测精度高、效率高的优点。
7.本发明的目的是通过下述技术方案实现的。
8.本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，根据热膨胀实验数据预设具有自适应性的滤波器窗口阈值；对热膨胀实验数据进行滤波，使数据平滑；为了避免对膨胀过程中拐点的误识别，筛除掉实验准备阶段的热膨胀实验数据；通过对热膨胀实验数据进行向前、向后差分，准确判断出用于滤除热膨胀实验中收缩阶段的拐点；根据判断出的收缩阶段的拐点滤出完整的膨胀阶段的数据；根据差分公式求膨胀过程中各个实验点的长度变化量对时间的偏差分；根据预设的具有自适应性的滤波器窗口阈值对各个实验点的长度变化量对时间的偏差分值再次进行滤波；根据偏差分平滑值，判定膨胀过程中奥氏体化过程的拐点位置；建立用于求得奥氏体体积分数的杠杆；求得奥氏体体积分数；根据奥氏体体积分数以及相应的温度数据，拟合、标定出jmak方程中的参数n、c、δh并求奥氏体体积分数计算值；求奥氏体体积分数计算值与实际实验数据的奥氏体体积分数的方差；重复迭代过程，测算出高瞬时性固态相变jmak方程的最佳参数。本发明具有标定精度高、效率高的优点，适合于任何可通过jmak方程预测其相变过程的金属材料。
9.本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，包括如下步骤：
10.步骤1、根据热膨胀实验数据预设具有自适应性的滤波器窗口阈值。
11.步骤1的具体实现方法为：
12.步骤1.1、计算热膨胀实验数据中各个数据点的温度变化速率，公式为式中i为当前计算的实验数据的编号，t为当前实验数据的温度，t为当前实验数据的时间，v
t
为当前实验数据的温度变化速率。
13.步骤1.2、基于步骤1.1求得的各个热膨胀数据点的温度变化速率，求出最大的温度变化速率v
tmax
。
14.步骤1.3、根据公式δt＝t
i+1-ti计算各个数据点的时间间隔δt，式中i为当前计算的实验数据的编号，t为当前实验数据的时间。
15.步骤1.4、基于步骤1.3得到的时间间隔δt，求得出现次数最多的时间间隔δt
multiple
。
16.步骤1.5、基于步骤1.2和步骤1.4的计算结果，求得用于定义滤波器窗口大小阈值的自适应算子operator，公式为operator＝v
tmax
×
δt
multiple
。
17.步骤1.6、基于步骤1.5求得的自适应算子operator，根据指数函数求得滤波器窗口阈值threshold。
18.作为优选，步骤1.6根据指数函数求得滤波器窗口阈值threshold，选用如下公式实现：
19.threshold＝0.2732
×
operator-2.667
20.步骤1.7、基于步骤1.6求得的滤波器窗口阈值threshold，对滤波器窗口阈值threshold进行向上取整。
21.步骤1.8、对步骤1.7向上取整后的滤波器窗口阈值threshold，进行分段函数修正。
22.作为优选，对步骤1.7向上取整后的滤波器窗口阈值threshold进行分段函数修正，选用如下公式实现：
[0023][0024]
步骤2、基于步骤1预设的具有自适应性的滤波器窗口阈值，对热膨胀实验数据进行滤波，使数据平滑。
[0025]
步骤3、为了避免对膨胀过程中拐点的误识别，筛除掉实验准备阶段的热膨胀实验数据。
[0026]
步骤3的具体实现方法为：
[0027]
步骤3.1、取实验数据的长度变化量并与预设初始值对比，筛除掉实验准备阶段的热膨胀实验数据。
[0028]
作为优选，取实验数据的长度变化量并与0.05初始值对比。
[0029]
步骤3.2、重复步骤3.1直至发现数据大于预设初始值，记录当前数据的位置。
[0030]
作为优选，取试验数据的长度变化量大于0.05，记录当前数据的位置。
[0031]
步骤3.3、舍弃掉记录的位置之前的实验数据，即筛除掉实验准备阶段的热膨胀实验数据以避免对膨胀过程中拐点的误识别。
[0032]
步骤4、基于步骤3整理、筛选过的实验数据，通过对热膨胀实验数据进行向前、向后差分，准确判断出用于滤除热膨胀实验中收缩阶段的拐点。
[0033]
步骤4的具体实现方法为：
[0034]
步骤4.1、取热膨胀实验数据点，并计算该点的长度变化对时间的向前偏差分以及向后偏差分，具体计算公式为：
[0035][0036][0037]
式中为长度变化对时间的向前偏差分，为长度变化对时间的向后差分，δl为实验中试件长度变化量，t为对应的实验时间，i为当前计算的实验数据的编号。
[0038]
步骤4.2、基于步骤4.1的计算结果，计算的值，若乘积小于0且实验数据的时间在整体实验时间的前一部分内，记录该点为拐点。即定义tj为目前计算数据的实验时间，δlj为该点的长度变化量，当该点满足
[0039][0040]
时，则记录该点为一个拐点。上述公式中，为长度变化对时间的向前偏差分；
为长度变化对时间的向后差分；t
max
为热膨胀试验数据的最大持续时间。
[0041]
步骤4.3、重复步骤4.1、步骤4.2，直至所有实验数据判断完成，准确判断出用于滤除热膨胀实验中收缩阶段的拐点。
[0042]
步骤5、根据步骤4准确的判断出热膨胀实验中试样长度变化的拐点，确定完整的试样膨胀过程的数据。
[0043]
步骤5的具体实现方法为：
[0044]
步骤5.1、统计步骤4.3中用于滤除热膨胀试验中收缩阶段的拐点中实验时间最靠后的拐点。
[0045]
步骤5.2、以步骤5.1中试验时间最靠后的拐点为分界线，舍去该拐点后的所有实验数据，仅保留该拐点前的数据，此时保留的数据即为完整的试样膨胀过程的实验数据。
[0046]
步骤6、根据差分公式求膨胀过程中各个实验点的长度变化量对时间的偏差分，得到各个实验点的长度变化量对时间的偏差分值。
[0047]
作为优选，步骤6所述差分公式为其中为长度变化量对时间的偏差分；δl为各个实验点的长度变化量；ti为实验点的时间。
[0048]
步骤7、根据预设的具有自适应性的滤波器窗口阈值，对步骤6得到各个实验点的长度变化量对时间的偏差分值再次进行滤波，使各个实验点的长度变化量对时间的偏差分值平滑。
[0049]
步骤8、根据步骤7中各个试验点的长度变化量对时间的偏差分平滑值，判定膨胀过程中奥氏体化过程的拐点位置。
[0050]
步骤8的具体实现方法为：
[0051]
步骤8.1、根据公式进行计算，若求得值为负数，则视为一个可能的拐点。
[0052]
步骤8.2、针对这个可能的拐点的位置进行判断，是否处于实验数据的最前段或最后段，即是否为实验数据的前4个或后4个，如果不是则可继续视为一个可能的拐点。
[0053]
步骤8.3、以三个数值构成一个数据窗口，若所述窗口内的三个数值均为负数，则判定该可能的拐点位于振荡中，不再视为拐点，否则视为真实的拐点。
[0054]
步骤8.4、重复步骤8.1、步骤8.2、步骤8.3直至所有实验数据判定完毕，直至准确判定出实验中膨胀过程中奥氏体化过程的拐点位置。
[0055]
步骤9、根据步骤8判定的膨胀过程中奥氏体化过程的拐点位置，建立用于步骤10求得奥氏体体积分数的杠杆。
[0056]
步骤9的具体实现方法为：
[0057]
步骤9.1、根据步骤8.4准确判定出的实验中膨胀过程中奥氏体化过程的拐点位置中的第一个拐点以及最后一个拐点将实验数据分成三段。设第一个拐点的位置为g1，最后一个拐点的位置为gn，倒数第二个拐点为g
n-j
，实验中膨胀过程的第一个数据点为g
min
，最后
一个数据点为g
max
，则其分段表达式为：
[0058][0059]
步骤9.2、以步骤9.1将实验数据分成三段后的第一段和第三段数据为基础，以长度变化为因变量，实验温度为自变量进行最小二乘拟合，得到两个函数，分别记为
[0060][0061]
步骤9.3、基于步骤9.2的最小二乘回归结果，求拟合优度r2。若回归结果的r2值过小，则舍弃gn，使用g
n-j
作为拐点，重复步骤9.1、步骤9.2的计算过程，直至求得比较好的回归结果。
[0062]
作为优选，应将r2＞0.75作为筛选回归结果的判定标准。
[0063]
步骤9.4、基于步骤9.3得到的比较好的回归结果，将第二段数据的实验温度作为自变量带入至y1、y2中，求得虚拟长度变化量δl
’1、δl
’2。
[0064]
步骤9.5、以第二段数据的实际长度变化δl为杠杆中心支点，步骤9.4中求得的虚拟长度变化量δl
’1、δl
’2分别作为杠杆上支点与下支点，构建用于步骤10求得奥氏体体积分数的杠杆。
[0065]
步骤10、基于步骤9建立的奥氏体体积分数的杠杆，求得奥氏体体积分数，公式为式中f
α
为求得的奥氏体体积分数，为百分数形式；δl
’1、δl
’2分别为杠杆的上支点与下支点；δl为杠杆的中心支点。
[0066]
步骤11、重复步骤9、步骤10直至所有奥氏体体积分数数据计算完成，获得所有奥氏体体积分数以及相应的温度数据。
[0067]
步骤12、根据步骤11获得的所有奥氏体体积分数以及相应的温度数据，标定出jmak方程中的参数n、c、δh，即基于热膨胀实验标定出高瞬时性固态相变jmak方程参数，并求得在标定参数下的奥氏体体积分数计算值。
[0068]
步骤12的具体实现方法为：
[0069]
步骤12.1、对传统jmak方程进行差分变化，所得方程为：
[0070][0071]
式中n,c,δh为带求解的参数，k为玻尔兹曼常数，ti为当前计算点的温度，f
α(i)
为当前计算点的奥氏体体积分数。
[0072]
步骤12.2、对差分后的方程两边同时取对数并进行整理，可得：
[0073][0074]
步骤12.3、根据步骤11获得的所有奥氏体体积分数以及相应的温度数据，随机任
取3组数据带入至步骤12.2中的等式中，构建矩阵形式的方程组。即将数据带入后，将等式变换为的形式。
[0075]
步骤12.4、求解广义逆矩阵
[0076]
步骤12.5、通过求解求得
[0077]
步骤12.6、以步骤12.5求得的为基础，通过将实验数据带入至原始jmak方程组
[0078][0079]
求得在标定参数下的奥氏体体积分数计算值。
[0080]
为了进一步提高jmak方程的预测精度，还包括步骤13、步骤14：
[0081]
步骤13、基于步骤12标定出的jmak方程中的参数n、c、δh以及奥氏体体积分数计算值和实际实验数据的奥氏体体积分数求得当前拟合参数下的方差。
[0082]
步骤14、重复步骤12、步骤13迭代计算过程直至测算出具有最小方差的一组参数，即测算出高瞬时性固态相变jmak方程的最佳参数。
[0083]
作为优选，步骤14可以通过如下方式实现：
[0084]
当新计算出的方差值小于当前方差值时，则舍弃当前参数组，记录方差更小的参数组，一直循环此过程，直至当前参数组的迭代次数超过总允许迭代次数的一半，输出所记录的参数组，完成优化过程。此时即为jmak方程针对于当前实验数据的最佳标定参数，也是测算出的高瞬时性固态相变jmak方程最佳参数。
[0085]
还包括步骤15：基于步骤14测算出的高瞬时性固态相变jmak方程最佳参数，能够预测金属材料在受到高瞬时性加热时，奥氏体体积分数的变化情况，进而预测材料在高瞬时性加热过程中的弹性模量、泊松比、抗拉强度，为疲劳寿命计算提供支撑数据，从而达到延长材料使用寿命的目的。另外，得到的最佳参数还能为受热过程中工件危险点的测算提供支持。
[0086]
有益效果：
[0087]
1、本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，根据热膨胀实验数据预设具有自适应性的滤波器窗口阈值，对热膨胀实验数据进行滤波，使数据平滑；通过对热膨胀实验数据进行向前、向后差分，准确判断出用于滤除热膨胀实验中收缩阶段的拐点；根据判断出的收缩阶段的拐点滤出完整的膨胀阶段的数据；根据预设的具有自适应性的滤波器窗口阈值对滤出的完整的膨胀阶段的数据各个实验点的长度变化量对时间的偏差分值再次进行滤波；根据滤波后的偏差分值准确判定膨胀过程中的长度变化拐点；根据准确判定的长度变化拐点建立用于求得奥氏体体积分数的杠杆，并求得奥氏体体积分数；根据求得的奥氏体体积分数以及相应的温度数据，拟合、标定出jmak方程中的参数n、c、δh。为了进一步提高jmak方程的预测精度，采用求奥氏体体积分数计算值与实际实验数据的方差并以方差最小为目标进行迭代的方法，测算出高瞬时性固态相变jmak方程的最佳参数，具有标定精度高、效率高、适合于任何可通过jmak方程预测其相变过程的金属材料的优点。
[0088]
2、本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，采用自适应滤波器窗口阈值的方式对实际实验数据进行平滑处理，并以偏差分乘积作为判断实际拐点的依据，能够准确的识别收缩阶段的拐点，提高测算精度和数据处理效率。
[0089]
3、本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，采用随机取实验数据带入至方程中的计算的方式，相较于传统的优化方式具有更高的计算效率。
[0090]
4、本发明公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，采用广义逆矩阵的方式求解拟合参数，避免矩阵奇异的问题，使得求解结果更加精确，计算效率更高。
附图说明
[0091]
图1为本发明的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法流程图，其中图1(a)为数据处理部分，图1(b)为参数测算部分；
[0092]
图2为45crnimova钢在热膨胀实验中的实际实验数据曲线图；
[0093]
图3为对45crnimova钢的实验数据进行第一次拐点计算的计算结果曲线图；
[0094]
图4为45crnimova钢在实验过程中的纯膨胀阶段曲线图；
[0095]
图5为45crnimova钢实验数据对时间偏差分后滤波的曲线图；
[0096]
图6为对45crnimova钢的纯膨胀实验数据进行第二次拐点计算的计算结果曲线图；
[0097]
图7为杠杆建立示意图；
[0098]
图8为实验计算结果与数值拟合计算结果的对比曲线图。
具体实施方式
[0099]
为了更好地说明本发明的目的和优点，下面结合附图和实例对发明内容做进一步说明。
[0100]
实施例1：
[0101]
采用本发明方法对45crnimova钢材的热膨胀实验数据进行处理与拟合。
[0102]
如图1(a)、(b)所示，本实施例公开的基于热膨胀标定高瞬时性固态相变jmak方程参数的方法，具体实现步骤如下：
[0103]
步骤1、根据热膨胀实验数据预设具有自适应性的滤波器窗口阈值。计算结果如表1所示。
[0104]
表1自适应算子确定过程
[0105]
最大温度变化速率6.5℃/s最高频次出现时间间隔1s自适应算子7确定的滤波器窗口阈值8
[0106]
步骤2、基于步骤1预设的具有自适应性的滤波器窗口阈值，对热膨胀实验数据进行滤波，使数据平滑。
[0107]
步骤3、为了避免对膨胀过程中拐点的误识别，筛除掉实验准备阶段的热膨胀实验数据。最终处理后的数据曲线图如图2所示。
[0108]
步骤4、基于步骤3整理、筛选过的实验数据，通过对热膨胀实验数据进行向前、向后差分，准确判断出用于滤除热膨胀实验中收缩阶段的拐点。最终通过计算得到的拐点位置如图3所示。
[0109]
步骤5、根据步骤4准确的判断出热膨胀实验中试样长度变化的拐点，确定完整的试样膨胀过程的数据。筛选后的试样实验过程的完整膨胀过程如图4所示。
[0110]
步骤6、根据差分公式求膨胀过程中各个实验点的长度变化量对时间的偏差分，得到各个实验点的长度变化量对时间的偏差分值。
[0111]
步骤7、根据预设的具有自适应性的滤波器窗口阈值，对步骤6得到各个实验点的长度变化量对时间的偏差分值再次进行滤波，使各个实验点的长度变化量对时间的偏差分值平滑。滤波平滑后的偏差分数值如图5所示。
[0112]
步骤8、根据步骤7中各个试验点的长度变化量对时间的偏差分平滑值，判定膨胀过程中奥氏体化过程的拐点位置。最终确定的拐点位置如图6所示。
[0113]
步骤9、根据步骤8判定的膨胀过程中奥氏体化过程的拐点位置，建立用于步骤10求得奥氏体体积分数的杠杆。杠杆的建立示意图如图7所示。
[0114]
步骤10、基于步骤9建立的奥氏体体积分数的杠杆，求得奥氏体体积分数。根据图7中所示，奥氏体体积分数的计算公式为eg/ef
[0115]
步骤11、重复步骤9、步骤10直至所有奥氏体体积分数数据计算完成，获得所有奥氏体体积分数以及相应的温度数据。最终计算完成的奥氏体体积分数以及其对用温度如表2所示。
[0116]
表2奥氏体体积分数实际实验值
[0117][0118][0119]
步骤12、根据步骤11获得的所有奥氏体体积分数以及相应的温度数据，标定出jmak方程中的参数n、c、δh，即基于热膨胀实验标定出高瞬时性固态相变jmak方程参数，并求得在标定参数下的奥氏体体积分数计算值。最终得到的参数为n＝1.7483，c＝0.7483，δh＝2.957
×
10-20
，在该组参数下的奥氏体体积分数计算值如表3所示。
[0120]
表3奥氏体体积分数计算值
[0121][0122][0123]
步骤13、基于步骤12标定出的jmak方程中的参数n、c、δh以及奥氏体体积分数计算值和实际实验数据的奥氏体体积分数求得当前拟合参数下的方差。根据计算结果，当前参数情况下的方差为0.1226。
[0124]
步骤14、重复步骤12、步骤13迭代计算过程直至测算出具有最小方差的一组参数，即测算出高瞬时性固态相变jmak方程的最佳参数。根据计算结果，n＝1.7483，c＝0.7483，
δh＝2.957
×
10-20
即为具有最小方差的参数组，其方差为0.1226即为最小值。该组参数组的计算结果与实际实验值的对比如图8所示。
[0125]
步骤15、基于步骤14测算出的高瞬时性固态相变jmak方程最佳参数，可以预测金属材料在受到高瞬时性加热时，奥氏体体积分数的变化情况，进而预测材料在高瞬时性加热过程中的抗拉强度，为疲劳寿命计算提供支撑数据。根据步骤14的计算结果，当升温速率6.5℃/s、实际温度为767.8℃时，45crnimova钢材的奥氏体含量为89.96％，当钢材在该温度下完成淬火后，其奥氏体全部转变为马氏体，即完成淬火后，该钢材的马氏体含量为89.96％，根据束德林所发表的结果，金属抗拉强度与马氏体含量之间存在如下关系：
[0126]
σb(mpa)＝365+16.3
×
m(％)
[0127]
式中σb为抗拉强度，单位为mpa，m为马氏体百分含量。
[0128]
根据该公式求得45crnimova钢材在经过高瞬时性加热过程后的抗拉强度为1831.348mpa，同时普通的未经过热处理的45crnimova钢材的抗拉强度为1553mpa，钢材的性能有所提升。
[0129]
以上所述的具体描述，对发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。