ABAQUS_材料本构模型与编程
基于Abaqus子程序的高分子材料本构关系实现

基于Abaqus子程序的高分子材料本构关系实现作者:汪品红来源:《计算机辅助工程》2013年第05期摘要:对于高分子材料的仿真,业界一般使用经典的弹塑性本构模型来描述其应力应变关系,但其真实的应力应变关系与经典的弹塑性本构模型存在一定差异,从而导致仿真与实际测试之间的差异.Abaqus提供UMAT/VUMAT子程序接口,让用户可以自己构建新的材料本构模型.通过撰写新的材料本构子程序实现高分子材料应力应变关系在仿真中的准确描述,减少仿真与测试之间的差异.同时,在卸载段可以通过卸载标志符的选择定义不同的卸载路径,方便用户使用.关键词:高分子材料;本构关系; Abaqus; UMAT; VUMAT中图分类号: TB324; TB115.1文献标志码: B引言高分子材料在日常生活中有着广泛的应用,因此其不可避免地出现在仿真分析中.当前没有一种商业软件具有适合高分子材料的材料本构模型.Abaqus是一款优秀的商业软件,其提供的子程序接口UMAT/VUMAT允许用户根据使用需求自定义材料本构.[1]使用该方法,可有效解决在仿真中由于材料本构不适用而导致的仿真与实际测试差异过大的问题.1高分子材料本构一般描述方法业界通常使用弹塑性本构定义高分子材料的材料属性.屈服强度一般取材料曲线上第一个峰值点.弹性模量的取法有2种不同的方式:对于应力应变关系曲线有明显直线段的,以第一段直线的斜率作为材料的弹性模量(切线法);对于曲线没有明显直线段的材料,则使用原点与屈服点连成的直线的斜率作为弹性模量(割线法).2种方式与真实应力应变曲线的比较见图1.图 1高分子材料测试材料曲线与仿真曲线比较由图1可知,无论使用何种方式,仿真使用的应力应变曲线都与实际材料的应力应变曲线有较大差异.将切线法获得材料数据代入到手机电池盖三点弯曲中进行仿真,见图2,其仿真与测试力位移曲线在最高点的差异约为23%,见图3.对于手机等一些电子类产品,高分子材料的仿真非常重要.在跌落或弯折测试中,高分子材料的应力应变关系与弹塑性本构的差异造成仿真预测不准确,必须定义正确的高分子材料本构.2Abaqus VUMAT子程序Abaqus提供丰富的材料本构模型库,能够满足绝大多数仿真材料模型的需要;同时,还提供UMAT/VUMAT子程序接口,让用户可以用FORTRAN语言编程,自己定义需要的材料本构模型,对Abaqus材料库中没有包含的材料进行计算.几乎可以把用户材料属性赋予Abaqus 中的任何单元,其中UMAT用在隐式仿真计算中,VUMAT用在显式仿真计算中.由于隐式计算与显式计算的差别,导致UMAT与VUMAT也有一定的差异,但是经过简单的改写即可完成它们之间的转换.本文使用准静态仿真分析方法,属于显式求解,所以只介绍VUMAT.3高分子材料VUMAT本构介绍由图1可知,高分子材料的本构与弹塑性本构最大的差异在于弹性段是直线还是曲线.弹性段的路径也直接影响到卸载的路径.因此,对高分子材料本构的定义关键在于非线性弹性段的实现,即要根据当前的应力值实时获取下一增量步所用的弹性模量值.程序整体流程见图4.图 4程序整体流程3.1弹性段多段线性的实现在弹性段,程序根据弹性模量和泊松比计算应力增量.由于弹性段为非线性,需要根据应力或应变更新用于计算的弹性模量值,直至达到屈服点,因此需要在输入文件中输入材料真实应力应变曲线,通过查表计算的函数,根据当前应力σ所在的位置,计算当前的弹性模量.应力应变曲线输入时,输入格式为:用查表的方法,直到σn3.2卸载路径的选择屈服发生后,需要选择弹性模量参与相关计算,有2个作用:一是用来计算屈服后加载段的应力试探值(不对该增量步真实应力产生影响,只起对比判断的作用);二是用来作为屈服后卸载的路径(为实现不同卸载路径,在程序中设置一个flag位,其值由用户自己输入),用户可以根据实际的需要选择卸载的路径.如图4中,共设置3种卸载路径:沿切线卸载、沿割线卸载以及沿曲线卸载等.用户也可以根据需要增加其他的卸载方式.4子程序的验证为验证子程序是否能实现设计的功能,取一个1/8的网格模型进行单轴拉伸仿真,单元类型为C3D8R.输出其应力应变曲线,与材料真实应力应变曲线比较,见图5.图 5使用VUMAT后加载应力应变曲线与材料曲线对比使用VUMAT后,加载的应力应变曲线与材料测试得到的真实应力应变曲线完全重合,说明VUMAT可以完全反映材料在加载过程中的力学行为.在卸载过程中,分别实现沿弹性段的切线、割线以及曲线卸载.为进一步验证,将VUMAT用于图2所示的手机电池盖三点弯模型中进行仿真与试验对比.在使用弹塑性本构模型时,仿真与测试力位移曲线的最大差异约为23%,而引入使用VUMAT 编写的高分子材料本构后,其仿真与测试的差异减少到4.5%,见图6.从实际项目的验证结果看,使用VUMAT后电池盖测试的力位移曲线与仿真的力位移曲线基本重合,仿真与测试的差异也明显减小.将该本构应用于其他高分子材料和实际案例,其仿真精度均明显改善,也说明该子程序在实际工程中的适用性.图 6使用VUMAT后电池盖力位移曲线对比5结束语使用VUMAT子程序后,高分子材料在加载段的力学特性与测试的真实应力应变曲线一致,同时将其应用在工程实际问题上,也与测试曲线基本一致,验证该程序的适用性.由于高分子材料的卸载特性较为复杂,还需进一步研究,所以程序只给出3种方式供用户按照实际需求进行选择.参考文献:[1]庄茁,张帆,岑松,等. 基于Abaqus的有限元分析和应用[M]. 北京:清华大学出版社, 2009: 509512.(编辑武晓英)。
abaqus操作流程

