2.18 正交各向异性材料的常见问题

LS-PrePost里加入了“Fiber Dir”，从代码中查看可以得出：只能用于mat36、mat234和mat235，需要一起读入k文件和d3plot。单元的前两个历史变量包含两个角度（用度度量），意味着需要输出至少两个历史变量到d3plot，才能显示纤维方向。

纤维方向一是通过材料x轴和绕着壳单元法向旋转向量histvar[0]得到的。纤维方向二用同样的方法，但是是绕着壳单元法向旋转向量histvar[1]得到的。

LS-PrePost 中显示的材料坐标系是直角坐标系。它的方向基于输入文件在 t=0 时刻被初始化，在计算过程中，根据单元坐标系的旋转而旋转。设置不同的*control_accuracy的INN值，旋转的算法有区别。这样的材料坐标系用于LS-DYNA中所有的亚弹材料模型。

然而，如果使用了一个超弹材料模型，需要计算总应变/变形，材料坐标系会跟随材料（而不是单元坐标系），这意味着在变形状态中材料轴不需要是直角的。在LS-PrePost中显示材料坐标系时，没有考虑这种情况。另外，如果应力张量按照材料坐标系写到d3plot中（CMPFLG），LS-PrePost也会在材料坐标系中显示应力，对于超弹性材料模型也是这样的。

LS-DYNA 中的壳单元应力更新计算是在单元局部坐标系中，不是全局坐标系。默认情况下，局部坐标系的 x 轴沿着 1-2 边。因此意味着，默认的单元坐标系会随着 1-2边的旋转而旋转，材料坐标系也是绕着单元的 1-2 边旋转的。对于正交各向异性材料，这会引起单元往某些非物理的低能模态变形。

试想一个沿着试件方向是材料强轴的试件进行单轴拉伸实验。如果网格划分成单元的1-2边是沿着材料的弱轴方向（垂直于加载方向），很可能单元变形引起1-2边往加载方向旋转，以至于单元的能量更低。这当然会扭转单元从而产生应力和单元能量，因此1-2 边不会旋转很大，但是由于任何有限的旋转都是错的，因此我们会得到一个错的分析和丑陋的网格。

不变的节点编号选项（*control_accuracy）用了一个不同的平均方法来定义单元坐标系的x轴。这个方法是：两个向量eta和mu，在壳的平面分别由相对的壳边的中点连线组成。换句话说，eta在一个方向上把一个单元一分为二，mu在另一个方向上把这个单元一分为二。eta和mu的一半是一个向量，叫作eta+mu，单元坐标系的x轴是45°指向eta+mu的一边，而y轴方向是45°指向eta+mu的另一边。在一个矩形单元中，单元x轴平行于N1-N2边，在非矩形单元中，单元x轴不平行于N1-N2边。

不变节点编号选项定义了一个局部单元坐标系，当单元重编号时，单元坐标系正好旋转了90°，因此节点号1、2、3、4变成了2、3、4、1。

这种方法的一个明显优势是网格和材料变形是独立用于单元的编号的。对于正交各向异性材料，在面内剪切和沙漏变形时，单元坐标系旋转一致于变形。更重要的是，材料往非物理的低能模态变形的趋势被降低了。虽然不变节点编号选项对于各向异性材料给出了更好的结果，但也不是完美的，可以用于壳单元公式1、2、5、7、9、10、11和16，不能用于壳单元公式3、4、6和8。

在LS-PrePost中“Local”项指的是壳单元的局部坐标系，这个坐标系由单元的连接性决定，而不是输入文件中*define_coordinate 定义的坐标系。节点 1 到节点 2 的向量（N1-N2），是局部坐标系的x方向，局部坐标系的z方向是壳单元的法向（N1-N2和N1-N4的叉乘），局部坐标系的y方向是z和x的叉乘。

下面是定义材料坐标系的相关知识。

对于使用各向异性材料的壳单元，有 3 个选项可以用来定义材料轴的初始方向，AOPT=0、2和3（可见手册中*mat_optiontropic的AOPT部分）。在求解过程中，单元坐标系随着单元旋转和平动，所以可以认为单元坐标系和材料坐标系之间的角度不变。换句话说，材料坐标系始终是随着单元的旋转和变形更新的。因此在变形前的几何上定义材料坐标系是可以的。

对于这个讨论，材料坐标系被称为 a-b-c 坐标系，这个名字和用户手册里说的名字一致。对于壳单元，c轴与单元法向共轴，a轴在壳单元的平面上，b轴由叉乘决定，b=c×a。事实上，对于翘曲单元，a轴不在单元平面上，是在沿着c轴的投影上垂直于c轴。这个投影也是为了局部单元坐标系。

