2. 应力状态与失效模式判据

在LS-DYNA的失效建模中，仅凭等效应力或等效塑性应变往往无法解释同一材料在不同载荷下呈现的截然不同破坏形式。例如，同一种铝合金在光滑圆棒拉伸时表现为韧性颈缩断裂，在缺口试样中却可能因高静水压力而更晚开裂，而在剪切试样中又会沿剪切带迅速失效。造成这种差异的核心，是应力状态本身。描述应力状态对失效影响的两个关键无量纲量，就是应力三轴度与Lode参数。

应力三轴度：静水压力与剪切的相对强弱

应力三轴度定义为平均应力（静水应力）与von Mises等效应力之比：

η = σ_m / σ_eq
σ_m = (σ_1 + σ_2 + σ_3) / 3
σ_eq = sqrt( 0.5 * [ (σ_1-σ_2)^2 + (σ_2-σ_3)^2 + (σ_3-σ_1)^2 ] )

其中σ_1、σ_2、σ_3为三个主应力，符号约定为拉正压负。η的物理含义是“静水拉伸（或压缩）相对于剪切驱动力的比例”。η越大，材料内部越倾向于通过微孔洞扩张而破坏；η越小甚至为负，则孔洞被压合，失效更多由剪切滑移主导。

典型应力状态下的η量级如下（示例值，用于建立数量级直觉）：

Lode参数：区分拉伸型与压缩型剪切

应力三轴度只描述了“静水部分”，却无法区分同样是低三轴度下的两种不同剪切状态：一种是“拉伸型剪切”（如纯剪切中一个主应力为正、一个为负、一个为零），另一种是“压缩型剪切”。Lode参数正是用来刻画偏应力张量在π平面上的位置：

L = (2σ_2 - σ_1 - σ_3) / (σ_1 - σ_3)   （σ_1 ≥ σ_2 ≥ σ_3）

L的取值范围为[-1, 1]。L = -1对应广义拉伸（σ_2 = σ_3），L = 0对应纯剪切，L = +1对应广义压缩（σ_1 = σ_2）。在韧性断裂中，Lode参数影响孔洞与剪切带的竞争：L接近-1时孔洞扩张占优，L接近0时剪切局部化更易发生。因此现代失效模型（如Modified Mohr-Coulomb、Hosford-Coulomb、DF2016等）通常同时以(η, L)作为断裂应变的输入变量。

应力状态如何改变失效模式

把η和L放在一起，可以形成一个“失效模式地图”的直觉：

• 高η、L≈-1：孔洞扩张主导，断裂应变随η升高而下降，断口呈韧窝。
• 低η、L≈0：剪切主导，断裂应变对η不敏感，但对L敏感，断口呈剪切唇或平坦剪切面。
• 负η：孔洞被压合，损伤累积被抑制，材料可承受很大压缩应变而不失效。
• L由-1向+1变化：在相同η下，压缩型剪切（L>0）通常比拉伸型剪切（L<0）更难触发孔洞型损伤，但可能更易形成绝热剪切带。

这意味着，如果只用一条“等效塑性应变—应力”曲线标定失效，就会在剪切主导的工况中严重高估或低估寿命。正确做法是设计多组不同(η, L)的试样，分别标定断裂应变曲面。

在LS-DYNA中输出与解读η和L

LS-DYNA本身不会自动把η和L作为历史变量输出，需要借助*DATABASE_EXTENT_BINARY与用户自定义输出，或使用*INITIAL_STRESS、*MAT_ADD_EROSION配合状态变量。更常用的做法是通过*DATABASE_HISTORY_SOLID与*DATABASE_BINARY_D3PLOT输出应力张量，再在后处理中计算。以下给出一个典型的输出配置片段：

*DATABASE_EXTENT_BINARY
$   NEIPH   NEIPS   MAXINT   STRFLG   SIGFLG   EPSFLG   RLTFLG   ENGFLG
        0       0         8        1        1        1        1        1
$   CMPFLG   IEVERP   BEAMIP   DCOMP    SHGE     STSSZ
        1        0        0        0        0        0

*DATABASE_HISTORY_SOLID
$      ID1      ID2      ID3
         1        2        3

*DATABASE_BINARY_D3PLOT
$       DT      LCDT
     1.0E-4

其中STRFLG=1表示输出应力张量分量，SIGFLG=1表示输出应力。在LS-PrePost中，可以通过“Fringe Component”选择von Mises应力，但η和L需要手动计算。一个实用的做法是：在LS-PrePost的“ASCII”中导出单元历史变量，用Python或Excel计算每个时间步的η和L，再绘制η–等效塑性应变曲线。

下面给出一个可直接参考的Python计算片段（示例，假设已从LS-PrePost导出sx, sy, sz, sxy, syz, szx）：

import numpy as np

def triax_lode(sx, sy, sz, sxy, syz, szx):
    # 构造应力张量
    S = np.array([[sx, sxy, szx],
                  [sxy, sy, syz],
                  [szx, syz, sz]])
    # 主应力（降序）
    p = np.sort(np.linalg.eigvalsh(S))[::-1]
    s1, s2, s3 = p
    sm = (s1 + s2 + s3) / 3.0
    seq = np.sqrt(0.5*((s1-s2)**2 + (s2-s3)**2 + (s3-s1)**2))
    eta = sm / seq if seq > 1e-12 else 0.0
    denom = (s1 - s3)
    L = (2*s2 - s1 - s3) / denom if abs(denom) > 1e-12 else 0.0
    return eta, L

解读时重点关注三点：第一，η是否在加载过程中显著变化，若变化剧烈，说明单一η标定不适用；第二，L是否接近0，若接近0则必须考虑剪切型失效判据；第三，η与L的组合是否落在已标定的断裂应变曲面范围内，若超出范围则外推需谨慎。

把上述流程落到具体操作上，可以按以下顺序推进：先在LS-PrePost中打开d3plot，选择“Fringe Component”确认应力输出已开启；再通过“History > Element”选取关注单元，导出ASCII数据；然后用上面的Python函数逐时间步计算η和L，并与等效塑性应变对齐；最后绘制η–ε_p曲线，观察失效点处的η和L值。若失效点落在高η、L≈-1区域，应优先采用孔洞扩张型判据；若落在低η、L≈0区域，则应改用剪切型判据或同时包含η和L的耦合模型。

需要强调的是，η和L都是当前应力状态的瞬时量，而失效是历史累积过程。因此在实际建模中，通常将η和L作为损伤演化方程的输入变量，而不是直接用瞬时值判断失效。例如GISSMO模型通过一个与η和L相关的等效塑性应变阈值来累积损伤，当损伤变量达到1时触发失效。理解这一点，才能正确使用LS-DYNA中的*MAT_ADD_EROSION、*DAMAGE_INITIATION等关键字。