abaqus操作流程Abaqus操作流程Abaqus是一款强大的有限元分析软件,广泛应用于工程、科学和研究领域。
本文将介绍Abaqus的操作流程,包括软件安装、模型建立、材料定义、边界条件设置、求解和后处理等步骤。
一、软件安装需要从官方网站下载Abaqus软件,并按照安装向导进行安装。
安装完成后,需要激活软件,通常需要输入许可证文件或者许可证服务器地址。
如果是学术版或者试用版,可以直接使用。
二、模型建立在Abaqus中,可以通过几何建模、导入CAD模型或者手动输入节点和单元来建立模型。
几何建模是最常用的方法,可以使用Abaqus CAE中的几何建模工具,例如绘制线、面、体等基本几何体,然后进行布尔运算、切割、倒角等操作,最终生成复杂的几何模型。
导入CAD模型需要将CAD文件转换为Abaqus支持的格式,例如STEP、IGES、ACIS等。
手动输入节点和单元需要了解节点和单元的类型、编号、坐标等信息,比较繁琐,不建议使用。
三、材料定义在Abaqus中,需要定义材料的力学性质,例如弹性模量、泊松比、屈服强度等。
可以选择预定义的材料模型,例如线弹性、非线性弹性、塑性等,也可以自定义材料模型。
自定义材料模型需要了解材料的本构关系,例如应力-应变曲线,可以通过实验或者理论计算得到。
四、边界条件设置在Abaqus中,需要设置边界条件,包括约束和载荷。
约束是指模型的某些部分不能移动或者旋转,例如固定支座、铰链等。
载荷是指模型受到的外部力或者压力,例如重力、风荷载、温度载荷等。
可以选择预定义的边界条件,例如固定支座、均布载荷等,也可以自定义边界条件。
自定义边界条件需要了解模型的物理特性和边界条件的作用方式。
五、求解在Abaqus中,需要进行求解,即求解模型的应力、应变、位移等物理量。
可以选择不同的求解器,例如标准求解器、隐式求解器、动态求解器等,也可以选择不同的求解方法,例如直接法、迭代法等。
求解过程中需要注意模型的收敛性和稳定性,如果模型不收敛或者不稳定,需要调整求解器和求解参数。
abaqus输入应力应变曲线本构

