浅析有限元法中的沙漏现象

本文转自-公众号有限元仿真分析
原文链接：https://mp.weixin.qq.com/s/mHe25cOXFLV1s-bULEZvPg
沙漏 (Hourglass) 模式是一种非物理的零能变形模式，产生零应变和应力。
在有限单元法的力学分析中，一般以节点的位移作为基本变量，单元内节点的位移以及应变均采用形函数对各点位移进行插值计算得到。应力根据本构方程由应变计算得到，之后就可以计算单元的内能了。如果采用单点积分（积分点在等参元中心），在某些情况下节点位移不为零（即单元有形变），但插值得到的应变却为零。比如，一个正方体单元变形为等腰梯形，节点位移相等却方向相反，各点的形函数为零，所以插值结果为零，这样内能计算结果也为零（单元没有变形）。在这种情况下，一对单元叠在一起有点像沙漏，所以这种模式被称之为沙漏模式或者沙漏。如果单元变成交替出现的梯形形状（两两在一起类似沙漏以及Windows系统中的鼠标动画图标），这时就需要小心了。
为了说明问题，首先假定选择一个弯矩作用来模拟纯弯曲荷载的一小块材料。在弯矩作用下，材料中轴线处的长度没有改变，与纵向轴线的夹角也没有改变。这意味着单元单个积分点上的所有应力分量均为零。由于单元变形，没有产生应变能，因此，这种变形的弯曲模式是一个零能量模式。由于单元在此模式下没有刚度，所以，单元不能抵抗这种形式的变形。在粗划的网格中，这种零能量模式会通过网格扩展，从而产生无意义的结果。
一般来说，如果从变形的网格中看不出沙漏效应的话，就认为它造成的影响不大。一个更为量化的途径就是研究伪应变能。它是控制沙漏变形所耗散的主要能量。如果伪应变能过高，说明过多的应变能可能被用来控制沙漏变形了。判断过高伪应变能的来源，最有效的途径是比较伪应变能和其他内部能量的值。一般而言，伪应变能与实际应变能的能量耗散比率应低于5%。处于完全弹性范围阶段的固体，由变形而存储的能量称为应变能，当外力消除时，应变能将释放做功，变形体恢复原状。
总的来说，作用在单独节点上的荷载、边界条件或接触易于产生沙漏现象。而将荷载或约束分布在两个或更多节点上，则大大减轻了沙漏问题，如对角部进行圆角处理等。一般的规律是，两点或多点接触发生得越快，沙漏将减弱的越快。
总能量=内能+动能+滑移界面能
能量之间是可以转化的。但是，对于动力学问题，总能量一般是不变的，也就是能量守恒原理。沙漏模式也就是零能模式在理论上是存在的，大多数实际的模型中是不可能的。零能模式是指有变形，但不消耗能量。显然，这是一种伪变形模式，若不加以控制，计算模型会变得不稳定，并且计算出来的结果是没有意义的。要抵制这种变形模式就需要消耗一定的能量，也就是沙漏能。如果，这个比值太多，就说明计算模型与实际模型的变形有很大差距，当然结果也就是不正确的。这也是使用缩减积分所付出的代价。用完全积分单元可以解决这个问题，但是计算效率不高，还有可能导致单元锁死，过刚度等问题。
由于沙漏问题与网格质量息息相关，虽然几何模型的划分方法至今没有一个统一的标准可以遵循，但是利用经验的积累，以下推荐几项改善网格质量的原则与方法：
01
对于多大数问题而言，采用线性缩减积分单元的细划网格产生的误差可在一个可接受的范围之内。建议当采用这类单元模拟承受弯曲荷载的任何结构时，沿厚度方向上至少应采用四个单元。线性缩减积分单元能很好地承受扭曲变形，因此，在任何扭曲变形很大的模拟中可以采用网格细划的这类单元。
02
标准单元的边长通常以几何模型的最下尺寸确定，即如果几何模型的厚度是结构的最小尺寸，那么标准单元的边长至少应与此厚度相当。在高应力梯度区域的单元应细分，单元大小取决于计算精度与规模等。
03
在高应力梯度区，要进行网格细分的应力稳定性计算。既采用多次，多类型的局部细化网格进行计算，当前后两次计算结果满足计算精度要求时，通常差别不大于5%，用以确定合适的网格。
04
网格划分时，单元各边之间的比例不能太大。对于线性单元，如四节点四边形单元，八节点六面体单元等，要求小于3；对于二次单元，例如十节点四边形单元，20节点六面体单元等，要求小于10等。对于梁结构，在两个节点之间可根据需要划分多个单元。