ABAQUS_子程序VUMAT_的坐标系的问题的讨论

SimWe仿真论坛---(论坛实行邀请码注册)'s Archiver SimWe仿真论坛---(论坛实行邀请码注册) » A04：ABAQUS » ABAQUS 子程序VUMAT 的坐标系的问题的讨论 billowriver 发表于 2009-4-6 05:50 ABAQUS 子程序VUMAT 的坐标系的问题的讨论 在ABAQUS/Explicit里面，Vumat里的所有变量都是相对于一个corotational 坐标系 来说的，而在ABAQUS的odb 里面显示的是cauchy stress, 是相对于当前坐标系来讲的。 因为变形梯度 F=RU，意思是材料先变形U，再做一个刚体旋转R，而这个corotational坐标系是在材料做刚体旋转 前的坐标系，ABAQUS上面说由于做刚体旋转不会影响应力在相应坐标系里的大小，因此刚体旋转前后应力大小应该 一样，即VUMAT里面的StressNew和odb里面显示的cauchy stress大小应该一样。 我在做单轴拉伸模拟得到的结果二者一样，但是如果做剪切模拟，如果剪切力和正应力同时存在的话，StressNew 和odb里面的应力不一样，即使把StressNew设成常数，在odb里得到的应力也不是常数。很奇怪。 这种影响在小变形下差别很小，但是大变形时候会产生很大误差。如何用VUMAT计算出当前的cauchy stress, 哪位 高手知道或者有兴趣的可以一起讨论下。 S_true=RSR'，但是ABAQUS里面没有直接提供R，而是提供了F_old，F_new, U_old, U_new，但是Vumat里 面应变也是基于corotational坐标系的，而且是lnV,在这里=lnU shawn2008 发表于 2009-4-6 06:05 [i=s] 本帖最后由 shawn2008 于 2009-4-6 06:07 编辑 [/i] 太多概念错误，我都不知道从哪开始说起了。 billowriver 发表于 2009-4-9 01:02 请多指教，这个都把我搞糊涂了 billowriver 发表于 2009-4-9 01:09 上面我讲的都是出自ABAQUS document，而且从连续介质力学角度讲也应该是这样，只是ABAQUS处理VUMAT 时多了个corotational坐标系，得到的结果有些问题 不知道哪写概念有问题，还请shawn2008 说说 shawn2008 发表于 2009-4-9 01:40 呵呵，别误会，我以前也是像你写的这样理解的；document的东西是有前提的。 ok, 我们先说第一句的前半句话，’在ABAQUS/Explicit里面，Vumat里的所有变量都是相对于一个corotational 坐 标系 来说的‘。 这里有几个概念请你要先说明一下， ‘所有变量’指的是什么变量？ corotational坐标系是怎么定义的？ 相对于是指 的什么含义？ billowriver 发表于 2009-4-9 02:29 [i=s] 本帖最后由 billowriver 于 2009-4-9 02:43 编辑 [/i] 谢谢shawn2008的答复，所有变量就是应力，应变，应变率这些，当然ABAQUS给出的变形梯度F应该还是度量 current坐标系和initial坐标系的。 corotational坐标系是材料只做变形，在发生刚体转动前的坐标系，这样我们就有了三个坐标系，initial coordinates, corotational coordinates, current coordinates， 有限元 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 核心增值经销商 北京盛维提供3D CAD/CAE/CAM/PDM/CA... 整 体解决方案。010-68250396/68253809 “相对于corotational坐标系”, 位移u_corotational=x_corotational-X_initial, 定义应变要用u_corotational 和 x_corotational，定义应力时，微元面积和垂直方向也要用x_corotational，此时的变形梯度其实就是U shawn2008 发表于 2009-4-9 09:51 [i=s] 本帖最后由 shawn2008 于 2009-4-9 09:54 编辑 [/i] 好的，所有的变量可以归纳为材料相关的张量，比如应力，应变，应力率，应变率，变形梯度F，速度梯度L，变形率 张量D等等。 张量的表达是先选择一个坐标系，或一组基，然后就能写出张量的各个分量。 现在先不讨论umat中的这些量，abaqus中所有的这些材料相关张量的定义都是在一个local system下的。 比如你 在后处理中见到的应力和umat中的变形梯度F的11， 22，33，12，23，31分量的值都是在local system下给出 的数值。abaqus中使用了3种方式定义这个local system，这个document中有详细介绍。 这里要区别张量的表达，张量在初始构型和当前构型中定义，张量的共旋部分这些概念。 我们现在只来讨论张量的表达，现在抛开abaqus，可以这样总结一下，要描述一个材料的本构行为，需要先选择一 个全局的坐标系（这个是最基础的坐标系，剩下的其它坐标系的基都是表达在这个最基础坐标系下的），有了这个全 局坐标系就可以给出应力，应变，变形梯度，速度梯度这些张量的分量了。 张量分量的脚标和张量分量的值是相对 于你选择坐标系的。 现在回到你的 ‘当然ABAQUS给出的变形梯度F应该还是度量current坐标系和initial坐标系的’， 变形梯度是来描述 当前构型和初始构型的形状变化，但它的分量的值是在全局坐标系下给出的（现这样理解）。 在谈张量的表达时，我并没有区分，如Cauchy 应力, 第一 Piola-Kirchhoff应力，第二Piola-Kirchhoff应力， Kirchhoff应力，还有Green应变的种类，因为这些概念的区别不是由你选择的坐标系决定的。 今天就先写这么多，明天还要上班了。 总结一句话，要描述本构特性，需要先选择一个全局坐标系，材料相关张量的表达就是在这个坐标系下给出的。 billowriver 发表于 2009-4-9 21:57 [i=s] 本帖最后由 billowriver 于 2009-4-9 22:22 编辑 [/i] 谢谢shawn2008说了这么多，这些都是连续介质力学里的概念，我也知道。如果没有这个corotational坐标系，没 有任何问题。如你所说，无论是度量哪个坐标系里的物理量都要跟未变形时候的初始坐标系比较，就是你说的全局坐 标系。 我的问题是这个corotational到底是哪个坐标系，我最开始以为它就是你说的local坐标系，但是ABAQUS里面有句 话说当考虑hyperelastic的时候，不能使用F，而要使用U，因为使用F的话你还需要再把应力转动回来。 还有，你说的Cauchy和1st,2nd PK应力跟坐标系选取没有关系应该是对固定时间来讲的，也就是说对于当前坐标系 x(X,t)，时间确定了，这些物理量也就确定了，但是如果你把整个材料变形过程分成两部分，第一步先做纯变形，x1 （X,t1)然后再做刚体转动x2=(x1,t2),这样得到的坐标系x1和x2是对不同时间点来讲的，因此如果你在x1坐标系下 计算得到的应力应变是不同的。 可以肯定的一点是，ABAQUS后处理输出的是Cauchy应力，他与Vumat里面得到的应力是不一样的。 Vumat和Umat的坐标系使用应该还是有点不一样的，ABAQUS认为Vumat得到的应力不需要做任何转动，但是 Umat有时候需要 shawn2008 发表于 2009-4-10 01:23 [i=s] 本帖最后由 shawn2008 于 2009-4-10 01:33 编辑 [/i] 1. corotational坐标系是与某个旋率相对应的坐标系，该坐标系剔除了变形过程的旋转影响，是满足本构模型的客观 性要求的坐标系。 可能有的人对corotational坐标系的定义与我的理解不一样，但我习惯称这样的坐标系为 corotational坐标系。 我这样回答，可能你还是没有明白什么是corotational坐标系。 2. 'Cauchy和1st,2nd PK应力跟坐标系选取没有关系应该是对固定时间来讲的，也就是说对于当前坐标系x(X,t)， 时间确定了，这些物理量也就确定了，但是如果你把整个材料变形过程分成两部分，第一步先做纯变形，x1（X,t1) 然后再做刚体转动x2=(x1,t2),这样得到的坐标系x1和x2是对不同时间点来讲的，因此如果你在x1坐标系下计算得 到的应力应变是不同的。可以肯定的一点是，ABAQUS后处理输出的是Cauchy应力，他与Vumat里面得到的应力是 不一样的。 这些张量的分量的值当然与你选择的坐标系有关，但当你选定了一个坐标系，怎么样区分这些应力？也就是它们的定 义不是根据坐标系来确定的。还有这里要区分configuration 和坐标系的概念，这两个概念对理解这些应力的定义非 常的重要，你的理解是什么？ 3. Vumat和Umat的坐标系使用应该还是有点不一样的，ABAQUS认为Vumat得到的应力不需要做任何转动，但是 Umat有时候需要. 这两个里面的坐标系是有区别的，一个使用的是Green-Naghdi corotational rate一个使用的是Jaumann corotational rate。 4. 你经常理解的先变形后旋转，或先旋转后变形，这当然是可能的变形路径。但在实际上的成形过程中，通常都是边 旋转边变形，我们只能观察到物体产生了怎样的形状变化，里面有多少是与本构无关的旋转，有多少是与本构有关的 变形我们是不能观察到的。正因如此，不同的人提出了不同的理论来定义旋转，所以就出现了不同的corotational rate。 billowriver 发表于 2009-4-10 05:06 1，你说的应该是对的， ABAQUS在介绍Vumat时有这么一句话说： Hyperelastic constitutive models in VUMAT should be defined in a corotational coordinate system in which the basis system rotates with the material. This is most effectively accomplished by formulating the hyperelastic constitutive model in terms of the stretch tensor, U, instead of in terms of the deformation gradient,F . Using the deformation gradient can present some difficulties because the deformation gradient includes the rotation tensor and the resulting stresses would need to be rotated back to the corotational basis 我得理解是这样的，其实变形梯度的分解无碍乎两种F=RU和F=VR，这里R是unitary matrix，就是R*R'=I，只是 让矩阵做了个纯旋转，U和V是right和left Green Deformation tensor，是纯变形对称阵，当然材料是一边变形一 边旋转的，但是你可以把它考虑成先变形后旋转或者先旋转后变形，得到的结果都是一致的。 上面所说的是指一般意义上的连续介质力学里的定义，是把initial 和current两个坐标系固定在一起，包括原点和轴 向，但是既然ABAQUS里面存在corotational坐标系，就是两个坐标系之间存在一个旋转，这样我们再采用变形梯度 F的话就相当于多旋转了一个R，所以要旋转回来。 2，对应力的定义，我得理解是traction force不因材料变形影响的，但是这个force作用的面积dA和单位正方向n是 受材料变形影响的，这就是为什么会有1st和2nd PK应力，至于configuration和坐标系，我们常说initial configuration 和 deformed configuration 但是说坐标系只能说initial coordinates 和current coordinates，是 因为configuration是指材料的构型，或者说是 geometrical state of the body。 但是这里既然存在一个corotational坐标系，而且是一个非惯性系，那么这个traction force应该也是跟坐标系的选 取有关了，需要和rotation rate 相关。 3，这个我同意，对于solid element来讲无论嵌入的模型还是Umat采用的都是Jaumann rate，但是VUMAT采用 Green-Naghdi rate，而Explicit本身又是采用Jaumann rate的，这可能就是为什么Umat的结果没问题，而 Vumat里面定义的应力和odb里面的应力不同的原因。 还是谢谢shawn2008 的热心答复，我一直都没想到率相关这个问题，而且坐标系和configuration也没有很好的区 分，不知道现在我这样理解是不是正确。 shawn2008 发表于 2009-4-10 06:40 在继续讨论之前，我先对本构的理论做个简单的总结。 本构理论建立的基础是热力学的三大定律，它的根本任务是建立‘力’与‘变形’之间的关系，这里需要注意的是‘变形’是 一个可以观测的量，‘力’是不可观测的，尽管这样这里我们先接受这两个概念，至于如何得到的这两个概念我们不在 这讨论了。 要定量化描述‘力’和‘变形’之间的关系，我们需要先选定一个坐标系，然后才是给出这两个概念的具体定 义。 假设我们已经选定了这样的全局坐标系，并且在经过人们的不断探索后，发现‘力’的描述用‘应力’比较合适，‘变 形’的描述用‘变形梯度’比较合适。 ‘应力’和‘变形梯度’给出的是静态的描述， 就是在时空中，某一个瞬间的变形体变 形情况。其实，这时我们如果能给出每个时刻，‘应力’全量与‘变形梯度’全量之间的一一映射关系，本构模型就建立完 了，但实际上对于不可逆变形过程，存在加载的历史相关性，所以我们很难给出全量之间的本构关系（不过线弹性材 料就可以给出这样的关系，所以线性弹性问题的本构关系比较简单）。 为了描述这类材料的本构关系，人们想到了率 （时间的导数），如果能够发现‘应力’的率和‘变形梯度’的率之间的关系，然后通过时间积分，就可以得到变形体变形 后的‘应力’和‘变形梯度’。 所以就有了‘应力率’和‘速度梯度’这两个描述动态性质的概念。也就是‘先静后动’， 材料的 变形是在一个四维的时空中， 用静态的两个概念‘应力’和‘变形梯度’描述三维空间的属性，用动态的两个概念‘应力率’ 和‘速度梯度’描述另外的一个时间维度， 这样四个定义就完全覆盖了这个四维时空。 本构理论的客观性原理要求我们建立 旋转无关的那部分‘应力率’和旋转无关的那部分‘速度梯度’之间的关系。旋转相 关的‘应力率’和‘速度梯度’里面都涉及了一个旋率的概念，即旋转对时间的导数， 旋率与旋转的关系是， \Omega=\dot{Q}Q^T 这个描述转动的Q就是corotational坐标系在材料变形过程中对应的转动。 注意这里我用 的是Q来描述corotational坐标系的转动，而不是用的刚体转动R，因为这两个是有区别的。 如何求出Q，只需我们 给出Omega，我一直强调总的变形是可以观察到得，但里面有多少是旋转，有多少是与本构相关的变形，不同的人 有不同的定义，比如Jaumann定义了旋率W，这个W是速度梯度的非对称部分，还有其他的旋率，比如Mandel旋。 到现在我没有强调任何坐标系，这里我们所有的量都是在全局坐标系下描述的。 这个分解过程可以简单理解一下，比 如速度梯度L是在全局坐标系下表达的，Jaumann定义了， L=D+W, 当然D和W也是在全局坐标系下表示的了，其 中W是与旋转相关的旋率定义，D是与本构特性相关的变形。 经常你在有的论文中会看到，有的本构建立时符号用的 是D，有的用的是\dot{\varepsilon}，这两个是有区别的，D的表达我已经解释了， \dot{\varepsilon}是在共旋 坐标系下表达的，全局坐标系与共旋坐标系之间的转动是Q，因此，\dot{\varepsilon}=Q^T D Q。 vumat和 umat中说的应变和应力是已经旋转后的，就是像\dot{\varepsilon}这种在corotational坐标系下的表达。 上面的就是本构模型建立思路，如果想与abaqus软件中的很多结果对应起来，还需要澄清一些基本的概念，比如 configuration和坐标系，变形梯度的极分解和乘法分解，拉格朗日坐标系与欧拉坐标系等等。 以后有时间在写了。 shawn2008 发表于 2009-4-10 06:50 [i=s] 本帖最后由 shawn2008 于 2009-4-10 06:53 编辑 [/i] 或许我的理解也有偏差，但我一直觉得本构模型里面最重要的就是基本概念。比如“应力”这个概念，如果我们再深入 的探讨下去，里面也有无数的新认识，而不像“单位面积上所受的力”这么简单。 billowriver 发表于 2009-4-11 06:21 跟shawn2008 讨论真是受益匪浅，你讲的非常好， 不过我记得连续介质和本构理论是建立在9个基本守恒方程基础上的。 质量守恒(1)，动量守恒（3），角动量守恒（3），能量守恒（1）和热力学第二定律也就是熵增原理（1） 采用Jaumann的定义，Q和R确实不同，而且实际上Jaumann rate是通过把速度梯度分解成对称的变形相关矩阵D 和反对称转动相关阵W，这样W之和Q有关，与R确实没有直接关系。但是Green-Naghdi rate里面确实直接定义 Q=R，这时旋率就是dot(R)*R^T。 我也觉得采用Jaumann rate定义物理意义更合理一些，因为材料是一直一边变形一边转动的，如果是简单的拉压变 形，二者应该没任何区别，但是如果有剪切变形，因为剪切变形会产生额外的旋转角和旋转加速度，因此二者就可能 会产生很大差别。 但是考虑到剪切也会造成一部分旋转，所以本构的提法是不是应该改成“与刚性旋转无关的那部分‘应力率’和与刚性旋 转无关的那部分‘速度梯度’之间的关系”？ 不过对于梁板壳这些结构单元，Green-Naghdi rates应该更容易求解，但是为什么对于solid element，ABAQUS 单独把Vumat单独拿出来采用Green-Naghdi rates？ 不管怎样，这个Q在一般情况下是不能等同于R的，但是在VUMAT里面，因为采用了Green-Naghdi rates这个Q就 是R。因此corotational coordinates与global coordinates相互相差了一个R. 我觉得采用corotational坐标系的好处应该是，当你计算应力变化率时，你只需要考虑应力对时间的导数，而不需要 考虑坐标系对时间导数的影响，不知道是否准确。 还有，在Vumat里面我们是已知应变增量直接计算出应力，而不是像Umat里面那样需要给出Jacobi，即 d_stress/d_strain，或者说是stress_rate/strain_rate。VUMAT里面的应变增量很好计算，直接对D从tn到 t_n+1进行积分就可以，因此D是受到旋转角度Q的影响，这一步ABAQUS已经自动完成，但是Vumat里面计算得到 的应力为什么又会和odb里面给出的应力不同？难道是因为这个应力是在corotational 坐标系下得到的，而这个 corotational 是对应Green-Naghdi rates，即旋转张量R的，而在odb输出时ABAQUS需要先把它旋转回global 的 坐标系，然后再旋转Q，采用Jaumann定义的另一个corotational 坐标系？ shawn2008 发表于 2009-4-11 16:42 [i=s] 本帖最后由 shawn2008 于 2009-4-11 23:11 编辑 [/i] 1.不过我记得连续介质和本构理论是建立在9个基本守恒方程基础上的。 质量守恒(1)，动量守恒（3），角动量守恒 （3），能量守恒（1）和热力学第二定律也就是熵增原理（1） 质量守恒(1)，动量守恒（3），角动量守恒（3）， 这几个是介质力学的基本原理；热力学第一，二原理是与本构方 程的直接表达有关的，状态势和耗散势； 我个人认为绝对零度的零熵原理是第二原理的一个补充。 以上的说法是我 个人的理解与习惯，你的表述更经典一些。 2. 采用Jaumann的定义，Q和R确实不同，而且实际上Jaumann rate是通过把速度梯度分解成对称的变形相关矩 阵D和反对称转动相关阵W，这样W之和Q有关，与R确实没有直接关系。但是Green-Naghdi rate里面确实直接定 义Q=R，这时旋率就是dot(R)*R^T。 不同人有不同的定义，没有什么是很绝对的。 3. 我也觉得采用Jaumann rate定义物理意义更合理一些，因为材料是一直一边变形一边转动的，如果是简单的拉压 变形，二者应该没任何区别，但是如果有剪切变形，因为剪切变形会产生额外的旋转角和旋转加速度，因此二者就可 能会产生很大差别。 通常来说，材料的变形过程都是即有旋转又有变形的。在剪切变形的时候，如果剪应变大到一定程度，这些率之间会 存在一定的差异。 不过 ‘因为剪切变形会产生额外的旋转角和旋转加速度’， 你这个指的是什么我就不太明白了。 4. 但是考虑到剪切也会造成一部分旋转，所以本构的提法是不是应该改成“与刚性旋转无关的那部分‘应力率’和与刚性 旋转无关的那部分‘速度梯度’之间的关系”？ 剪切会造成旋转，我们在客观率中不是已经消除掉这种旋转效果了吗？ 我笼统的称为旋转，是因为有些人认为塑性 旋也应该考虑到对最后旋转的影响。 不过我不知道你这里的刚体旋转指的是什么？ 这里有个问题需要你确认， 一 个单元只作剪切变形，corotational坐标系发生变化了吗？ polar decomposition 中F＝RU， R是零张量吗？ 5. 不过对于梁板壳这些结构单元，Green-Naghdi rates应该更容易求解，但是为什么对于solid element， ABAQUS单独把Vumat单独拿出来采用Green-Naghdi rates？ 不管怎样，这个Q在一般情况下是不能等同于R的， 但是在VUMAT里面，因为采用了Green-Naghdi rates这个Q就是R。因此corotational coordinates与global coordinates相互相差了一个R. 在剪应变不是非常大的情况下，这几个率的结果不会有太大差异。 为什么你说'Green-Naghdi rates应该更容易求 解'？ 6. 我觉得采用corotational坐标系的好处应该是，当你计算应力变化率时，你只需要考虑应力对时间的导数，而不需 要考虑坐标系对时间导数的影响，不知道是否准确。 本构相关的应力率都是客观的，都不用考虑应力对时间的导数，这是客观性原理的要求。 为什么说选corotational坐 标系的好处？ 7. 还有，在Vumat里面我们是已知应变增量直接计算出应力，而不是像Umat里面那样需要给出Jacobi，即 d_stress/d_strain，或者说是stress_rate/strain_rate。 Jacobi矩阵是平衡方程迭代的需要，用来保证收敛速度的。 除此以外，umat和vumat在应力更新上都是一样的。 8. VUMAT里面的应变增量很好计算，直接对D从tn到t_n+1进行积分就可以，因此D是受到旋转角度Q的影响，这 一步ABAQUS已经自动完成，但是Vumat里面计算得到的应力为什么又会和odb里面给出的应力不同？ odb中的分量是表达在local system下，abaqus用了3种方式的local system。 umat和vumat中的应力如果你转 换到与odb表达相同的local system下，结果是一样的。 9. 难道是因为这个应力是在 corotational 坐标系下得到的，而这个corotational 是对应Green-Naghdi rates， 即旋转张量R的，而在odb输出时ABAQUS需要先把它旋转回global 的坐标系，然后再旋转Q，采用Jaumann定义 的另一个corotational 坐标系？ 忘掉极分解的R， 本构方程建立时，只有Q， 只不过在特定的旋率定义下，Q＝R而已。 umat和vumat的应力\sigma_n的确是已经旋转过的了，但这个应力只是corotational的，其表达还是在local system 下， 而不是在corotational 坐标系下。 对不同的单元abaqus采用了不同的方式计算变形梯度。 你可以在 材料子程序里面对DFGRD0和DFGRD1分别进行一下F=RU＝VR分解，看看变形梯度里面包含了多少R，特别是 shell单元。 最后，如果能从DFGRD0, DFGRD1还有DTIME计算得到与dstran一样的结果，你的所有问题就完全清晰了。 billowriver 发表于 2009-4-11 23:54 3，和4，我是觉得对于刚体转动而言，就不存Q和R的不同，但是对于变形体，哪部分算作旋转哪部分算作剪切变形 就不好说了，所以即使对于纯剪切，如果采用F=UR，此时R是单位矩阵，就是材料不做任何刚性旋转，根据Green- Naghdi rate的定义Q就是R，但是如果用Jauhman的定义，Q和R是有很大差别的，所以即使是纯剪切也会有旋转角 度。可能我的表达也有问题。 5，但是为什么ABAQUS里面采用了这个corotational坐标系，在global坐标系下表达不好么。 7，我是记得有限元里面有个D矩阵，即d_stress/d_strain，所以实际上我们定义材料本构子程序就是定义这个D， 和Umat是一致的 8，你说的这个local system是指什么，是跟单元属性有关么，比如solid element，它的local system是什么呢 不过我Umat没怎么用过，我去试试看 shawn2008 发表于 2009-4-12 01:07 [i=s] 本帖最后由 shawn2008 于 2009-4-12 20:11 编辑 [/i] 一，3，和4，我是觉得对于刚体转动而言，就不存Q和R的不同，但是对于变形体，哪部分算作旋转哪部分算作剪切 变形就不好说了，所以即使对于纯剪切，如果采用 F=UR，此时R是单位矩阵，就是材料不做任何刚性旋转，根据 Green-Naghdi rate的定义Q就是R，但是如果用Jauhman的定义，Q和R是有很大差别的，所以即使是纯剪切也会 有旋转角度。可能我的表达也有问题。 Q 和 R的不同是由于来自于不同的旋率的定义。 对于纯剪切变形体， 由于我们能够观察到的只是总变形，至于里面 有多少是剪切变形，有多少是旋转，我们无法观察，就是因为这个原因，才出现了很多种旋率的定义。 现在被大家接受的有两种建立唯象本构的体系： （1） 材料当前的状态由可观测变量的当前值和历史值决定；（2） 材料当前的状态由可观测变量的当前值和一组内变量决定。 我们现在多数使用的都是第二种。 可观测变量指的是应 变和温度这类， 至于弹性应变，塑性应变，背应力，等效塑性应变这些都是不可观测的，所以被称为内变量。 这些 不可观察的东西，我们就很难判断是否是正确的，虽然我们可以通过一些实验间接验证。 我说这个的意思就是说， 力学里面的很多东西都是人们根据实验结果和数学演绎得到的这么一个自洽的理论体系，没有绝对正确的东西。 就像 一个变形过程，我们能看到的只是变形体的运动，至于里面的旋转，纯变形都是人为的分出来的。 本构要求的是剔 除旋转效果后的那部分变形，对剪切变形，这里面当然包括了旋转和纯剪切，极分解的本质是把一个特殊变形路径下 的变形过程分解为旋转和纯剪切。 F=RU分解后R不是单位张量。 附件中的剪切变形过程，你可以计算一下R。 二，5，但是为什么ABAQUS里面采用了这个corotational坐标系，在global坐标系下表达不好么。 如果没有材料用户子程序，是表达的global下比较好，只是为了方便。 对于金属板料而言，其有各向异性的对称中 心，比如RD,TD,ND，一般我们的各向异性屈服函数都是在这个标架下表示的，但这个标架通常和corotational坐标 系不重合。 还有比如shell单元，其位移场的描述在单元的几何坐标系下表达更方便一些，这个几何坐标系和local system常常也不一样，也要进行一些量的变换。 总之，软件里面考虑的实际情况很复杂的。 三，7，我是记得有限元里面有个D矩阵，即d_stress/d_strain，所以实际上我们定义材料本构子程序就是定义这个 D，和Umat是一致的 本构模型的根本任务是根据当前值和应变增量更新应力。 umat中的一致性切模量(算法切模量)是为了保证整体平衡 方程迭代的收敛性的。 四，8i，你说的这个local system是指什么，是跟单元属性有关么，比如solid element，它的local system是什么 呢。 你可以到手册中找找具体的解释， local system就是一个坐标系，是我们表达应力，应变这些张量的坐标系，odb 中应力分量的11，22，33......就是相对local system而言的。 五，vumat也是一样的，接口中也传入了DFGRD0，DFGRD1，DTIME和dstran。 billowriver 发表于 2009-4-13 08:15 多谢shawn2008 这么精彩的解释， 但是这里有个问题， 在ABAQUS的user manual1.2.2 Conventions里面提到 The 1-, 2-, and 3-directions depend on the element type chosen. For solid elements the defaults for these directions are the global spatial directions. For shell and membrane elements the defaults for the 1- and 2-directions are local directions in the surface of the shell or membrane, 但是又特别强调对于ABAQUS/Explicit， For geometrically nonlinear analysis（Explicit里面默认的就是geometrically nonlinear ）, the internally stored components rotate with the material, regardless of whether or not a user-defined orientation is used. 在ABAQUS document user subroutine里面的1.2.5 VUMAT.: All stresses, strains, stretches, and state variables are in the orientation of the local material axes. These local material axes form a basis system in which stress and strain components are stored. This represents a corotational coordinate system in which the basis system rotates with the material. 所以以上应该说明solid element在Abaqus/explicit 里面，local coordinates 应该就是corotational coordinates，但是对于Abaqus/Standad，local coordinates就是global coordinates. 我用Standard比较少， 不知道在Standard里面是不是这样的。而对于shell或 membrane element等的local coordinates是由几何形状 决定的。 所以对于odb里面的local coordinates 应该是基于Jaumann rate 的local coordinate？所以才会与Vumat里面 的应力不同？ABAQUS里面也提到在剪切变形很大的时候，Jaumann和Green-Naghdi rates的结果相差很大。 附件里面的剪切确实不是纯剪切，F=[1, sin(alpha),0 ; 0, 1, 0; 0, 0, 1]，而对于纯剪切 F=[cos(alpha), sin(alpha),0 ; sin(alpha), cos(alpha) 0; 0, 0, 1]，R是单位张量，F=U shawn2008 发表于 2009-4-13 19:09 [i=s] 本帖最后由 shawn2008 于 2009-4-13 22:00 编辑 [/i] corotational坐标系在t=0时刻是与local system重合的，如果local system相对于global 坐标系固定不动，odb 与材料子程序中的应力就很好对应了。 如果local sytem 也以一定得旋率在转动，在这种情况下，如果local system旋率与corotational坐标系的旋率一样，那么在变形过程中，每一时刻local system都与corotational坐标 系重合；当旋率不一样时，两个坐标系就不重合了。我倾向于把这几个坐标系区分开，即使有的时候它们重合，但物 理意义是完全不同的。 你这个剪切变形的确是纯剪切，可以拿这个来对应一样，global坐标系，local system和corotational坐标系。 F=[1, sin(alpha),0 ; 0, 1, 0; 0, 0, 1]，这个变形梯度应该是F=[1, -tg(alpha),0 ; 0, 1, 0; 0, 0, 1]，如果 alpha定义为正方形与y轴重合边的逆时针转角。 shawn2008 发表于 2009-4-13 19:14 我们讨论了半天，“在ABAQUS/Explicit里面，Vumat里的所有变量都是相对于一个corotational 坐标系 来说的， 而在ABAQUS的odb里面显示的是cauchy stress, 是相对于当前坐标系来讲的。”这句话的前半句算是基本上讨论清 楚了。 :) linrn 发表于 2009-4-13 20:33 写的真多 看来看去就是连续介质力学里面的基本概念 billowriver 发表于 2009-4-13 22:06 是啊，再次佩服shawn2008 的理论功底，很多东西在一般的连续介质力学 书 关于书的成语关于读书的排比句社区图书漂流公约怎么写关于读书的小报汉书pdf 里都是找不到的。 这次弄清楚了很多问题。 看来如果真要完全弄清楚ABAQUS里面这个local coordinates 是什么只能以后有时间的时候做些相关的数值模拟 了。在此之前我就先认为是由于旋率不同造成的吧。 那个是 tg(alpha)，是我看错了，以为边长是不变的 shawn2008 发表于 2009-4-13 22:25 对，就是要不断的动手和思考才能搞清楚。再过一段时间，你在重新读一下这个帖子，你会发现我们很多表述可能是 错的或者说是不严谨的。 很多基本概念就是在不断的认识，推翻，再认识中透彻的。我只不过多看了一些资料，然后 自己做了一些加工，其实，在我们的讨论中让我对很多概念也加深了印象。 你可以看一下，[url]http://imechanica.org/[/url] 这个论坛上的第二个帖子也是搞不清楚坐标系的问题。 billowriver 发表于 2009-4-14 06:06 [i=s] 本帖最后由 billowriver 于 2009-4-14 06:07 编辑 [/i] 不知道前一个帖子为什么不在了，重贴一下！ 佩服shawn2008 的理论功底，有些概念连一般的连续介质力学书里面都没有详细介绍， 看来要想彻底弄清楚ABAQUS，odb里面的local coordinates，只能有时间进行一些数值模拟了，在此之前我先假 定Vumat和odb里面应力的不同是由于旋率不同造成的吧。 那个变形梯度里面的是-tg(alpha），我以为边长没变呢，看错了 页: [1] 查看完整版本: ABAQUS 子程序VUMAT 的坐标系的问题的讨论 Powered by Discuz! Archiver 7.0.0 © 2001-2009 Comsenz Inc. 设计 领导形象设计圆作业设计ao工艺污水处理厂设计附属工程施工组织设计清扫机器人结构设计 仿真核心增值经销商 北京盛维提供3D CAD/CAE/CAM/PDM/CA... 整体解决方案。010-68250396/68253809 www.shengweitech.com 买输油臂请选江苏德威 我公司是"中国鹤管生产基地"欢迎前来 选购. 全国免费热线:400-699-1909 www.dmco.com.cn 选择热喷涂-认准伟俊喷涂 我们是一家专业从事喷焊技术生 公司,公司坚持用心去做,以质赢 www.wjrpt.com