一、概述1.1 abaqus是什么Abaqus是一款用于进行有限元分析的软件,它可以模拟各种工程材料在不同应力应变条件下的力学性能,为产品设计和工程分析提供了重要的支持。
1.2 应力应变曲线本构应力应变曲线本构是描述材料在不同应变条件下的应力应变关系的数学模型,用于分析材料在不同加载条件下的强度和变形性能。
在abaqus中,输入应力应变曲线本构是十分重要的一步,它决定了模拟分析的准确性和可靠性。
二、 abaous输入应力应变曲线本构的方法2.1 材料的本构模型在abaqus中,材料的本构模型是用于描述材料在应力应变条件下的力学性能的数学模型。
常见的本构模型包括线性弹性本构、非线性弹性本构和塑性本构等。
对于不同的材料和工程条件,需要选择合适的本构模型来准确模拟材料的力学性能。
2.2 输入应力应变曲线在abaqus中,输入应力应变曲线本构的方法主要包括以下步骤:(1)确定材料的本构模型,选择合适的本构模型对材料的应力应变曲线进行描述。
(2)测定材料的应力应变曲线,通过实验或者理论计算等方法获取材料在不同应变条件下的应力应变关系数据。
(3)将实验或理论得到的应力应变数据输入到abaqus软件中,通过合适的参数设置和插值方法,建立起材料的应力应变曲线本构模型。
2.3 应力应变曲线本构的应用输入应力应变曲线本构后,abaqus可以进行模拟分析,对材料在不同工况下的力学性能进行计算和预测。
利用输入的应力应变曲线本构,abaqus可以准确模拟材料的强度、刚度、变形性能等,并为工程设计和分析提供重要的参考。
三、输入应力应变曲线本构的注意事项3.1 数据的准确性输入应力应变曲线本构前,需要确保所使用的应力应变数据具有较高的准确性和可靠性。
对于实验数据,需要进行充分的测试和验证;对于理论计算数据,需要保证所使用的材料模型和参数的准确性。
3.2 本构模型的选择在进行输入应力应变曲线本构时,需要根据材料的力学性能和所处的工程条件,选择合适的本构模型。
基于Abaqus子程序的高分子材料本构关系实现

基于 A b a q u s 子程 序 的高 分 子 材 料本 构 关 系实 现
汪 品红
abaqus 模拟材料力学扭转试验方法

abaqus 模拟材料力学扭转试验方法Abaqus模拟材料力学扭转试验方法引言:材料力学扭转试验是评估材料在受到扭转载荷下的变形和破坏特性的一种重要方法。
在工程实践中,了解材料的扭转性能对于设计和分析各种结构和部件具有重要意义。
本文将介绍如何使用ABAQUS软件进行材料力学扭转试验的模拟,以便更好地理解材料的扭转行为。
一、试验背景材料力学扭转试验是通过施加一个扭转载荷来研究材料的变形和破坏行为。
在试验中,材料样本通常采用圆柱形状,通过扭转机构施加一个扭转力矩,测量样本的扭转角度和扭转力矩,从而得到材料的扭转应力和扭转应变。
二、建立模型1. 创建几何模型:根据实际材料样本的几何形状,使用ABAQUS提供的几何建模工具创建一个圆柱形样本模型。
2. 定义材料属性:根据实际材料的力学性能指标,定义材料的弹性模量、泊松比和屈服强度等参数。
3. 网格划分:对样本进行网格划分,确保网格密度适中,以便准确模拟材料的变形和破坏行为。
三、施加边界条件在模拟中,需要施加适当的边界条件来模拟真实试验中的约束和加载条件。
根据实际情况,可以采用固定边界条件约束样本的底部,并施加一个扭转力矩在样本的顶部。
四、定义材料本构模型ABAQUS提供了多种本构模型来描述材料的力学行为。
根据实际材料的性质和试验要求,选择适当的本构模型,如线性弹性模型、非线性弹性模型或塑性本构模型。
五、指定加载条件根据试验要求,指定加载条件,即施加扭转力矩的大小和加载速率。
可以根据实际试验数据或已有的材料力学性能曲线来指定加载条件。
六、运行模拟并分析结果在模型设置完成后,运行ABAQUS模拟,并通过后处理功能来分析结果。
可以得到材料样本的扭转角度、扭转力矩和应力应变分布等信息,以评估材料的扭转性能。
七、结果验证与优化将模拟结果与实验数据进行对比,并进行结果验证与优化。
如果模拟结果与实验结果吻合较好,则可以认为模拟是准确可靠的。
如果存在差异,需要对模型和参数进行调整,以更好地模拟材料的扭转行为。
ABAQUS显式分析梁单元的混凝土、钢筋本构模型共3篇

