ABAQUS-材料本构模型及编程
- 格式:docx
- 大小:20.90 KB
- 文档页数:8
材料本构模型及编程-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 implementation of any realis tic constitutive model requires extensive development and testing. Initial testing on a single element model with prescribed traction l oading is strongly recommended.但这并不意味着非力学专业,或者力学基础知识不很丰富者就只能望洋兴叹,因为我们的任务不是开发一套完整的有限元软件,而只是提供一个描述材料力学性能的本构方程(Constitutive equation)而已。
当然,最基本的一些概念和知识还是要具备的,比如应力(stress),应变(strain)及其分量; volumetric part和deviatoric part;模量(modulus)、泊松比(Poisson’s ratio)、拉美常数(Lame constant);矩阵的加减乘除甚至求逆;还有一些高等数学知识如积分、微分等。
abaquscdp本构原理
ABAQUS的CDP(Concrete Damaged Plasticity)模型是一种混凝土本
构关系模型,用于描述混凝土的非弹性行为。
该模型通过将各向同性下损伤弹性与拉伸和压缩塑性相结合的方式来描述混凝土的非弹性行为,适用于模拟混凝土在任意荷载作用下的受力情况。
CDP模型考虑了由于拉、压塑性
应变导致的弹性刚度的退化以及循环荷载作用下刚度的恢复,具有较好的收敛性。
CDP模型采用混凝土在单轴受力状态下的应力和非弹性应变,这里的非弹
性应变是根据混凝土的单轴应力-应变关系(混凝土本构关系)换算出来的。
混凝土本构关系有3种:GB《混凝土结构设计规范》欧洲规范、Kent-Park 模型。
CDP模型中,混凝土材料的弹性模量E c 可通过结构试验进行实测,也可以查表,也可以根据下式进行计算:E c = 10^5 × + ( / f cu , k)。
其中,fcu,k为混凝土的峰值抗压强度。
此外,CDP模型本构曲线末尾段的选取,对滞回曲线下降段的影响较大。
为了验证所编子程序的合理性与正确性,可以选用具体的有限元模型进行验证。
以上内容仅供参考,如需更多信息,建议查阅ABAQUS软件相关书籍或咨询软件专家。
一张图掌握Abaqus复合材料层板结构基础建模对于Abaqus复合材料分析初学者,会遇到各类奇怪的错误,其中最常见的一类问题就是由于材料模型、截面属性、网格类型等设置不匹配造成的错误,尤其是显式分析中使用实体单元时,还需要特殊对截面属性进行设置。
前面的文章里介绍了Abaqus复合材料层压板结构的主要建模方法有经典建模方法和layup快速建模方法,本文将针对经典建模方法,用一张图介绍一下几何类型、材料模型、截面属性、网格类型以及适用的求解器等多个因素的匹配关系。
复合材料经典建模方法建模顺序先回顾一下Abaqus复合材料层压结构的经典建模方法建模顺序,如下:Part→Material→Section→Assignsection→Assign material orientation →Createmesh→Assign mesh type→Assembly→Step→output→Interaction → Load → Job→ Visualization红色标注部分为与金属等各向同性材料建模有区别的步骤。
网格划分步骤建议提前到Part之后,避免因几何切分造成的几何信息丢失或错乱等错误。
匹配图几点说明▪几何类型:指part 模块下建立的几何形状,面指平面或曲面,无厚度。
体指三维实体几何。
▪本表格仅列举了2D lamina及3D Engineering constants两类材料模型,可基本满足绝大部分层压板建模。
▪不排除表格以外的其他建模方法,但是按照本表格中的匹配关系建模可保证无误。
▪特别需要指出的是,许多人误以为显式分析中无法使用三维材料模型及三维实体单元,其实是完全可以的,只是建模方法和普通的壳单元/连续壳元有较大差异。
显式求解其中使用Engineering constants这种材料模型时,截面属性要用Solid homogeneous这种截面属性,一层层切出来,逐层赋属性。
基于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中johnson-cook本构模型理解-回复什么是Johnson-Cook本构模型?Johnson-Cook本构模型是一种常用的冲击/爆炸动力学材料本构模型,用于描述高速冲击、爆炸等极端应变速率条件下材料的动态力学行为。
它是由Johnson、Cook等人在1983年基于强化流动材料在宇宙空间中引起的一系列实验数据提出的。
该模型可用于预测金属材料在高温、高速率变形下的应力-应变关系,以及材料的生命损伤和破坏。
Johnson-Cook本构模型的原理是将应力-应变关系分解成四个部分:弹性变形、塑性变形、强化流动和软化流动。
其中,弹性变形是指材料在受力后能恢复到原始形状的变形行为;塑性变形是指材料在超过屈服点后出现的不可逆塑性变形行为;强化流动是指材料在高速率下发生的塑性变形行为,它是由高速率塑性变形引起的材料强度增加;软化流动是指材料在高温下、高应变速率下出现的塑性变形行为,它是由高温下的材料软化引起的材料强度降低。
通过将以上四个部分结合起来,Johnson-Cook本构模型可以较为准确地描述金属材料在极端应变速率条件下的应力-应变关系。
具体来说,Johnson-Cook本构模型的形式为:σ= ε^(n) * (A + B * ε^(m)) * (1 + (C * ln(ε^(˙)/ε̇_0 ))^(α-1) )其中,σ为材料应力,ε为材料应变,ε˙为应变速率,ε̇_0为参考应变速率,A、B、C、n、m和α为经验参数。
其中A和B可表示材料的初始硬化性能,C可表示高速率下流动机制,n可表示材料的应变硬化效应,m可表示材料的应变硬化率,α可表示材料的软化速度。
Johnson-Cook本构模型的应用范围非常广泛,特别适合用于金属材料在高速率和高温环境下的动态加载计算。
该模型的参数可以通过实验数据拟合获得,使得模型的预测结果能够与实际测试结果较好地吻合。
同时,Johnson-Cook本构模型还可以用于模拟金属材料的损伤与破坏过程,包括疲劳寿命预测、裂纹扩展等。
ABAQUS显式分析梁单元的混凝土、钢筋本构模型共3篇ABAQUS显式分析梁单元的混凝土、钢筋本构模型1在ABAQUS中,梁单元是一种经常用于模拟混凝土和钢筋梁的元素。
它使用线性或非线性混凝土本构模型和钢筋本构模型来描述材料的行为,并考虑梁单元在三个方向上的应力和应变。
混凝土本构模型:ABAQUS提供了多个混凝土本构模型,它们可以用于描述混凝土的本构行为。
其中一个常用的模型是Mander本构模型,它考虑了混凝土的三个不同阶段的行为:1. 压缩阶段: 混凝土在受到压缩时会逐渐变硬,所以Mander模型使用一个非线性的应力-应变关系来描述混凝土的压缩行为。
该模型使用三个参数来描述混凝土在不同应变范围内的硬化行为。
2. 弯曲-拉伸阶段: 当混凝土受到弯曲或拉伸时,会发生一些微小的裂缝,导致其变得更容易受到破坏。
因此,Mander模型采用一个渐进应力-应变关系来描述混凝土的弯曲和拉伸行为。
该模型也使用三个参数来描述不同应变范围内的弯曲和拉伸行为。
3. 破坏阶段: 当混凝土受到极大应力时,会发生破坏。
为了模拟破坏行为,Mander模型使用两个参数来描述混凝土的弹性模量和极限应变。
当混凝土受到超过极限应变的应变时,该模型将输出一个非常大的应力值,这意味着梁单元已经破坏。
钢筋本构模型:ABAQUS也提供了多个钢筋本构模型。
其中一个常用的模型是多屈服弹塑性模型,它考虑了钢筋的应力-应变关系的多个拐点:1. 弹性阶段: 在应力小于屈服强度时,钢筋的行为是弹性的。
因此,多屈服弹塑性模型使用一个线性应力-应变关系来描述弹性阶段的行为。
2. 屈服阶段: 当钢筋的应力达到屈服强度时,它的行为将开始变得非线性。
因此,多屈服弹塑性模型使用一个拐点来描述屈服后的应力-应变关系。
该模型使用一组参数来描述每个拐点的应力和应变差。
3. 再次弹性阶段: 当钢筋的应变超过屈服点后,它的应变-应力关系将再次变得线性。
多屈服弹塑性模型也考虑了这个阶段的行为。
Abaqus混凝土材料塑性损伤模型浅析与参数设置【壹讲壹插件】欢迎转载,作者:星辰-北极星,QQ群:431603427Abaqus混凝土材料塑性损伤模型浅析与参数设置 (1)第一部分:Abaqus自带混凝土材料的塑性损伤模型 (2)1.1概要 (2)1.2学习笔记 (2)1.3 参数定义与说明 (3)1.3.1材料模型选择:Concrete Damaged Plasticity (3)1.3.2 混凝土塑性参数定义 (3)1.3.3 混凝土损伤参数定义: (4)1.3.4 损伤参数定义与输出损伤之间的关系 (4)1.3.5 输出参数: (4)第二部分:根据GB50010-2010定义材料损伤值 (5)第三部分:星辰-北极星插件介绍:POLARIS-CONCRETE (6)3.1 概要 (6)3.2 插件的主要功能 (6)3.3 插件使用方法: (6)3.3.1 插件界面: (6)3.3.2 生成结果 (7)3.4、算例: (9)3.4.1三维实体简支梁模型说明 (9)3.4.2 计算结果: (9)第一部分:Abaqus自带混凝土材料的塑性损伤模型1.1概要首先我要了解Abaqus内自带的参数模型是怎样的,了解其塑性模型,进而了解其损伤模型,其帮助文档Abaqus Theory Manual 4.5.1 An inelastic constitutive model for concrete讲述的是其非弹性本构,4.5.2 Damaged plasticity model for concrete and other quasi-brittle materials则讲述的塑性损伤模型,同时在Abaqus Analysis User's Manual 22.6 Concrete也讲述了相应的内容。
1.2学习笔记1、混凝土塑性损伤本构模型中的损伤是一标量值,数值范围为(0无损伤~1完全失效[对于混凝土塑性损伤一般不存在]);2、仅适用于脆性材料在中等围压条件(为围压小于轴抗压强度1/4);3、拉压强度可设置成不同数值;4、可实现交变载荷下的刚度恢复;默认条件下,由拉转压刚度恢复,由压转拉刚度不变;5、强度与应变率相关;6、使用的是非相关联流动法则,刚度矩阵为非对称,因此在隐式分析步设置时,需在分析定义other-》Matrix storate-》Unsymmetric。
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混凝土材料塑性损伤模型浅析与参数设置【壹讲壹插件】欢迎转载,作者:星辰-北极星,QQ群:431603427Abaqus混凝土材料塑性损伤模型浅析与参数设置 (1)第一部分:Abaqus自带混凝土材料的塑性损伤模型 (2)1.1概要 (2)1.2学习笔记 (2)1.3 参数定义与说明 (3)1.3.1材料模型选择:Concrete Damaged Plasticity (3)1.3.2 混凝土塑性参数定义 (3)1.3.3 混凝土损伤参数定义: (4)1.3.4 损伤参数定义与输出损伤之间的关系 (4)1.3.5 输出参数: (4)第二部分:根据GB50010-2010定义材料损伤值 (5)第三部分:星辰-北极星插件介绍:POLARIS-CONCRETE (6)3.1 概要 (6)3.2 插件的主要功能 (6)3.3 插件使用方法: (6)3.3.1 插件界面: (6)3.3.2 生成结果 (7)3.4、算例: (9)3.4.1三维实体简支梁模型说明 (9)3.4.2 计算结果: (9)第一部分:Abaqus自带混凝土材料的塑性损伤模型1.1概要首先我要了解Abaqus内自带的参数模型是怎样的,了解其塑性模型,进而了解其损伤模型,其帮助文档Abaqus Theory Manual 4.5.1 An inelastic constitutive model for concrete讲述的是其非弹性本构,4.5.2 Damaged plasticity model for concrete and other quasi-brittle materials则讲述的塑性损伤模型,同时在Abaqus Analysis User's Manual 22.6 Concrete也讲述了相应的内容。
1.2学习笔记1、混凝土塑性损伤本构模型中的损伤是一标量值,数值范围为(0无损伤~1完全失效[对于混凝土塑性损伤一般不存在]);2、仅适用于脆性材料在中等围压条件(为围压小于轴抗压强度1/4);3、拉压强度可设置成不同数值;4、可实现交变载荷下的刚度恢复;默认条件下,由拉转压刚度恢复,由压转拉刚度不变;5、强度与应变率相关;6、使用的是非相关联流动法则,刚度矩阵为非对称,因此在隐式分析步设置时,需在分析定义other-》Matrix storate-》Unsymmetric。
材料本构模型及编程-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 implementation of any realis tic constitutive model requires extensive development and testing. Initial testing on a single element model with prescribed traction l oading 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。