含损伤型 Johnson-Cook 本构模型迭代算法在 LS-DYNA 中的嵌入

本文利用商业软件 LS-DYNA 提供的自定义子程序接口，对 LS-DYNA 进行二次开发，形成新的求解器，文中提出的本构迭代算法基于 Johnson-Cook 本构模型，原有的模型可描述金属材料的应变率效应、温度效应和破坏效应。之后开展穿甲弹贯穿金属靶板的数值模拟研究，新求解器中嵌入了含损伤型本构模型，采用本构迭代算法，揭示侵彻过程中的主要物理图象，计算结果和实验结果符合较好。

1 本构方程描述及材料参数

1.1 含损伤型本构方程

含损伤型 Johnson-Cook 模型将具有如下形式

\[\bar \sigma = \left( A + B{\left( {\bar \varepsilon}^{\rm{p}} \right)}^n \right)\left({1 + C\ln {\dot \varepsilon}^ \ast } \right)\left(1 - {T^{ \ast m}} \right)\left(1 - D\right)\]

\[{\rm{d}}{{\bar \varepsilon}^{\rm{p}}} = {\phi ^{n + 1}}/\left({3G + \frac{{\partial \bar \sigma}}{{\partial {{\bar \varepsilon}^{\rm{p}}}}}} \right)\]

其中 \(\bar \sigma\) 为应力强度，\(\theta _T=(T-T_0)/(T_m-T_0)\)，\(T_0\) 是室温，\(T_M\) 是材料的熔化温度。该本构模型考虑了应变率硬化行为、其温度升高导致材料的热软化和损伤软化行为。

由于是绝热过程，由能量守恒方程可得温升方程

\[\dot T = \alpha \frac{v}{c_v}\mathit{\boldsymbol{s}}:{{\mathit{\boldsymbol{\dot e}}}^{\rm{I}}}\]

其中 \(c_v\) 是 \(\alpha \) 是塑性功转换成热的比例因子。\(s\) 为偏应力张量，\(e^{I}\) 为不可逆偏应变率张量。

至于压力的计算，在低压下，可以采用如下状态方程

\[P = K\theta /3\]

其中 \(\theta\) 为体应变。

1.2 损伤演化方程

损伤演化方程可定义为

\[\dot D =\left\{\begin{matrix} 0 &p < {p_{\rm{d}}} \\ \frac{D_c}{p_f- p_d} \dot p & p \ge {p_{\rm{d}}} \end{matrix}\right.\]

其中

\[{p_{\rm{f}}} = \left[{{D_1} + {D_2}\exp \left({{D_3}{\sigma ^ \ast}} \right)} \right]{\left( {1 + {{\dot p}^ \ast }} \right)^{{D_4}}}\left({1 + {D_5}{T^ \ast}} \right)\]

式中 \(\sigma^ \ast =\sigma_mm / \sigma_{eq}\)，\(\dot p^\ast =\dot p/ \dot \varepsilon_0 \)，\(D_1 \sim D_2\) 为材料参数，\(\sigma_m\) 为静水压，\(\sigma_{eq}\) 为等效应力，p 为塑性应变，\(\dot p\) 为塑性应变率，\(p_d\) 为损伤阈值，\(D_c\) 为损伤阈值。

当 \(D=D_c≤1\) 时材料破坏，\(D_c\) 为小于等于 1 的材料参数，一般取 0.3 左右。

2. 计算流程

如图所示，已知 n 时刻的内外变量 \((\sigma_n, \xi_n)\)，速度场 \(u\) 在主程序中由运动方程求得，dε 在主程序中由连续方程求得，其余由本构算法在自定义子程序中求得 n + 1 时刻的内外变量，从而可进入下一个循环。