GISSMO损伤失效模型对金属材料失效的预测
摘要：本文主要阐述一种基于应力状态的损伤失效模型（GISSMO）对金属材料失效的预测，GISSMO在预测材料成形和碰撞分析中显著提高了预测失效的精度。该模型具有以下的特点：依赖于当前应力状态的失效曲线进行损伤的累计，损伤变量是一个考虑了非线性加载的标量值，当损伤变量达到1时材料失效；同时考虑了材料局部的不稳定性，通过一条与当前应力状态相关的曲线对材料的局部不稳定性进行预测，当不稳定因子为1时，应力张量开始于损伤变量进行耦合。在宏观上反映出超过临界区域后（即应变局部化到断裂产生）对应力的影响。在模拟材料的行为时，网格的虚假依赖性也是需要考虑的，因此引入了与网格相关的缩放曲线（正则化曲线），对不同的网格尺寸时的失效曲线进行缩放，以弥补网格尺寸对失效的影响。本文对汽车中使用的钢板进行了不同应力状态的实验，并使用GISSMO损伤失效模型对该材料的失效进行预测，展现出GISSMO在预测塑性材料失效时的可靠性。
关键字：应力状态、损伤失效模型、损伤变量、正则化曲线
1、介绍
在过去的很长一段时间内，精确的预测材料的失效行为一直是一个困难复杂的问题，塑性材料的失效需要考虑到很多的因素，同时考虑到这些因素并不容易。
Bridgman[1]等人通过实验发现，金属材料的失效直接与加载在试样上的静水压力有关。引入标量化的数值描述一点处的应力状态有着重要的意义。通过定义应力三轴度和lode角参数即可以确定任意一点处的空间应力状态，在平面应力状态下，应力三轴度和Lode角参数不是相互独立的，在很多基于现象学的失效准则中引入了应力三轴度这一概念，例如Johnson-Cook模型。
屈服轨迹和塑性流动对金属材料失效也有很大的影响，大的塑性变形对塑性材料的退化有着明显的影响，塑性变形是路径相关的。非比例的应变加载路径对失效的行为也有很大的影响，特别在汽车碰撞应用中，碰撞过程高度非线性，有必要掌握整个非比例加载的过程。
网格的虚假依赖性是失效分析时要考虑的问题，在局部紧缩之后，材料的变形往往集中在很小的区域，大的网格尺寸会导致结构偏刚，而在碰撞分析时，所需的网格尺寸较大，需要将失效曲线进行缩放。
基于失效准则的这些需求，Dynamore的NeuKamm[2-10]等人开发出了GISSMO（Generalized Incremental Stress State Dependent Damage Model），尝试考虑到上述提出的所有的在现象学中涉及的方面。
2、应力状态介绍
Bridgman通过实验发现静水压力与失效应变相关，应力状态对失效的预测有很大的影响，因此将应力张量分解为影响材料体积变化的静水应力和影响材料形状改变的应力偏张量两项：
式2.1即为应力偏张量，
为克罗内克函数，
为平均应力。
Bridgman通过实验证明了在塑性过程中体积不会改变（即等容变换条件），失效一般是和压力相关的，学者Rice和Tracey等人提出了表征应力状态的标量公式“三轴度”，定义为平均应力和等效应力的比值：
(2.4)
式中，
为静水压力，
为等效应力。
冯.米塞斯的屈服准则和应力-应变关系在弹性和塑性区域基于应力偏张量的第二不变量，其中，的表达式如下：
(2.5)
对于大多数金属材料，屈服只和偏应力的第二不变量有关，屈服方程表示为：
(2.6)
当
时，材料处于弹性状态；当时，满足屈服条件，材料处于塑性状态；应力状态
是不能实现的。
我们可以使用单轴拉伸试验中的屈服应力
来确定k，即=
。得到：
应力三轴度在处理金属材料失效时有着重要的意义，当材料在三个相互垂直的方向上施加相同强度的拉力，这种情况下只会改变材料的体积，材料的形状不会发生改变，等效应力为0，平均应力等于施加在三个方向上的应力，应力三轴度趋于+∞；如果在三个方向上施加相同强度的压力，应力三轴度趋于-∞。其他的中间应力状态都在这个边界内。要确定空间一点处的应力状态，lode角参数也是需要，lode角参数的表达式为：
在平面应力状态下，各向同性材料，lode角参数和应力三轴度可以相互转换，并非独立变量，转换关系为：
平面应力在主应力空间可以表示为：
带入公式2.3、2.4、2.8得：
可以得到纯剪切时的应力三轴度为0，单轴拉伸时为1/3，等双轴拉伸时为2/3。如果为负，三轴度的下限为-2/3。
3、GISSMO模型描述
金属材料失效通常依赖于应力状态，可以定义失效为的曲线，对于比例加载工况，失效应变在这条曲线上；对于非比例加载工况，即应变状态在加载过程中发生变化的情况，需要采用增量累计损伤的形式记录损伤变量随应变路径的演变。
比例加载，损伤可以以全量的形式表示：
式3.1对时间微分，应力状态定义为恒定时，可以得到下面的表达式：
对式3.2进行调整，可以得到：
公式3.4与具体的应变值无关，只与当前的损伤和应力状态有关，当损伤变量D越接近1时，损伤演化越快，可以用来对非比例加载进行损伤的累计预测。
类似于损伤变量，定义一个和不稳定性测量相关的附加的变量F，是一个于临界塑性应变相关的方程，在比例加载时，类似于损伤D：
进行微分，增量形式的不稳定变量：
当F=1时，应力和损伤变量进行耦合，耦合公式为：
式中，为未损伤的应力，为当F=1时的损伤变量值。
在材料达到局部不稳定阶段后，不同网格尺寸的模拟结果差异比较明显，主要体现在小尺寸网格模型的局部塑性应变比大尺寸网格模型的局部塑性应变大，大尺寸网格模型的失效会靠后，引入网格正则化曲线以最小尺寸的失效曲线为基准，使用其他网格尺寸时，对这个基准的失效曲线进行缩放，使得不同尺寸模型失效前吸收的能量基本一致。
4、相关试验及参数获取
在汽车行业中，高强钢的使用特别广泛，本文针对一种复相高强钢开展相应的力学性能研究，以提高零部件仿真的精确性。复杂应力状态的失效行为很难通过单一的试验表征出来，需要进行多种不同应力状态的试验拟合得到表征损伤失效的失效曲线，这就意味着要开展不同应力状态下的实验，对材料的失效判据进行标定。在LS-DYNA软件中，通过与材料本构模型耦合使用的关键字*MAT_ADD_EROSION卡片中的GISSMO失效模型，能够灵活的反映出复杂应力状态下损伤的累计过程，以提高在零部件对标式的精度。
4.1、实验样件与实验数据介绍
假定材料为各向同性材料，通过准静态的拉伸试验得到材料的本构曲性，试验样件为标准样件，为了反映复杂的应力状态，制定了多种不同形状的缺口样件，如下图所示：
a)标准样件
b)缺口D2.0mm样件 c）45°剪切样件
图1 不同形状样件尺寸
试验数据如下：
a)单轴拉伸力-位移曲线 b)缺口D2.0mm力-位移曲线
c)45°剪切力-位移曲线
图2 试验数据
4.2、本构曲线拟合
为了得到LS-DYNA中*MAT_024材料需要输入的有效应力-有效应变曲线，首先需要对试验数据工程应力-应变曲线进行处理，将其转化为真实应力应变曲线，并且在颈缩后，通过不同的硬化公式进行插值外推，得到完整的真实应力应变曲线，进而代入LS-DYNA中进行计算，模拟拉伸试验过程，通过仿真计算，得到力-位移曲线或者工程应力-应变曲线，并与试验数据进行对比，看结果是否一致，如果一致，说明输入曲线能够反应材料的真实性能，如果不一致，就需要调整输入曲线，直到结果与试验数据一致。
可以将仿真过程具体化为三步：
第一步，材料在弹性变形阶段，其应力和应变成正比例关系（即符合胡克定律），其比例系数称为弹性模量。因此对工程应力应变曲线的线性段进行线性拟合，得到的斜率即为弹性模量；
第二步，将工程应力-应变曲线转换为真实应力-应变曲线，找到紧缩点，删除紧缩后的数据，使用swift硬化公式进行外推，得到材料紧缩后的真实应力应变曲线，公式具体为：
第三步，根据第一步和第二步计算的弹性模量和真实应力应变曲线，采用LS-DYNA中的24号材料模型，并建立与试验一致的样件模型，设置与试验一样的边界条件（一端固定，一端加载），通过对计算结果设置，得到工程应力应变曲线，并与试验结果进行对比，如果一致，则说明输入的材料曲线（真实应力-应变曲线）能够反应材料性能，可用于工程实际，如果不一致，调整输入的材料曲线，重复计算，直至得到计算曲线与试验曲线一致的结果。最终得到的输入曲线以及模拟和试验的对标曲线如图3所示。
a)输入曲线 b)试验和模拟对比曲线
图3 单轴拉伸输入曲线和对标曲线
4.3 失效参数拟合
GISSMO模型中需要确定的参数包括：失效预测曲线，不稳定性曲线，损伤累计指数n和应力衰减指数m。根据多组不同应力状态的试验拟合出这些失效参数是一个复杂的过程，工作量大，使用不断试错的方法很难得到与实验结果匹配较好的模拟结果。随着计算机科学的发展，使用数值优化的方法进行参数识别是一个比较方便的方法。LS-OPT能够与LS-DYNA进行交互，操作简单，因此本文采用LS-OPT进行失效参数的识别，流程图如下：
图4 GISSMO材料失效参数识别流程图
在拟合过程中，可以首先设置较多的变量，进行灵敏度分析，将对结果影响较小的变量去掉，然后再对这些变量进行域缩减的优化，最终得到和试验数据匹配较好的一组失效曲线：
图5 失效曲线和不稳定性曲线
图6 模拟和试验结果的对比
在工程应用中，模拟材料成形的网格尺寸为0.5mm，而进行碰撞的网格尺寸为5mm，要将成形过程汇总的塑性应变、损伤等映射到碰撞的模型中，需要研究网格尺寸对结果的影响，不同网格尺寸对应的失效结果不一致，引入正则化曲线，建立不同网格尺寸的模型，得到匹配较好的结果。
a) 不同网格尺寸的模型 b) 不同网格尺寸的结果
c) 正则化曲线 b) 正则化处理后不同网格尺寸的结果
图7 正则化曲线和对应的模拟结果
5、结论
GISSMO能够精确的预测复杂应力状态下金属材料的失效行为，以及在非比例加载的路径下损伤的累计过程，以及处理网格尺寸的虚假依赖性，通过实际的例子，仿真的结果匹配较好，失效参数能够为后续的分析使用，从而为 工程师选取合适的材料进行结构的轻量化设计提供有力的支撑。
参考文献：
[1] Hill R . Studies in large plastic flow and fracture[J]. Journal of the Mechanics & Physics of Solids, 1953, 1(2):146-146.
[2]Gurson, A. L . Continuum Theory of Ductile Rupture by Void Nucleation and Growth: Part I—Yield Criteria and Flow Rules for Porous Ductile Media[J]. Journal of Engineering Materials and Technology, 1977, 99(1):297-300.
[3] Lemaitre J . A continuous damage mechanics model for ductile fracture[J]. Transactions of the Asme Journal of Engineering Materials & Technology, 1985, 107(107):83-89.
[4]Neukamm, F, Andrade, et al. An incremental stress state dependent damage model for ductile failure prediction[J]. International Journal of Fracture, 2016.
[5] Neukamm F , Feucht M , Roll K , et al. On Closing the Constitutive Gap Between Forming and Crash Simulation[C]// International Ls-dyna® Users Conference. 2008.
[6] Basaran M , SD Wölkerling, Feucht M , et al. An Extension of the GISSMO Damage Model Based on Lode Angle Dependence. 2010.
[7] Lee B W , Jang J I , Kwon D . Evaluation of fracture toughness using small notched specimens[J]. Materials Science and Engineering A, 2002, 334(1-2):207-214.
[8] D Hörling. Parameter identification of GISSMO damage model for DOCOL 1200M : A study on crash simulation for high strength steel sheet components. 2015.
[9] Falkinger G , Sotirov N , Simon P . An investigation of modelling approaches for material instability of aluminum sheet metal using the GISSMO-model. 2015.
[10] Anderson D , Butcher C , Pathak N , et al. Failure parameter identification and validation for a dual-phase 780 steel sheet[J]. International Journal of Solids & Structures, 2017:S0020768317302846.
