本发明涉及摩擦焊工艺模拟技术，具体是摩擦焊工艺热力流微观组织多物理场数值计算方法。
背景技术：
摩擦焊是利用搅拌头与工件表面相互运动和摩擦所产生的热量，融化工件表面，使其达到热塑性状态，然后搅拌头迅速顶锻，使两块板焊接在一起的工艺方式。
摩擦焊工艺的计算机模拟主要集中在四个方面：温度场模拟，研究焊接过程工件的温度和传热过程；流程模拟，研究焊接过程中工件焊缝附近塑性区域的流动状况；应力应变模拟，研究焊接过程中工件整体的应力和应变分布或研究搅拌头的应力与寿命分析；组织性能模拟，研究焊接过程中工件微观组织转变过程。
目前摩擦焊工艺的计算机模拟只能独立的针对摩擦焊工艺的某一或两个方面进行数值模拟计算，借用的是通用有限元的计算算法，在计算过程中需做很多简化和假设，并不能完全真实的反映实际的摩擦焊工艺和摩擦焊接过程的真实物理现象，也不能建立完整的摩擦焊物理与数学耦合模型，相应的计算效率和精度也较低，且存在因二次开发导致的操作难度较大、运行不稳定等问题。
技术实现要素：
本发明的目的在于解决现有摩擦焊工艺模拟存在的效率和精度低、与实际的摩擦焊工艺吻合度低的问题，提供了一种摩擦焊工艺热力流微观组织多物理场数值计算方法，其用于摩擦焊工艺模拟时可以高效、准确地得到模拟计算结果，且能提升与实际摩擦焊工艺的吻合度。
本发明解决上述问题主要通过以下技术方案实现：摩擦焊工艺热力流微观组织多物理场数值计算方法，包括以下步骤：
步骤一、建立热-流耦合数值模型，利用考虑温度的欧拉算法进行温度和流场完全耦合计算，得到工件焊接过程任意时刻的温度和流场大小及分布；
步骤二、利用Leblond相变模型与算法进行微观组织转变计算，得到相变组织分布；
步骤三、基于得到的温度、流场及相变组织分布，建立热-力-组织耦合数值模型，计算焊接过程工件塑性区域和应力应变。其中，Leblond相变模型是一种对扩散型相变和马氏体相变都适用的连续冷却和连续加热相变动力学模型，可以应用于两相或多相相变的计算。其可以把HAZ所经历的热循环和发生的相变过程相联系，因此可以用于HAZ相变的计算。
进一步的，所述步骤一中获取温度时包括工件内部热传导计算和工件与外界的热交换计算，其中，热传导是根据傅里叶传热定律来计算，傅里叶传热定律的公式为：
其中，k(x,y,z)为材料的热传导系数，A是传热面积，T为温度，x为在导热面上的坐标；
工件与外界的热交换利用欧拉方程来进行计算，欧拉方程的公式为：
；
其中ρ是材料密度，C是材料比热，Q是材料内部总热量，V和为算子。
进一步的，所述步骤一中获取流场大小及分布时先通过雷洛系数来判断是否属于紊流，然后基于N-S方程和欧拉方程进行计算，其中N-S方程为：
其中，△是拉普拉斯算子，μ是动力粘性系数，F是载荷力，p为压力，v是速度分量，t是时间。其中，N-S方程为Navier-Stocks方程。
进一步的，所述步骤一～步骤三计算过程中根据实际模型采用拉格朗日-欧拉网格准备摩擦焊数值模拟网格，其中，焊缝区域采用拉格朗日网格，非焊接区域采用欧拉网格；所述步骤一中针对每个网格节点先计算温度，然后以温度作为该节点流体计算的边界条件之一，计算使得每个网格节点上都保存有温度和流场结果。传统焊接模拟中，一般采用欧拉网格或拉格朗日网格。欧拉网格可以方便的定义计算网格的初始状态，描述流体质点的变化，但很难跟踪流体质点，即精确定义液体形状和位置；拉格拉日网格可以很容易跟踪流体质点运动，描述液体形状，但若流体质点运动剧烈时，会导致计算网格的扭曲，网格奇异问题，造成不收敛。本发明计算过程中，首次采用拉格朗日-欧拉耦合网格，即，充分利用欧拉和拉格拉日网格的优点，网格节点运动方式自由灵活，网格节点可以跟随流体质点一起运动，也可以固定不变，甚至可以在一个方向固定而在其他方向上随流体一起运动变化，因此该网格特别适合大变形大位移非线性过程的摩擦焊模拟过程。本发明的每个网格节点上都保存有温度和流场结果，这样依次进行计算迭代即可实现温度流程完全耦合。
进一步的，所述步骤二中微观组织转变计算根据Leblond Model和K-M模型来计算，Leblond Model方程为：
其中，T为温度，为平衡时刻的相比例，为相变过程单一相比例，为相变在温度T时的延迟时间，为温度T时相的个数，Leblond Model方程中各种参数通过材料的CCT曲线获取；
K-M方程为：
其中，T为温度，为温度T时马氏体的比例，为温度T时奥氏体的比例，k为试验参数，Ms为马氏体开始转变温度。
进一步的，所述步骤三中应力应变计算过程中焊缝流体粘性力的计算公式为：
其中，是动力学粘性系数或摩擦系数，u是速度在z方向的分量。
进一步的，所述热-流完全耦合数值模型和热-力-组织耦合数值模型中搅拌头与工件间应力应变计算模型为：
其中，为粘性耗散，为热流密度，为温度，是材料密度，C是材料比热，是速度分量。
进一步的，所述摩擦焊过程中应变的计算公式如下：
其中，为弹性应变，为塑性应变，是蠕变应变，为相应变，为由相转变的体积应变，为相变过程的塑性应变；
相应变的方程为：
式中，T是温度，是初始温度，是膨胀系数，是指定温度下的单一相应变，是相的体积分数；
相变过程的塑性应变方程为：
式中，是应力张量，是由于相转变引起的体积的改变，s是屈服应力，是已经转变的相的百分数，是相转变增量。
进一步的，所述步骤一至步骤三的计算过程中材料的屈服准则包括特斯卡准则和米塞斯准则，特斯卡准则的公式为：
，其中，是纯剪切时的屈服剪切力，和分别是主应力；
米塞斯准则的公式为：，其中，是纯剪切时的屈服剪切力，、和分别是x、y、z三个方向主应力。其中，特斯卡准则（Tresca准则）和米塞斯准则（Mises准则）为焊接领域广泛应用准则，本发明同时应用两个准则以实现互补，使得本发明应用时更与工程实际相匹配。
进一步的，所述热-流完全耦合数值模型和热-力-组织耦合数值模型建立时采用Johndon-Cook模型来描述摩擦焊材料，Johndon-Cook模型的表达式为：
其中，A是材料的屈服强度，B是材料的抗拉强度，n是材料的硬化指数，C是应变率敏感指数，m是材料的温度软化指数，是等效塑性应变，是参考应变率，是材料熔点，是参考温度，T是材料当前温度。
与现有技术相比，本发明具有以下有益效果：本发明建立了摩擦焊完整的物理耦合模型，本发明应用时能同时计算摩擦焊温度、流体、微观组织及应力应变，进而可实现摩擦焊工艺的多物理场直接耦合数值模拟计算。综上所述，本发明用于摩擦焊工艺模拟时可以高效、准确地得到模拟计算结果，且能提升与实际摩擦焊工艺的吻合度。
附图说明
此处所说明的附图用来提供对本发明实施例的进一步理解，构成本申请的一部分，并不构成对本发明实施例的限定。在附图中：
图1为本发明一个具体实施例的流程图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚明白，下面结合实施例和附图，对本发明作进一步的详细说明，本发明的示意性实施方式及其说明仅用于解释本发明，并不作为对本发明的限定。
实施例：
如图1所示，摩擦焊工艺热力流微观组织多物理场数值计算方法，包括依次进行的以下步骤：步骤1、建立热-流耦合数值模型；步骤2、进行温度和流场完全耦合计算，得到工件焊接过程任意时刻的温度和流场大小及分布；步骤3、进行微观组织转变计算，得到相变组织分布；步骤4、建立热-力-组织耦合数值模型；步骤5、计算焊接过程工件塑性区域和应力应变。
本实施例的热-流完全耦合数值模型和热-力-组织耦合数值模型建立时采用Johndon-Cook模型来描述摩擦焊材料，Johndon-Cook模型的表达式为：
其中，A是材料的屈服强度，B是材料的抗拉强度，n是材料的硬化指数，C是应变率敏感指数，m是材料的温度软化指数，是等效塑性应变，是参考应变率，是材料熔点，是参考温度，T是材料当前温度。
本实施例摩擦焊过程中，搅拌头与工件摩擦加热功率为：
式中，n为转速，P是搅拌头压力，R是搅拌头半径，是摩擦系数，M是转矩，是角速度。
本实施例步骤2中进行温度和流场完全耦合计算时包括热传导中热计算和流体计算两部分内容，其中，热传导中热计算包括工件内部热传导计算和工件与外界的热交换计算，热传导是根据傅里叶传热定律来计算，傅里叶传热定律的公式为：；其中，k(x,y,z)为材料的热传导系数，A是传热面积，T为温度，x为在导热面上的坐标。
本实施例中工件与外界的热交换利用欧拉方程来进行计算，欧拉方程的公式为：；其中ρ是材料密度，C是材料比热，Q是材料内部总热量，V和为算子。
本实施例中流体的计算先通过雷洛系数来判断是否属于紊流，然后基于N-S方程和欧拉方程进行计算，其中N-S方程为：；其中，△是拉普拉斯算子，μ是动力粘性系数，F是载荷力，p为压力，v是速度分量，t是时间。
本实施例应用时根据实际模型采用拉格朗日-欧拉网格准备摩擦焊数值模拟网格，其中，焊缝区域采用拉格朗日网格，非焊接区域采用欧拉网格，两种网格交接的地方采用过渡网格进行匹配。温度和流场的完全耦合计算是通过选择合适的网格（拉格朗日-欧拉耦合网格），针对每个网格节点先计算温度，然后温度作为该节点流体计算的边界条件之一，利用N-S方程进行计算，每个网格节点上都保存有温度和流场结果，这样依次进行计算迭代实现温度流程完全耦合。本实施例利用欧拉-拉格朗日网格形成的热-流网格在保证温度计算精度的同时考虑流场的计算，组织-力学网格可以是欧拉-拉格朗日网格，也可以是普通网格，主要考虑组织应力耦合计算过程收敛性。
本实施例步骤3中微观组织转变计算属于材料冶金的计算范畴，本实施例步骤3在具体实施时通过耦合步骤2中的热根据Leblond Model(JMA)和K-M模型来计算，Leblond Model方程为：
其中，T为温度，为平衡时刻的相比例，为相变过程单一相比例，为相变在温度T时的延迟时间，为温度T时相的个数，Leblond Model方程中各种参数通过材料的CCT曲线获取；
K-M方程为：
其中，T为温度，为温度T时马氏体的比例，为温度T时奥氏体的比例，k为试验参数，Ms为马氏体开始转变温度。
本实施例中应力应变的计算，需要考虑流场对其影响，在利用N-S方程计算出焊接过程液体的流动情况下，得到焊接塑性区域，应力和应变会引起材料的非线性行为，从而影响材料组织转变。其中，应力应变计算过程中焊缝流体粘性力的计算公式为：；其中，是动力学粘性系数或摩擦系数，u是速度在z方向的分量。
本实施例中热-流完全耦合数值模型和热-力-组织耦合数值模型中搅拌头与工件间应力应变计算模型为：
其中，为粘性耗散，为热流密度，为温度，是材料密度，C是材料比热，是速度分量。
摩擦焊过程中应变的计算公式如下：
其中，为弹性应变，为塑性应变，是蠕变应变，为相应变，为由相转变的体积应变，为相变过程的塑性应变；
相应变的方程为：
式中，T是温度，是初始温度，是膨胀系数，是指定温度下的单一相应变，是相的体积分数。
相变过程的塑性应变方程为：
式中，是应力张量，是由于相转变引起的体积的改变，s是屈服应力，是已经转变的相的百分数，是相转变增量。
本实施例计算过程中材料的屈服准则包括特斯卡准则和米塞斯准则，特斯卡准则的公式为：，其中，是纯剪切时的屈服剪切力，和分别是主应力；米塞斯准则的公式为：，其中，是纯剪切时的屈服剪切力，、和分别是x、y、z三个方向主应力。
本实施例利用考虑温度的欧拉算法进行温度和流场完全耦合计算，得到工件焊接过程任意时刻的温度和流场大小及分布，利用Leblond相变模型与算法进行微观组织转变计算，得到相变组织分布，并基于得到的温度、流场及相变组织分布，计算焊接过程工件塑性区域和应力应变。
本实施例应用时，利用欧拉-拉格朗日网格可以很准确的描述摩擦焊热源形状并与实际结果比较吻合，可以精确地描述摩擦焊接过程搅拌头流场的形状和大小及方向，为耦合温度计算提供了准确的流体结果。本实施例能很好的描述热源处相变情况，在整个摩擦焊接过程中，考虑相变的应变计算是直观、准确的。综上所述，本实施例能很高效地计算出考虑相变的应变结果，并给予直观的显示。
以上所述的具体实施方式，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施方式而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。