ABAQUS显式分析梁单元的混凝土、钢筋本构模型共3篇ABAQUS显式分析梁单元的混凝土、钢筋本构模型1在ABAQUS中,梁单元是一种经常用于模拟混凝土和钢筋梁的元素。
它使用线性或非线性混凝土本构模型和钢筋本构模型来描述材料的行为,并考虑梁单元在三个方向上的应力和应变。
混凝土本构模型:ABAQUS提供了多个混凝土本构模型,它们可以用于描述混凝土的本构行为。
其中一个常用的模型是Mander本构模型,它考虑了混凝土的三个不同阶段的行为:1. 压缩阶段: 混凝土在受到压缩时会逐渐变硬,所以Mander模型使用一个非线性的应力-应变关系来描述混凝土的压缩行为。
该模型使用三个参数来描述混凝土在不同应变范围内的硬化行为。
2. 弯曲-拉伸阶段: 当混凝土受到弯曲或拉伸时,会发生一些微小的裂缝,导致其变得更容易受到破坏。
因此,Mander模型采用一个渐进应力-应变关系来描述混凝土的弯曲和拉伸行为。
该模型也使用三个参数来描述不同应变范围内的弯曲和拉伸行为。
3. 破坏阶段: 当混凝土受到极大应力时,会发生破坏。
为了模拟破坏行为,Mander模型使用两个参数来描述混凝土的弹性模量和极限应变。
当混凝土受到超过极限应变的应变时,该模型将输出一个非常大的应力值,这意味着梁单元已经破坏。
钢筋本构模型:ABAQUS也提供了多个钢筋本构模型。
其中一个常用的模型是多屈服弹塑性模型,它考虑了钢筋的应力-应变关系的多个拐点:1. 弹性阶段: 在应力小于屈服强度时,钢筋的行为是弹性的。
因此,多屈服弹塑性模型使用一个线性应力-应变关系来描述弹性阶段的行为。
2. 屈服阶段: 当钢筋的应力达到屈服强度时,它的行为将开始变得非线性。
因此,多屈服弹塑性模型使用一个拐点来描述屈服后的应力-应变关系。
该模型使用一组参数来描述每个拐点的应力和应变差。
3. 再次弹性阶段: 当钢筋的应变超过屈服点后,它的应变-应力关系将再次变得线性。
多屈服弹塑性模型也考虑了这个阶段的行为。
(完整版)ABAQUS中Cohesive单元建模方法讲解
复合材料模型建模与分析1. Cohesive单元建模方法1。
1 几何模型使用内聚力模型(cohesive zone)模拟裂纹的产生和扩展,需要在预计产生裂纹的区域加入cohesive 层。
建立cohesive层的方法主要有:方法一、建立完整的结构(如图1(a)所示),然后在上面切割出一个薄层来模拟cohesive单元,用这种方法建立的cohesive单元与其他单元公用节点,并以此传递力和位移.方法二、分别建立cohesive层和其他结构部件的实体模型,通过“tie”绑定约束,使得cohesive单元两侧的单元位移和应力协调,如图1(b)所示。
(a)cohesive单元与其他单元公用节点 (b)独立的网格通过“tie"绑定图1.建模方法上述两种方法都可以用来模拟复合材料的分层失效,第一种方法划分网格比较复杂;第二种方法赋材料属性简单,划分网格也方便,但是装配及“tie"很繁琐;因此在实际建模中我们应根据实际结构选取较简单的方法。
1.2 材料属性应用cohesive单元模拟复合材料失效,包括两种模型:一种是基于traction-separation描述;另一种是基于连续体描述。
其中基于traction-separation描述的方法应用更加广泛.而在基于traction—separation描述的方法中,最常用的本构模型为图2所示的双线性本构模型。
它给出了材料达到强度极限前的线弹性段和材料达到强度极限后的刚度线性降低软化阶段. 注意图中纵坐标为应力,而横坐标为位移,因此线弹性段的斜率代表的实际是cohesive单元的刚度。
曲线下的面积即为材料断裂时的能量释放率。
因此在定义cohesive的力学性能时,实际就是要确定上述本构模型的具体形状:包括刚度、极限强度、以及临界断裂能量释放率,或者最终失效时单元的位移。
常用的定义方法是给定上述参数中的前三项,也就确定了cohesive的本构模型。
ABAQUS材料用户子程序UMAT学习报告
具有友好的用户 界面和易用的操 作流程,方便用 户进行学习和使 用
UMAT子程序简介
UMAT子程序是 ABAQUS材料用户 自定义模块,允许 用户根据实际需求 编写材料本构模型。
UMAT子程序采用C 语言编写,用户需要 具备一定的编程基础。
UMAT子程序可以实 现多种材料本构模型 ,如弹性、塑性、蠕 变等。
UMAT子程序实现细节
编程语言和接口
A B A Q U S 材 料 用 户 子 程 序 U M AT 使 用Fortran语言编写
U M AT 子 程 序 中 可 以 定 义 材 料 属 性 、 本构关系等
添加标题
添加标题
添加标题
添加标题
U M AT 子 程 序 通 过 A B A Q U S 提 供 的 接口与主程序进行交互
不足:使用门槛较高,需要用户具备一定的编程基础
未来展望:期待更多的用户参与开发,不断完善子程序功能
总结:UMAT子程序为用户提供了强大的材料模型描述能力,但使用过程中需要注意其局 限性
在ABAQUS中的未来发展方向
开发更高效的材料模型 集成人工智能和机器学习技术 增强与CAD软件的集成 扩展对多物理场模拟的支持
适用于金属材料
适用于复合材料
适用于橡胶材料
适用于陶瓷材料
参数的合理选择
参数选择需符合实际物理模型 参数选择需考虑材料特性 参数选择需经过实验验证 参数选择需注意收敛性和稳定性
收敛性和稳定性问题
收 敛 性 : U M AT 子 程 序 在 迭 代 过 程 中应满足收敛条件,否则可能导致 计算失败或结果不准确。
边界条件和初始条件
边界条件:描述模型在边界上的行为,如位移、速度等 初始条件:描述模型在初始时刻的状态,如温度、压力等
abaqus中johnson-cook本构模型理解 -回复
abaqus中johnson-cook本构模型理解-回复什么是ABAQUS?ABAQUS是一种商业有限元分析软件,可用于模拟和分析复杂的结构、材料和物理过程。
它提供了广泛的分析功能和完整的建模工具,以帮助工程师和研究人员解决各种问题。
Johnson-Cook本构模型是什么?本构模型是材料力学中用于描述物质的应变-应力关系的数学模型。
Johnson-Cook本构模型是一种经验本构模型,适用于金属材料在高应变速率和高温条件下的塑性行为。
它是根据大量实验数据拟合而得,能够准确地预测材料在复杂工况下的力学行为。
Johnson-Cook本构模型的基本形式为:σ= (A + Bε^n)(1 + C(ln(ε_p))^m)(1 - (T - T_0)^p)其中,σ是应力,ε是总应变,ε_p是等效塑性应变,T是温度,T_0是参考温度。
A、B、C、n、m和p是拟合参数,可以根据材料的实验数据进行确定。
Johnson-Cook本构模型的理解和应用:1. 适用范围和限制Johnson-Cook本构模型适用于金属材料在高应变速率和高温条件下的塑性行为。
它广泛应用于冲击、爆炸、高速碰撞等动力学加载条件下的模拟分析。
然而,它也有一些局限性,如适用于单一相材料,对复杂的孪晶结构和多相材料的预测能力有限。
2. 模型参数的确定Johnson-Cook本构模型的参数包括A、B、C、n、m和p,这些参数可以通过材料的实验数据和拟合方法进行确定。
A和B参数分别表示材料的强度和硬化特性,C参数表示材料的温度敏感性,n参数表示材料的硬化指数,m参数表示材料的温度敏感指数,p参数表示材料的温度衰减指数。
通过对实验数据的统计分析和拟合,可以得到这些参数的值。
3. 模型的应用Johnson-Cook本构模型可以在ABAQUS软件中使用,用于模拟和分析金属材料在高应变速率和高温条件下的力学行为。
通过将该模型与适当的加载条件和边界条件结合使用,可以准确地预测材料在复杂工况下的力学响应。
ABAQUS子程序UMAT里弹塑本构的实现之欧阳科创编
前言有限元法是工程中广泛使用的一种数值计算方法。
它是力学、计算方法和计算机技术相结合的产物。
在工程应用中,有限元法比其它数值分析方法更流行的一个重要原因在于:相对与其它数值分析方法,有限元法对边界的模拟更灵活,近似程度更高。
所以,伴随着有限元理论以及计算机技术的发展,大有限元软件的应用证变得越来越普及。
ABAQUS软件一直以非线性有限元分析软件而闻名,这也是它和ANSYS,Nastran等软件的区别所在。
非线性有限元分析的用处越来越大,因为在所用材料非常复杂很多情况下,用线性分析来近似已不再有效。
比方说,一个复合材料就不能用传统的线性分析软件包进行分析。
任何与时间有关联,有较大位移量的情况都不能用线性分析法来处理。
多年前,虽然非线性分析能更适合、更准确的处理问题,但是由于当时计算设备的能力不够强大、非线性分析软件包线性分析功能不够健全,所以通常采用线性处理的方法。
这种情况已经得到了极大的改善,计算设备的能力变得更加强大、类似ABAQUS这样的产品功能日臻完善,应用日益广泛。
非线性有限元分析在各个制造行业得到了广泛应用,有不少大型用户。
航空航天业一直是非线性有限元分析的大客户,一个重要原因是大量使用复合材料。
新一代波音 787客机将全部采用复合材料。
只有像 ABAQUS这样的软件,才能分析包括多个子系统的产品耐久性能。
在汽车业,用线性有限元分析来做四轮耐久性分析不可能得到足够准确的结果。
分析汽车的整体和各个子系统的性能要求(如悬挂系统等)需要进行非线性分析。
在土木工程业, ABAQUS能处理包括混凝土静动力开裂分析以及沥青混凝土方面的静动力分析,还能处理高度复杂非线性材料的损伤和断裂问题,这对于大型桥梁结构,高层建筑的结构分析非常有效。
瞬态、大变形、高级材料的碰撞问题必须用非线性有限元分析来计算。
线性分析在这种情况下是不适用的。
以往有一些专门的软件来分析碰撞问题,但现在ABAQUS在通用有限元软件包就能解决这些问题。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
材料本构模型及编程-ABAQUS-UMAT材料本构模型及编程实现:简介1、什么时候用用户定义材料(User-defined material, UMAT)?很简单,当ABAQUS没有提供我们需要的材料模型时。
所以,在决定自己定义一种新的材料模型之前,最好对ABAQUS已经提供的模型心中有数,并且尽量使用现有的模型,因为这些模型已经经过详细的验证,并被广泛接受。
2、好学吗?需要哪些基础知识?先看一下ABAQUS手册(ABAQUS Analysis User's Manual)里的一段话:Warning: The use of this option generally requires considerable expertise. The user is cautioned that the imple mentation of any realistic constitutive model requires extensive development and testing. Initial testing on a s ingle element model with prescribed traction loading is strongly recommended.但这并不意味着非力学专业,或者力学基础知识不很丰富者就只能望洋兴叹,因为我们的任务不是开发一套完整的有限元软件,而只是提供一个描述材料力学性能的本构方程(Constitutive equation)而已。
当然,最基本的一些概念和知识还是要具备的,比如应力(stress),应变(strain)及其分量;volumetric part和deviatoric part;模量(modulus)、泊松比(Poisson’s ratio)、拉美常数(Lame constant);矩阵的加减乘除甚至求逆;还有一些高等数学知识如积分、微分等。
3、UMAT的基本任务?我们知道,有限元计算(增量方法)的基本问题是:已知第n步的结果(应力,应变等),;然后给出一个应变增量, 计算新的应力。
UMAT要完成这一计算,并要计算Jacobian矩阵DDSDDE(I,J) =。
是应力增量矩阵(张量或许更合适),是应变增量矩阵。
DDSDDE(I,J) 定义了第J个应变分量的微小变化对第I 个应力分量带来的变化。
该矩阵只影响收敛速度,不影响计算结果的准确性(当然,不收敛自然得不到结果)。
4、怎样建立自己的材料模型?本构方程就是描述材料应力应变(增量)关系的数学公式,不是凭空想象出来的,而是根据实验结果作出的合理归纳。
比如对弹性材料,实验发现应力和应变同步线性增长,所以用一个简单的数学公式描述。
为了解释弹塑性材料的实验现象,又提出了一些弹塑性模型,并用数学公式表示出来。
对各向同性材料(Isotropic material),经常采用的办法是先研究材料单向应力-应变规律(如单向拉伸、压缩试验),并用一数学公式加以描述,然后把讲该规律推广到各应力分量。
这叫做“泛化“(generalization)。
5、一个完整的例子及解释下面这个UMAT取自ABAQUS手册,是一个用于大变形下的弹塑性材料模型。
希望我的注释能帮助初学者理解。
需要了解J2理论。
SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,RPL,DDSDDT,1 DRPLDE,DRPLDT,STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,2 CMNAME,NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,3 PNEWDT,CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC)STRESS--应力矩阵,在增量步的开始,保存并作为已知量传入UMAT ;在增量步的结束应该保存更新的应力;STRAN--当前应变,已知。
DSTRAN—应变增量,已知。
STATEV--状态变量矩阵,用来保存用户自己定义的一些变量,如累计塑性应变,粘弹性应变等等。
增量步开始时作为已知量传入,增量步结束应该更新;DDSDDE=。
需要更新DTIME—时间增量dt。
已知。
NDI—正应力、应变个数,对三维问题、轴对称问题自然是3(11,22,33),平面问题是2(11,22);已知。
NSHR —剪应力、应变个数,三维问题时3(12,13,23),轴对称问题是1(12);已知。
NTENS=NTENS NSHR,已知。
PROPS材料常数矩阵,如模量啊,粘度系数啊等等;作为已知量传入,已知。
DROT—对finite strain问题,应变应该排除旋转部分,该矩阵提供了旋转矩阵,详见下面的解释。
已知。
PNEWDT—可用来控制时间步的变化。
如果设置为小于1的数,则程序放弃当前计算,并用新的时间增量DTIME X PNEWDT作为新的时间增量计算;这对时间相关的材料如聚合物等有用;如果设为大余1的数,则下一个增量步加大DTIME为DTIME X PNEWDT。
可以更新。
其他变量含义可参看手册,暂时用不到。
CINCLUDE 'ABA_PARAM.INC'定义了一些参数,变量什么的,不用管CCHARACTER*8 CMNAMECDIMENSION STRESS(NTENS),STATEV(NSTATV),DDSDDE(NTENS,NTENS),1 DDSDDT(NTENS),DRPLDE(NTENS),STRAN(NTENS),DSTRAN(NTENS),2 PREDEF(1),DPRED(1),PROPS(NPROPS),COORDS(3),DROT(3,3),3 DFGRD0(3,3),DFGRD1(3,3)矩阵的尺寸声明CC LOCAL ARRAYSC ----------------------------------------------------------------C EELAS - ELASTIC STRAINSC EPLAS - PLASTIC STRAINSC FLOW - DIRECTION OF PLASTIC FLOWC ----------------------------------------------------------------C局部变量,用来暂时保存弹性应变、塑性应变分量以及流动方向DIMENSION EELAS(6),EPLAS(6),FLOW(6)CPARAMETER(ZERO=0.D0,ONE=1.D0,TWO=2.D0,THREE=3.D0,SIX=6.D0,1 ENUMAX=.4999D0,NEWTON=10,TOLER=1.0D-6)CC ----------------------------------------------------------------C UMAT FOR ISOTROPIC ELASTICITY AND ISOTROPIC MISES PLASTICITYC CANNOT BE USED FOR PLANE STRESSC ----------------------------------------------------------------C PROPS(1) - EC PROPS(2) - NUC PROPS(3..) - SYIELD AN HARDENING DATAC CALLS HARDSUB FOR CURVE OF YIELD STRESS VS. PLASTIC STRAINC ----------------------------------------------------------------CC ELASTIC PROPERTIESC获取杨氏模量,泊松比,作为已知量由PROPS向量传入EMOD=PROPS(1) EENU=PROPS(2) νEBULK3=EMOD/(ONE-TWO*ENU) 3KEG2=EMOD/(ONE ENU) 2GEG=EG2/TWO GEG3=THREE*EG 3GELAM=(EBULK3-EG2)/THREE λDO K1=1,NTENSDO K2=1,NTENSDDSDDE(K1,K2)=ZEROEND DOEND DO弹性部分,Jacobian矩阵很容易计算注意,在ABAQUS中,剪切应变采用工程剪切应变的定义,所以剪切部分模量是G而不是2G!CC ELASTIC STIFFNESSCDO K1=1,NDIDO K2=1,NDIDDSDDE(K2,K1)=ELAMEND DODDSDDE(K1,K1)=EG2 ELAMEND DODO K1=NDI 1,NTENSDDSDDE(K1,K1)=EGEND DOCC RECOVER ELASTIC AND PLASTIC STRAINS AND ROTATE FORWARDC ALSO RECOVER EQUIVALENT PLASTIC STRAINC读取弹性应变分量,塑性应变分量,并旋转(调用了ROTSIG),分别保存在EELAS和EPLAS中;CALL ROTSIG(STATEV( 1),DROT,EELAS,2,NDI,NSHR)CALL ROTSIG(STATEV(NTENS 1),DROT,EPLAS,2,NDI,NSHR)读取等效塑性应变EQPLAS=STATEV(1 2*NTENS)先假设没有发生塑性流动,按完全弹性变形计算试算应力CC CALCULATE PREDICTOR STRESS AND ELASTIC STRAINCDO K1=1,NTENSDO K2=1,NTENSSTRESS(K2)=STRESS(K2) DDSDDE(K2,K1)*DSTRAN(K1)END DOEELAS(K1)=EELAS(K1) DSTRAN(K1)END DOC计算Mises应力C CALCULATE EQUIVALENT VON MISES STRESSCSMISES=(STRESS(1)-STRESS(2))**2 (STRESS(2)-STRESS(3))**21 (STRESS(3)-STRESS(1))**2DO K1=NDI 1,NTENSSMISES=SMISES SIX*STRESS(K1)**2END DOSMISES=SQRT(SMISES/TWO)C 根据当前等效塑性应变,调用HARDSUB得到对应的屈服应力C GET YIELD STRESS FROM THE SPECIFIED HARDENING CURVE CNVALUE=NPROPS/2-1CALL HARDSUB(SYIEL0,HARD,EQPLAS,PROPS(3),NVALUE)CC DETERMINE IF ACTIVELY YIELDINGC 如果Mises应力大余屈服应力,屈服发生,计算流动方向IF (SMISES.GT.(ONE TOLER)*SYIEL0) THENCC ACTIVELY YIELDINGC SEPARATE THE HYDROSTATIC FROM THE DEVIATORIC STRESS C CALCULATE THE FLOW DIRECTIONCSHYDRO=(STRESS(1) STRESS(2) STRESS(3))/THREEDO K1=1,NDIFLOW(K1)=(STRESS(K1)-SHYDRO)/SMISESEND DODO K1=NDI 1,NTENSFLOW(K1)=STRESS(K1)/SMISESEND DOC根据J2理论并应用Newton-Rampson方法求得等效塑性应变增量C SOLVE FOR EQUIVALENT VON MISES STRESSC AND EQUIVALENT PLASTIC STRAIN INCREMENT USING NEWTON ITERATION CSYIELD=SYIEL0DEQPL=ZERODO KEWTON=1,NEWTONRHS=SMISES-EG3*DEQPL-SYIELDDEQPL=DEQPL RHS/(EG3 HARD)CALL HARDSUB(SYIELD,HARD,EQPLAS DEQPL,PROPS(3),NVALUE)IF(ABS(RHS).LT.TOLER*SYIEL0) GOTO 10END DOCC WRITE WARNING MESSAGE TO THE .MSG FILECWRITE(7,2) NEWTON2 FORMAT(//,30X,'***WARNING - PLASTICITY ALGORITHM DID NOT ',1 'CONVERGE AFTER ',I3,' ITERATIONS')10 CONTINUEC更新应力,应变分量C UPDATE STRESS, ELASTIC AND PLASTIC STRAINS ANDC EQUIVALENT PLASTIC STRAINCDO K1=1,NDISTRESS(K1)=FLOW(K1)*SYIELD SHYDROEPLAS(K1)=EPLAS(K1) THREE/TWO*FLOW(K1)*DEQPLEELAS(K1)=EELAS(K1)-THREE/TWO*FLOW(K1)*DEQPLEND DODO K1=NDI 1,NTENSSTRESS(K1)=FLOW(K1)*SYIELDEPLAS(K1)=EPLAS(K1) THREE*FLOW(K1)*DEQPLEELAS(K1)=EELAS(K1)-THREE*FLOW(K1)*DEQPLEND DOEQPLAS=EQPLAS DEQPLCC CALCULATE PLASTIC DISSIPATIONCSPD=DEQPL*(SYIEL0 SYIELD)/TWOCC 计算塑性变形下的Jacobian矩阵FORMULATE THE JACOBIAN (MATERIAL TANGENT)C FIRST CALCULATE EFFECTIVE MODULICEFFG=EG*SYIELD/SMISESEFFG2=TWO*EFFGEFFG3=THREE/TWO*EFFG2EFFLAM=(EBULK3-EFFG2)/THREEEFFHRD=EG3*HARD/(EG3 HARD)-EFFG3c...if (props(7).lt..001) go to 99c...DO K1=1,NDIDO K2=1,NDIDDSDDE(K2,K1)=EFFLAMEND DODDSDDE(K1,K1)=EFFG2 EFFLAMEND DODO K1=NDI 1,NTENSDDSDDE(K1,K1)=EFFGEND DODO K1=1,NTENSDO K2=1,NTENSDDSDDE(K2,K1)=DDSDDE(K2,K1) EFFHRD*FLOW(K2)*FLOW(K1)END DOEND DOc...99 continuec...ENDIFC将弹性应变,塑性应变分量保存到状态变量中,并传到下一个增量步C STORE ELASTIC AND (EQUIVALENT) PLASTIC STRAINSC IN STATE VARIABLE ARRAYCDO K1=1,NTENSSTATEV(K1)=EELAS(K1)STATEV(K1 NTENS)=EPLAS(K1)END DOSTATEV(1 2*NTENS)=EQPLASCRETURNENDc...c...子程序,根据等效塑性应变,利用插值的方法得到对应的屈服应力SUBROUTINE HARDSUB(SYIELD,HARD,EQPLAS,TABLE,NVALUE)CINCLUDE 'ABA_PARAM.INC'CDIMENSION TABLE(2,NVALUE)CPARAMETER(ZERO=0.D0)CC SET YIELD STRESS TO LAST VALUE OF TABLE, HARDENING TO ZERO CSYIELD=TABLE(1,NVALUE)HARD=ZEROC IF MORE THAN ONE ENTRY, SEARCH TABLECIF(NVALUE.GT.1) THENDO K1=1,NVALUE-1EQPL1=TABLE(2,K1 1)IF(EQPLAS.LT.EQPL1) THENEQPL0=TABLE(2,K1)IF(EQPL1.LE.EQPL0) THENWRITE(7,1)1 FORMAT(//,30X,'***ERROR - PLASTIC STRAIN MUST BE `,1 `ENTERED IN ASCENDING ORDER')CALL XITENDIFCC CURRENT YIELD STRESS AND HARDENINGCDEQPL=EQPL1-EQPL0SYIEL0=TABLE(1,K1)SYIEL1=TABLE(1,K1 1)DSYIEL=SYIEL1-SYIEL0HARD=DSYIEL/DEQPLSYIELD=SYIEL0 (EQPLAS-EQPL0)*HARD GOTO 10ENDIFEND DO10 CONTINUEENDIFRETURNEND。