对于上一段讨论的原因，AOPT中3个可用于壳单元定义a-b-c坐标系的选项都是用来定义a轴的。

对于AOPT=0，a轴假定与局部单元坐标系x轴一致。

对于AOPT=2，a轴定义为用户输入的向量a，投影到单元表面。注意用户自定义的向量d不再使用；也注意用户自定义的向量a，不等同于材料坐标系的a轴，除非a正交于单元法向。

对于AOPT=3，a轴定义为用户定义的向量v和单元法向的叉乘，如a=v× c。给予相同的输入向量，AOPT=3定义的材料坐标系与AOPT=2定义的坐标系正好旋转90°。

当正交各向异性材料使用*element_shell_beta时，材料坐标轴由材料输入中的AOPT选项定义的方向旋转单元中的PSI角度得到单元的参考方向。单元积分点上的材料坐标轴再旋转*section_shell中的BETA角度。总结：AOPT、*element_shell_beta中的PSI和积分点的BETA角度共同定义积分点上的材料方向。

BETA 角度允许用户通过绕着法向量旋转一定角度来重定向材料坐标系，当使用材料 22、23、33、34、36、41～50、54、55、56、59 和 103 时，用户通过*section_shell中的ICOMP参数和B1、B2等参数可以为每一层单元定义一个BETA角（厚度方向的每一个积分点）。可以使用*element_shell中的BETA选项为每一个单元定义BETA角度。如果两个地方都定义了BETA角，BETA角会相加。

对于材料2、21、86和117，当AOPT=3时，对于所有的单元一个默认的BETA角可以使用材料模型中的BETA参数；但是当AOPT=0或2时不可以，*element_shell中的BETA选项值会覆盖这个默认值。*section_shell中的ICOMP参数不能用于这些材料模型，因此这些材料模型不能简单地用于层叠复合材料。

当使用*element_shell的BETA选项时，理论上任意几何的任意材料方向都可以定义。然而，没有一个通用方法自动产生正确的BETA角。这些角度大多数需要手动计算，或者可能通过用户自定义开发的程序嵌入一个公式来描述某个区域内的材料方向进行自动生成。为了在材料坐标系上输出各向异性材料的应力和应变值，设置 CMPFLG=1（见*database_extent_binary）。这个参数影响壳单元、实体单元、厚壳单元的应力输出，包括ASCII elout文件、二进制d3plot文件中的结果。

对于正交各向异性材料，AOPT 只用来建立初始的材料坐标系，它对材料坐标系如何随着时间更新没有影响。

mat2是个例外，这个材料模型用全Lagrangian公式，材料坐标系的更新基于单元坐标系的旋转。这对受单元连接影响比较大的实体单元（当实体单元发生剪切\扭曲时）可能是一个现实的问题（可以通过设置*control_accuracy中的INN=3或4解决）。对于壳单元，使用*control_accuracy中的INN消除材料坐标系受单元连接的影响。

默认情况下，局部单元坐标系的x轴是沿着节点1到节点2，y轴垂直于x轴和由节点1、节点2、节点4构成的平面，对于壳单元，z轴沿着壳单元的法向。为了使单元局部坐标系不敏感于单元连接的节点顺序，通过设置*control_accuracy 中的 INN=2（壳单元）或者INN=3（实体单元）引入不变节点编号。对于实体单元的不变节点编号，只在970版本或之后的版本中可用。

局部材料坐标系和局部单元坐标系是不一样的。这两个坐标系不需要重合，但是它们是用同样的方式来更新的。

对于126号材料，它的材料坐标系更新方式取决于所用的单元公式。单元公式0和9用了不同的方法。当单元发生大的剪切变形时，单元公式0更合适。

正交各向异性的弹性常数的说明可参见理论手册中的介绍，壳单元是平面应力单元，针对三维正交异性（实体单元）的理论不适用于壳单元。对于正交各向异性的壳单元，输入参数EC、PRCA、PRCB不使用；并且，在*mat_fabric中，GBC、GCA也不使用。

在R11和后续版本中加入了一个检查：实体单元的正交异性弹性参数对于壳单元无效。对于正交异性的壳单元，实现了一个不同的检查（Error 21413）。

运行 dev 版本时，在命令行中加上“msg=20397”可以显示出关于实体单元的正交异性参数更多的检查。

运行最新版本时，在命令行中加上“msg=21413”可以显示出关于壳单元的正交异性参数更多的检查。

对于2D壳单元的材料，为了保持数值稳定，材料常数EA、EB、PRBA必须输入，并且满足PRAB*PRBA＜1.0，这里，PRAB=PRBA*EA/EB。