免费访问
金博宝问题
金博宝
体积19日,数 3.,二千零一十八
文章编号 308
数量的页面(年代) 15
DOI https://doi.org/10.1051/meca/2018021
在线发布 2018年9月12日

©AFM,EDP科学2018

1引言

在过去的几年里,大量的研究集中在材料的动态力学性能上。构成法,内力与变形的关系,也受到了很多关注。在大变形和高变形率下,如成形,加工过程和影响分析,弹塑性本构关系对预测材料的力学性能具有重要意义。

虽然在商用非线性有限元软件ABAQUS中已经实现了多种本构关系。[1,在某些情况下,它们不能很好地适应材料的某些复杂行为,例如涉及大应变的变形过程,高应变率,温度软化,等。因此,ABAQUS为用户提供了仅实现本构流的能力,或者通过FORTRAN用户子例程自行构建完整的本构行为律:UHARD或VUHARD表示硬化流律,UMAT或VUMAT表示完整的行为律。虽然用户强化实现非常简单和容易,用户资料的实现需要更多的专业知识。VUMAT与UMAT的主要区别在于VUMAT子程序不需要计算所谓的一致雅可比矩阵,但只有应力张量的计算西格玛1在增量的末尾。通常用于编写Vumat子例程的方法,由于显式时间积分方案使得时间增量很小,是否使用等效塑性应变增量的直接计算不需要任何局部迭代算法[3.]。

在本文中,我们选择了一种不同的方法,在VUMAT子程序中,塑性应变增量采用隐式径向回归算法第二节。为了说明所建议的方法并能够验证它,约翰逊-库克本构流动定律[,介绍于第三节,选择通过VUHARD和VUMAT子例程实现到Abaqus/显式代码中。文中给出了一个VUMAT子程序对Johnson-Cook弹塑性本构律的完整实现第四节。最后研究了该算法的数值效率和精度,在里面第五节,通过一些动态基准测试,证明了所提出的集成方案的有效性和鲁棒性。第五节并给出了泰勒冲击试验的一些替代本构律结果。

的时间积分算法J可塑性

基于有限差分法的时间积分算法在力学有限元求解中得到了广泛的应用。在这些算法中,时间在有限网格上离散,网格上连续点之间的距离定义为时间步长Δt。加上位置和一些时间导数t,集成方案稍后计算相同的数量t+γδt。通过过程的迭代,这些量的时间演化是可以追踪的。

二点一J塑性模型

弹性应力应变方程通常用以下增量形式表示:(1)在哪里为考虑大变形客观性的客观应力率,:是双点积,为各向同性的线性四阶弹性张量,由:(2)二阶同一张量,四阶单位张量,⊗二价的产品,G剪切模量和K体积模量与杨氏模量有关e泊松比V通过以下:(3)是所谓的弹性部分的变形率张量一般假设[]次弹性材料应符合以下加性分解方程:(4)在哪里是变形率张量的塑性部分。根据文献,J可塑性,也叫做冯·米塞斯屈服准则,由Huber于1904年和von Mises于1913年提出[]。广泛应用于机械工程和金属成形领域。它假设存在一个标量屈服函数f由:(5)在哪里为von Mises等效应力,或有效应力,由偏应力定义为:(6)西格玛y是材料的当前屈服应力,在我们的应用中取决于等效塑性应变,等效塑性应变率和温度T由:(7) (8)在哪里ε是泰勒-昆尼吗[7]定义塑性功转化成热能的潜热系数,Cp是比热系数和π为材料的密度。所谓的屈服函数f定义屈服面(当f=0)弹性域(当f<0)和无访问域(当f> 0)。假设有一个联想塑性流动规律,塑性应变率可以用以下公式表示:(9)在哪里γ一个标量是否表示流强度和为二阶张量(流应力法向单位,仅由试验弹性应力确定[8])表示以下给出的流动方向:(10)

从方程(7)可以很容易地得到:(11)

上面描述的塑性模型必须与增量目标算法进行及时集成。柯西应力张量时间导数的不同表达式用于方程(1)可以使用,在最终结果中导致一些明显的差异。Abaqus/Explicit,例如,对于实体元素,对所有内置本构模型使用jaumann应力率,对于Vumat用户子例程使用green–naghdi应力率。[1]。当我们希望将VUMAT结果与内置模型进行比较时,这可能会导致一些困难,我们将在后面看到这一点。

2.2径向回归映射算法

一种应用广泛的时间积分方法J各向同性硬化塑性是威尔金斯[]还有麦肯和萨克斯[]。这种算法现在被广泛使用,并在许多书籍中详细介绍,如[11]或者像西蒙等人的论文。(8]或Ponthot [十二]。我们假设,由于有限元公式,时间之间的应变增量t0(最后一个平衡态)和t1(新的平衡配置t1=t0+γδt)是以及当时的偏应力张量t0。应力张量的试偏部分最后的静液压p1是从应变增量计算的吗假设整个步骤第一次完全有弹性,以便:(12)(13)

符合方程(5),在时间t1,屈服函数f因此:(14):(15)电流增量开始时最后平衡态的屈服应力(t=t0)为了知道预测的试验应力是否允许,必须对屈服函数的符号进行测试。f这两种选择如下:

  • 如果f小于0,整个步骤是完全弹性的,这意味着预测的应力是可以接受的,我们可以得出结论吗

  • 如果f>0,试验应力是不允许的,为了计算

利用离散广义一致性参数对塑性修正进行了计算f(Γ)= 0的当前增量(这时间t=t1):(16)(17)

方程组组合(16)(17)导致一致性参数的最终形式如下f(Γ):(18)其中唯一未知的值是标量_,已经在方程中定义(11)。如果我们假设在这段时间内材料呈线性硬化,由于在Abaqus/ explicit中使用了显式集成方案,方程的下列解析解(16)可以写,如Gao等人提出的。(]或《Abaqus手册》[3.](19),此时硬化律的斜率,在时间增量Δ假定为常数t。由于Vumat实现的简单性,通常采用这种方法,因为稍后不需要任何迭代来解决所提出的问题。不幸的是,因为它将被进一步介绍,当本构流规律的非线性项变得重要时,本构流方法所提出的近似会导致许多不稳定性。在我们的方法中,因此我们选择了解方程(18)使用类似于Zaera等人提出的方法。(13解方程的一种所谓的求根方法(11)将在2.3节

一旦正确的最终值的塑料校正Γ已经获得,应力张量偏分量的终值更新后的内部变量通过方程计算得到(17)。在算法的最后,增量末端的最终应力计算公式如下:(20)

2.3非线性方程的求根

几种求根方法可以用来求方程的数值解。(18)。在我们的应用中,在正确选择了搜索区间的边界后,我们知道请求的根被括在给定的区间内。一旦我们知道建议的区间包含一个根,有几个经典的程序可以用来改进它。在提出的方法中,其中一些已知有快速收敛速度,精度及/或稳健性[十四15]。不幸的是,保证收敛的方法是已知的较慢的方法,而那些表现出快速收敛速度的人,如果不采取措施来避免这种行为,也可以在没有警告的情况下迅速冲向无穷大。在众多可用方法中(弦,二分法,Regula-Falsi,骑手',牛顿-拉斐逊)由于效率和鲁棒性,最后只选择了其中的两个,下面将介绍:对分法。缓慢而确定的方法和牛顿-拉斐逊方法是快速但有时失败的方法。

2.3.1等分法

本文提出的第一种找根方法是对分法。这种方法是不能失败的,但是收敛速度很慢。这是一种迭代寻根方法,在这种方法中,根的预测间隔连续减半,直到足够小,也叫区间减半法。重复此过程,直到间隔足够小,满足:(21)在哪里γ国际清算银行为二等分法的误差容忍度。假设要找到的根最初被一个区间括起来[x0x1令人满意f(x0)f(x1)< 0,该方法线性收敛,相对来说比较慢,但是这种方法的主要优点是只需要对函数进行评估f(x)而不是它的导数f(x)和牛顿-拉斐逊等其他方法一样。因此,当计算函数的导数时,这种方法会变得很有竞争力。f(x)变得难以计算;即当f(x)是复杂的。

2.3.2 Newton-Raphson方法的安全版本

牛顿-拉弗逊法[十四15]是最著名的求根方法,因为它简单高效(收敛速度非常快)。这种方法唯一明显的缺点是需要对导数求值f(x)函数的f(x),所以它只在问题中有用f(x)容易计算,或者我们会进一步看到数值计算。从泰勒级数展开f(x关于x,我们得到如下表达式:(22)

牛顿-拉弗森是一个从第一猜测开始的迭代过程x0对于函数的根f(x),重复此过程,直到收敛准则为:(23)在哪里是牛顿-拉斐逊法的误差容限,是达到了。

在某些情况下,牛顿-拉弗逊方法的全局收敛性较差,因为切线并不总是函数的一个可接受的近似值。但更重要的是,当f(x)是一个分段函数,由给定点左右两个不同的表达式定义xSf(x)将不会是可微的xS,导致牛顿-拉弗森算法在求解时的数值困难和发散x穿过点xS。作为一个结果,所谓的牛顿-拉弗森的安全版本,描述于图1,本文将牛顿-拉斐逊方法与二分法相结合,提出了相应的计算方法。

如果区间中有一个根[x0x1令人满意f(x0f(x1)< 0,牛顿-拉斐逊方法的安全版本考虑到[x0x1当根和Newton-Raphson迭代的第一个猜测开始时。每一次迭代后,通过用新的解决方案替换这两个边界中的一个来更新间隔。如果迭代超出了间隔,它被忽略,取而代之的是一个二分法步骤来重新定位初始的猜测。牛顿-拉斐逊算法需要计算导数屈服函数f(Γ)对Γ参数求解方程(18),以便:(24)收益函数导数的定义如下:(25)

计算屈服函数导数的一种常用方法西格玛y关于T是使用分析方法,根据材料的硬化流动规律,计算出各偏导数的解析表达式。然而,对于大多数屈服函数,用解析法求导数可能比较困难。因此,在这种情况下,我们在这里提出一个数值解作为替代。这个数值解是通过在等效塑性应变上加一个小增量得到的。分别用塑性应变速率和温度来计算这三个偏导数在方程(25):(26) (27) (28)

当然,准确的结果取决于三个增量的正确选择和ΔT。在以下所有数值试验中,三个增量任意固定在相同的值Δ上。x=10−6

缩略图 无花果。1

安全nr算法流程图。

3硬化流动规律

3.1约翰逊-库克硬化流动定律

在本文中,约翰逊-库克硬化流动定律[16]被选择通过VUMAT和VUHARD子例程在Abaqus/Explicit中实现。虽然该方法可以应用于其它弹塑性本构定律的求解,约翰逊-库克法是验证所提议的方法的最佳选择,因为它已经在Abaqus/显式代码中实现,这意味着可以在Abaqus内置的本构法和我们的FORTRAN实现之间进行基准测试。

Johnson-Cook硬化流动定律可能是考虑塑性应变的高应变率变形过程模拟中应用最广泛的流动定律,塑性应变率和温度效应。由于过去对许多材料的本构流动规律参数的识别做了大量的工作,它在许多有限元代码中实现,如abaqus[1]。一般公式由下式给出:(29)在哪里为参考应变率,T0T分别是材料的参考温度和熔化温度,以及一个BCn是五个本构流动规律参数。该流动定律的乘法公式允许采用下列形式:(30)哪里,根据(1317,塑性应变率的依赖性只有在:(31)以及对温度的依赖定义为,如果T<T0屈服应力与温度无关,如果T±T假定该材料表现为液体:(32)

由方程定义的这两个条件(31)(32)导致硬化关系出现一些不连续,它的导数而屈服函数本身在迭代求解过程中产生了数值困难。因此,屈服函数在T0。尽管如此,文中提出了鲁棒牛顿-拉弗逊方法2.3.2节已被发现具有足够的健壮性来克服这些问题。约翰逊-库克流动定律导数的解析形式西格玛y关于T由以下三个方程给出:(33) (34) (35)

最后一个经常遇到的问题是,因为这个术语在方程(33)H当第一个塑性增量(即当)因此,在第一步塑性硬化参数的计算中,我们必须进行特殊处理。这将在后面进一步介绍第4.3节

3.2 Abaqus/Explicit中的VUHARD实现

VUHARD子程序是在Abaqus/Explicit中实现一种新的本构流动规律的简单方法,它只需实现一个FORTRAN子程序来计算材料的屈服应力及其衍生产品T。内建本构定律的主要部分用于应力的时间积分。对于给定的时间增量,并利用所提供的用户子程序计算硬化流动规律。在ABAQUS文档中很少给出有关此实现的详细信息,但Jansen van Rensburg等人提供了一些有用的信息。(十八]。

针对Johnson-Cook流的VUHARD子程序的数值实现是通过一个FORTRAN程序完成的,该程序定义了硬化流定律根据方程式(29)以及西格玛y关于T由方程的定义(33)- - - - - -(35)。这部分代码与本文下一节介绍的VUMAT子例程中使用的代码完全相同。

4 ABAQUS/Explicit中的Vumat实现

本节将详细介绍使用Fortran VUMAT子例程为Johnson-Cook流律实现径向返回映射算法的所有细节。详细的流程图算法如图所示图2。该算法的第一个块“VUMAT的开始”用于获取定义为用户材料常数的材料属性。abaqus/explicit为用户子程序vumat提供以下定义的数量:

  • 应变增量对于当前时间步长;

  • 应力张量西格玛0和温度T0在电流增量的开始处;

  • 时间增量Δt对应于当前时间步长;

  • 用于存储重要数据(如,并将它们从一个增量转移到另一个增量。

VUMAT子例程必须计算并返回应力值西格玛1以及每个积分点增量末尾的sdvs变量。为了计算模型中的温度,还必须评估内部能量和耗散能量。

缩略图 无花果。二

VUMAT实施流程图。

4.1 ABAQUS包子程序

ABAQUS软件中的包部分是用于计算起始值(计算时间增量Δ的值)的强制步骤。t例如,通过考虑材料的性能)作为参考时间t=0。为此,我们需要计算由虚拟应变增量引起的弹性应力。Abaqus使用以下表达式提供:(36)

计算应力西格玛1将在包步骤结束后丢弃。

4.2弹性预测

本小节的目的是计算预测弹性应力的偏差部分。和冯·米塞斯等效应力并与增量开始时的屈服应力进行了比较以测试当前步骤是完全弹性还是部分塑性。因此,我们需要首先分解初始应力西格玛0进入其静水压力(p0)及其偏差()部分:(37)

然后,我们计算新的压强(p1)和试应力张量的偏微分部分()使用:(38)并对增量开始时的屈服应力进行了比较对冯·米塞斯等效试验应力(见情商。(15)以后再说)。根据这个测试,然后我们会纠正或不纠正最后的重音。

  • 如果,塑料校正器为零,因此可以跳过塑料校正步骤。我们假设预测的应力是最后一个,塑料校正Γ= 0,屈服应力保持不变我们可以直接进行最终的计算4.4节

  • 如果,所述步骤至少是部分塑料和所述塑料校正器第4.3节必须进行计算,以便在电流增量结束时将预测应力拉回材料的屈服表面。

4.3塑料校正器

根据弹塑性本构行为规律,采用牛顿-拉弗森算法的安全版本进行应力恢复。这个牛顿迭代法用于计算Γ参数定义修正由于增加的压力。为了提高计算效率,我们初始化的值Γ参数的值年底前增加(Γ=Γ0)如果电流增量是第一个塑性增量(即如果),因为我们无法计算H(Γ)当塑性应变为零,我们初始化Γ一个初始值不同的零(Γ= 10−8)二等分部分的间隔初始化为。预测的等效塑性应变,塑性应变率和温度T1在增量结束时,通过方程组计算(17)。屈服应力,屈服函数f()及其衍生产品f(Γ)计算。然后,我们测试的牛顿迭代算法计算的收敛的增量ΔΓΓ参数:(39)并与牛顿-拉弗森精度进行了比较γNR用户定义:

  • 如果ΔΓ>γNR我们需要使用以下关系迭代计算_值的校正:

(40)重新计算。的值 f(Γ)和 f艾氏(γ);
  • 如果δ-γ-εγNR我们得到的最终价值Γ参数。

应力张量的最后一个偏微分部分根据预测值计算修正项_使用:(41),在径向回归算法的框架内,定义如下:(42)

4.4最终计算

最后计算的主要工作是更新状态变量(SDVs)和能量。最后一个应力张量在增量的最后西格玛1是由静水压力计算出来的吗p1以及应力张量的偏态部分由于方程(20)。等效塑性应变,等效塑性应变率和最终温度T1增量的末尾存储在各自的SDVs变量中,以便在下一个增量中重用。我们还需要计算新的比内能e1来自:(43)以及耗散的非弹性能来自:(44)

在这一点上,VUMAT子例程结束。在流程图的这一点上,最终的温度没有计算出来,由于软件使用后续的热步骤来评估由于非弹性能量耗散和传导部分而引起的温度升高。热力学能在这个热步骤中被修正,这似乎在这个额外的步骤中被考虑(在VUMAT子例程之外),因为对1个元素施加所有节点位移的测试给出了与内置例程完全相同的结果。

5建议实施的验证

在本节中,VUMAT子例程中编程的Johnson-Cook法的性能与Abaqus/显式内置Johnson-Cook法和VUHARD Johnson-Cook实现的性能进行了比较,以验证拟议的实施方法。基准测试包括两个不同的单元测试(拉伸测试和剪切测试),圆棒的颈缩和著名的泰勒冲击试验。同样的材料,一个42 crmo4钢铁、所有这些测试都被选中了,和材料性能的报告表1。建议的基准点使用“动态温度位移,Abaqus显式FEM代码的显式过程。非弹性热分数参数已设置为其默认值ε=0.9。这种热力耦合方式允许通过塑性耗散或粘弹性耗散产生热量。没有其他热边界条件适用于全球绝热情况。

所有基准测试都是在运行Ubuntu 16.04 64位、12 GB RAM和两个4核e5620 Intel Xeon处理器的Dell Precision T7500计算机上使用ABAQUS/Explicit V.6.14解决的。所有计算都是使用Abaqus的双精度选项完成的,with one CPU and the VUMAT and VUHARD subroutines have been compiled using the Intel Fortran 64 v.14 compiler.以下模型已分别为四个基准测试:

  • A-N-R模型:采用牛顿-拉斐逊程序的Vumat模型,利用方程(33)–(35)对导数进行分析计算;

  • N-N-R模型:采用牛顿-拉弗逊法建立VUMAT模型,利用式(26)-(28)对其导数进行数值计算;

  • 直接的模型:VUMAT模型的直接评价Γ使用方程(19)如其他文件所述[

  • VUHARD模型:仅由方程定义的约翰逊-库克本构流动定律(29)它的三个分析导数(33)–(35)已经使用fortran子程序实现;

  • 内置的模型:内置Johnson-Cook宪法的原生实现,以便将结果与参考溶液进行比较。

测试也采用了对分法,但是为了减少结果的数量,我们在这一节中选择省略它们,因为这些结果与A-N-R模型非常接近,但是计算时间比A-N-R和N-N-R模型长10-30%。在这一点上,必须注意的是,当使用内置实现和Vumat子例程时,abaqus/explicit不使用相同的目标应力率。如文献所述[1],并经我们基于单单元超弹性剪切试验的试验验证,Jaumann速率用于内置公式和VuHard子例程,而绿色Naghdi速率用于VuMat子例程。因此,当发生大旋转时,无法直接比较内置和Vumat结果。

表1

42CrMo4钢的材料性能。

5.1一个元素基准测试

单元素测试,也叫单元素检验,是一种非常简单和实用的方法来研究的准确性和灵敏度的一个元素的行为对外部负荷。在本节,4结点双线性位移和温度的变形,减少集成沙漏控制元素CPE4RT将使用建议的VUMAT子例程和Abaqus内置模型进行模拟。由于只使用了一个未集成的元素,我们只有一个积分点位于元素的中心。在这种测试中,该单元的所有节点都受到规定位移的约束。由于以下两个基准点的几何变化完全相同,很容易用塑性变形来比较结果,压力和温度。元件的原始尺寸为10 mm×10 mm。

5.1.1单元件拉伸试验

在第一个基准测试中,元素的两个左节点被封装,并具有指定的水平位移。d=10 mm应用于同一元素的两个右节点,如中所示。图3。当我们使用显式积分方案时,总模拟时间设置为t=0.01μs。一元拉伸试验的五个结果完全吻合。最终塑性应变值最后的温度T= 164.09°C。由于最终结果在不到1秒的时间内得到,因此本次测试的计算时间差异不明显。

缩略图 无花果。3.

一元拉伸和剪切试验的数值模型。

5.1.2单构件剪切试验

第二个基准与前一个类似,但是现在,元素的两个底部节点被封装,并具有规定的水平位移。d=10 mm应用于同一元素的两个顶部节点,如中所示。图3。再一次,在所有五个结果之间找到了完美的匹配,最终塑性应变值最后的温度T=192.22摄氏度。

5.2圆棒的颈缩

圆钢试件的颈缩试验有助于评价材料在塑性和大变形条件下的Vumat子程序性能。[十二8]。由于对称结构,建立了试件的轴对称四分之一模型。试样尺寸报告见图4。加载是通过施加沿…的总位移7mm来实现的试件左侧轴线,同一边缘的径向位移保持为零。在另一侧,轴向位移受到限制,而径向位移是自由的。网格由400个CAX4RT元素组成,在总高度的1/3右侧有200个元素的细化区域。再一次,总模拟时间设置为t= 0.01S由于采用显式积分方案。

图5分别为内置模型(左侧)和N-N-R模型(右侧)的变形杆等效塑性应变等值线图。最大等效塑性应变位于条的中心(因此,我们选择在无花果。四),所有模型给出的值都与所报告的值完全相同表2T图6显示了冯米塞斯应力的演变在不同模型的试件末端施加位移。如图所示,内置的模型,Vuhard模型和牛顿-拉斐逊模型的两个版本给出了几乎相同的结果。在直接模型中可以看到细微的差别,在模拟过程中对von Mises应力有一点高估。

表2报告计算结束时所有模型的增量总数和总计算时间的一些结果。如本表所述,两个版本的牛顿-拉斐逊模型的结果是相同的,这往往证明,在使用流动定律的数值或分析导数方面没有区别。这个T/列表示总计算时间与总增量数之比,即计算时间/增量。如报道表2,就计算时间而言,本机过程对应于最快的算法。这很容易解释,因为这个内部过程是在Abaqus代码中本地优化的。将所有数据传输到Vumat或VuHard子程序的简单事实,增加了计算时间,而与所实现的应力评估算法的优化无关。这就意味着由于使用VUHARD子例程而增加的CPU时间约为18%,而Vumat程序的不同版本在T/比率。

在这种情况下,值得注意的是,牛顿-拉斐逊模型执行整个模拟所需的增量总数比其他三个模型要少。这种差异在图7其中时间增量Δ的演变t在计算过程中报告。时间增量的平滑变化,和最大的值,为牛顿- raphson的两个版本,以最小的增量完成模拟。使用内置模型和VUHARD模型,从Δ减少稳定时间增量t=6.5×10−8年代Δt=5.2×10−8位移2.52 mm后发现S,使用直接法时,大量减少Δ稳定增加t=3.4×10−8s在位移2.03 mm后会产生一些残余振荡。使用牛顿-拉弗森模型的两个版本的时间增量都是恒定的——位移为3.36 mm,在这一点之后由于试件的收缩而平稳地减小。在最后三分之一的模拟时间内,N-R模型的增量几乎等于内置模型和Vuhard模型的增量。

从这些结果我们可以得出结论,本构方程的集成有一个影响稳定时间增量Δt。我们还可以注意到,根据报告的结果表2那使用VUMAT Newton-Raphson不会增加太多的总计算时间(在本例中大约增加10%),因为增量的总数减少了3.8%。我们还可以注意到VUMAT子例程比VUHARD子例程快。当然,优化FORTRAN例程可以降低VUMAT模型的计算成本(删除子例程中大量的测试,优化计算和数学表达式,等等)。还需要注意的是,Newton-Raphson子例程要求的精度γNR=10−8也会影响计算时间,减少它也减少了计算时间。

缩略图 无花果。四

圆棒颈缩的数值模型。

缩略图 无花果。五

等效塑性应变圆棒颈缩的轮廓图。

表2

对圆棒基准缩颈的结果进行了比较。

缩略图 无花果。六

•冯•米塞斯应力VS圆杆颈的位移。

缩略图 无花果。7

时间增量ΔtVS圆杆颈的位移。

5.3泰勒冲击试验

5.3.1轴对称泰勒冲击试验

最后,通过对Taylor冲击试验的仿真,验证了VUMAT子程序在高变形率下的性能[十九]。在泰勒冲击试验中,发射圆柱形试样以规定的初始速度撞击刚性目标。数值模型,报道图8建立为轴对称。高度为32.4 mm,半径为3.2 mm。当径向位移是自由的时,轴向位移被限制在试样的右侧(以使弹丸在目标上无摩擦的情况下实现完美接触)。预定的速度Vc对试样施加287 m/s。网格由250个CAx4RT元件(5×50元件)组成。泰勒冲击试验的总仿真时间为t= 8.0×105年代。

图9分别为内置模型(左侧)和N-N-R模型(右侧)的变形杆等效塑性应变等值线图。两种模型的等效塑性应变分布基本相同。最大等效塑性应变位于模型的中心元素(中的红色元素无花果。8),模型给出的值与文献报道的值相当一致表3T以及试样的最终尺寸Lf(最后的长度)Df(冲击面最终直径)图10显示了等效塑性应变的演变随着时间的推移,不同模型的元素在冲击面的中心(红色的元素在无花果。8)如图所示,根据表3,所有模型的结果几乎相同。

从这之后可以注意到,直接模型(即假设显式的时间积分方案允许使用对塑性应变的直接评估,那么编写VUMAT子程序的经典方法不会根据关于试件内部场和最终几何形状的其他模型提供结果。这往往证明,如果要根据泰勒冲击试验和Fortran Vumat子程序对动态应用中的材料本构参数进行反向表征,以实现奇异本构定律,这可能是本构方程参数识别的误差来源。

缩略图 无花果。8

泰勒冲击试验模型。

缩略图 无花果。九

等效塑性应变泰勒冲击试验的等值线图。

表3

泰勒冲击试验结果的比较。

缩略图 无花果。十

等效塑性应变VS泰勒冲击试验的时间到了。

5.3.2三维泰勒冲击试验

为了在泰勒冲击试验的数值模拟框架内保持较长的计算时间,选择了3D建模。因此,泰勒圆柱试样的三维四分之一模型,与泰勒二维模型的尺寸相同,如报告中所述图11。新的网格由37 422个C3D8RT单元组成,这是一个相当大的模拟基准。

有关CPU时间的结果(增量,时间和T/比率)变形和温度的最终值以及几何结果如所示。表4。我们可以注意到关于几何尺寸的结果有很好的相关性,然而,由于Vumat子程序更好地考虑了本构定律,变形和温度方面的结果(非常受约束元素上的非常局部化结果)不同(这已通过细化网格得到证实)。类似于5.2节关于圆棒的颈缩,的T/对于VUMAT子例程的不同版本,然而,由于模拟所需的增量总数减少了大约5%,CPU时间的总体增加仅为5%。关于Vuhard方法,增加T/比例约为10%,所需增量与内置方法的增量大致相同,导致总的CPU也增加了大约10%。这再次证明了VUMAT算法的有效性。

缩略图 无花果。11

三维泰勒冲击试验模型。

表4

三维泰勒冲击试验结果比较。

5.4 Vumat实施替代本构法

在本节中,在N-N-R Vumat和直接Vumat子程序中实现了Johnson-Cook本构律和其他三种本构律,以模拟泰勒压缩试验。为了验证该算法在替代弹塑性本构关系中的应用J塑性和各向同性硬化有一般形式。使用n-n-r方法的主要优点是不需要显式计算硬化定律的导数,它们是用公式(26)-(28)进行数值计算的。其他的本构法是坦赫本构法,修正tanh本构法和b_ker本构法。

  • 卡拉马什等。(20.]提出了所谓的TANH本构定律,在Johnson-Cook本构定律中加入了应变软化的模型项,由:

(45) (46)在这 西格玛 JC代表最初的约翰逊-库克宪法。常数 可以调节峰值应力对应的应变, p为附加本构律参数,和 T 矩形为应变软化现象的起始温度。
  • 霍尔等。(二十一]后来通过修改丹的本构法提出了本构法,由:

(47)在哪里 (48) (49)一个BC 1C 1 np为本构律参数。和Tanh模型一样,常数 能调节峰值应力对应的应变吗 T 矩形为应变软化现象的起始温度。
  • 我们使用VUMAT子例程实现的最后一个本构律是Baker [二十二二十一并给出:

(50) (51)在哪里 一个n 0是与温度有关的材料参数,和 C是常数。

仿真所用的二维模型与文中介绍的二维模型相同5.3.1节。在这种情况下使用的材料是42CrMo4钢和铁素体珍珠岩(称为42CrMo4 FP)。Hor等人提出了42CrMo4 FP钢的四种本构关系参数。(二十一]。将预定义的速度设置为Vc= 200 m / s。四种模型变形杆的等效塑性应变等值线图如所示。图12

塑性变形集中在底部,最大等效塑性应变位于试件的中心单元。有关增量总数的更多详细信息,总计算时间,最大等效塑性应变,最后的长度Lf,最后的半径Rf底部和最高温度T有报道表5。Johnson-Cook和TANH模型的结果几乎相同。与前两种模型相比,改进后的TANH模型在底部实现了较大的等效塑性应变。相应地,它也实现了更大的最终半径和更高的底部温度。然而,正如我们所看到的,另一部分的塑性变形较小。无论用最大等效塑性应变和温度还是几何响应,贝克模型的塑性变形最小。关于直接进近,可以看出,这一个通常高估了塑性变形和温度树第一模型。就全球结果而言,前三个模型给出了几乎相同的结果,而Baker模型给出的结果与之前的模型有很大的不同。

缩略图 无花果。十二

四种本构模型泰勒试件的等效塑性应变等值线图。

表5

泰勒压缩测试的结果。

5.5基准讨论

在所有拟议的基准测试中,值得注意的是,Vumat子程序和ABAQUS内置模型的结果非常接近,但可以指出一些细微的差异。这些差异主要可以用以下几点来解释:

  • 通过Vumat子程序实现的Johnson–Cook本构定律的表达式与ABAQUS/显式内置模型的表达式不完全相同,因为我们的实现中已经修改了对变形率的依赖性,如前所述。第三节。这可以在酒吧缩颈基准中看到,Vuhard和内置模型不能提供完全相同的结果;

  • 内建模型和VUHARD模型通过明确的中心差分时间积分规则进行积分,而径向回归法,属于隐式积分算法,在VUMAT子程序中使用;

  • 由于函数的根不是一个精确解,但是一个近似值,精度公差的选择对最终结果有影响;

  • 内置模型和VUHARD子例程使用Jaumann应力速率,VUMAT例程使用Green-Naghdi应力速率。当材料点的有限旋转伴随有限剪切时,这会导致大变形模型的结果不同。[1]。

6结论

在本文中,提出了一种通过Vumat子程序实现ABAQUS/显式弹塑性本构关系的方法。建议的实施是稳健和有效的,可用于以下任意弹塑性本构模型所定义的各种材料的行为模拟J塑性和各向同性硬化。在所提出的算法中,采用基于径向回归法的隐式时间积分方法,这是无条件稳定的,可以保证在计算过程中始终满足von Mises屈服准则。通过数值基准验证,一种鲁棒高效的寻根方法,所谓的牛顿-拉斐逊方法的安全版本,已经实现了计算塑料校正器项。

该算法现在可以重复使用,只需改变流定律的表达式,以实现替代本构定律,如Rule等人提出的Johnson-Cook法修订版。(16,或同一家庭的任何其他人。这类发展的应用主要是基于动态冲击试验的本构方程参数反演识别。或者在abaqus/explicit中编程替代本构定律。在直接评估塑性应变率的基础上,我们发现,对于其他类型的Vumat或Vuhard程序,建议的实施更为精确,以便为识别程序提供准确的结果。

确认

这项工作部分得到了中国奖学金委员会(CSC)第号奖学金的资助。201406290010。

工具书类

  1. Abaqus诉6.14文档,达索系统公司[谷歌学者]
  2. C。高Abaqus/显式程序中VUMAT热粘塑性本构模型的有限元实现在:计算力学,施普林格,2007,页。301—301(CrossRef) [谷歌学者]
  3. Abaqus诉6.14用户子程序参考指南,达索系统公司[谷歌学者]
  4. G.R.约翰逊,W.H.做饭,大应变金属的本构模型和数据,高应变率和高温,第七届国际弹道研讨会论文集,卷。21日,海牙荷兰,1983,页。541 - 547[谷歌学者]
  5. 年代。Nemat Nasser关于有限变形弹塑性,INTJ。固体结构。18(1982)857–872(CrossRef) [谷歌学者]
  6. M.-H。于广义塑性,斯普林格科学与商业媒体,二千零六[谷歌学者]
  7. 特种部队泰勒,H。Quinney,冷加工后金属中剩余的潜在能量,PROCR.SOC。朗德爵士。一个,143(1934) 307-326(包括数学和物理性质的论文)(CrossRef) [谷歌学者]
  8. J.C.SimoT.J.R.休斯计算不弹性,施普林格,一千九百九十八[谷歌学者]
  9. M.L.威尔金斯,弹塑性流动计算,技术。代表。,在:阿尔德(E.)计算物理学的方法,卷。三,学术出版社,1963,页。211 - 263[谷歌学者]
  10. G.Maenchen,年代。麻袋,张量的代码,在:桤木年代。FernbackMRotenberg(编辑)计算物理方法,卷。三,学术出版社,纽约,1964年,页。181—210[谷歌学者]
  11. F。唐恩N。Petrinic,计算塑性导论,牛津大学出版社,纽约,二千零五[谷歌学者]
  12. j。Ponthot,大变形弹塑性和弹粘塑性过程数值模拟的统一应力更新算法,INTJ。体。18(2002)91 - 126(CrossRef) [谷歌学者]
  13. R.Zaera,J。Fernandez-Saez,绝热和有限变形条件下热塑性本构方程积分的隐式一致算法,INTJ。固体结构。43 (2006)1594 - 1612(CrossRef) [谷歌学者]
  14. F.S.阿克顿有效的数值方法,共,1990[谷歌学者]
  15. W.H.出版社,数字配方,在:科学计算的艺术,第三版,剑桥大学出版社,2007[谷歌学者]
  16. W.K.规则,年代。琼斯,Johnson-Cook强度模型的修正形式,INTJ。影响Eng。21 (1998)609 - 624(CrossRef) [谷歌学者]
  17. l施沃,Johnson-Cook本构模型的可选应变速率形式及参数0的作用,在:第六届欧洲LS-Dyna用户会议记录中,Anwenderforum法兰肯塔尔,德国2007[谷歌学者]
  18. G.Jansen说伦,年代。角,基于状态变量的可塑性教程:一个abaqus UHARD子程序[谷歌学者]
  19. 特种部队泰勒,材料在高载荷下的试验,J。安特CIV。Eng。26 (1946)486 - 519(CrossRef) [谷歌学者]
  20. MCalamaz,D。CoupardF。Girot基于应变软化本构律的钛合金干切削过程数值模拟马赫。SCI。抛光工艺。14 (2010)244 - 257(CrossRef) [谷歌学者]
  21. 一个。何珥F。莫雷尔J.-L。勒布伦G.杰曼建模,在大应变率和温度范围内现象学本构定律的识别和应用,机器。板牙。64 (2013)91 - 110(CrossRef) [谷歌学者]
  22. M巴克,高速切削力的有限元模拟,J。板牙。PROC抛光工艺。176 (2006)117 - 126(CrossRef) [谷歌学者]

将本文引用为:L明OPantale,ABAQUS/显式有限元程序中弹塑性本构关系的有效和稳健Vumat实现,金博宝力学与行业十九,308 (2018)

所有表

表1

42CrMo4钢的材料性能。

表2

对圆棒基准缩颈的结果进行了比较。

表3

泰勒冲击试验结果的比较。

表4

三维泰勒冲击试验结果比较。

表5

泰勒压缩测试的结果。

所有的数据

缩略图 无花果。1

安全nr算法流程图。

在文本中
缩略图 无花果。二

VUMAT实施流程图。

在文本中
缩略图 无花果。3.

一元拉伸和剪切试验的数值模型。

在文本中
缩略图 无花果。四

圆棒颈缩的数值模型。

在文本中
缩略图 无花果。五

等效塑性应变圆棒颈缩的轮廓图。

在文本中
缩略图 无花果。六

•冯•米塞斯应力VS圆杆颈的位移。

在文本中
缩略图 无花果。7

时间增量ΔtVS圆杆颈的位移。

在文本中
缩略图 无花果。8

泰勒冲击试验模型。

在文本中
缩略图 无花果。九

等效塑性应变泰勒冲击试验的等值线图。

在文本中
缩略图 无花果。十

等效塑性应变VS泰勒冲击试验的时间到了。

在文本中
缩略图 无花果。11

三维泰勒冲击试验模型。

在文本中
缩略图 无花果。十二

四种本构模型泰勒试件的等效塑性应变等值线图。

在文本中

当前使用指标显示文章视图的累计计数(全文文章视图,包括HTML视图,PDF和ePub下载,并在Vision4Press平台上对视图进行抽象。

数据与2015年后在模板上的使用情况相对应。当前使用指标在在线发布后48-96小时内可用,并在每周更新。

度量标准的初始下载可能需要一段时间。