\[
f(\boldsymbol{\sigma},\bar{\varepsilon}^{p})=\sqrt{\frac{3}{2}}\|\mathbf{s}\|-\sigma_{y}(\bar{\varepsilon}^{p}),
\]
\[
\dot{\bar{\varepsilon}}^{p}=\sqrt{\frac{2}{3}}\|\dot{\boldsymbol{\varepsilon}}^{p}\|=\sqrt{\frac{2}{3}}\|\dot{\gamma}\mathbf{N}\| =\dot{\gamma},
\]
\[\begin{split}
\begin{equation}
\begin{aligned}
\boldsymbol{\varepsilon}^{e}_{n+1} &= \boldsymbol{\varepsilon}^{e, \text{trial}} - \Delta\gamma\mathbf{N}_{n+1},\\
\bar{\varepsilon}_{n+1}^p &= \bar{\varepsilon}_n^p + \Delta \gamma,\\
f_{n+1}&=\sqrt{\frac{3}{2}}\|\mathbf{s}_{n+1}\|-\sigma_{y}(\bar{\varepsilon}^{p}_{n+1})=0,
\end{aligned}
\end{equation}
\end{split}\]
\[
\boldsymbol{\varepsilon}_{n+1}^{e},\quad \bar{\varepsilon}_{n+1}^p,\quad \Delta \gamma
\]
方程组求解
这一非线性方程组可使用 Newton-Raphson 方法进行求解：
记 \(\mathbf{v} = \left[\boldsymbol{\varepsilon}=\boldsymbol{\varepsilon}_{n+1}^{e},\bar{\varepsilon}=\bar{\varepsilon}_{n+1}^p,\Delta \gamma\right]^{T}\)，
记 \(\mathbf{r}=\left[\mathbf{r}_{\boldsymbol{\varepsilon}},r_{\bar{\varepsilon}},r_{f}\right]^{T}\)，其中
\[\begin{split}
\begin{equation}
\begin{aligned}
\mathbf{r}_{\boldsymbol{\varepsilon}} &= \boldsymbol{\varepsilon}^{e}_{n+1} - \boldsymbol{\varepsilon}^{e, \text{trial}} + \Delta\gamma\mathbf{N}_{n+1},\\
r_{\bar{\varepsilon}}&=\bar{\varepsilon} - \bar{\varepsilon}_n^p - \Delta \gamma,\\
r_{f}&=\sqrt{\frac{3}{2}}\|\mathbf{s}_{n+1}\|-\sigma_{y}(\varepsilon),
\end{aligned}
\end{equation}
\end{split}\]
于是，Jacobian 矩阵为
\[\begin{split}
\begin{equation}
\mathbf{J} = \frac{\partial \mathbf{r}}{\partial \mathbf{v}}=
\begin{bmatrix}
\frac{\partial \mathbf{r}_{\boldsymbol{\varepsilon}}}{\partial \boldsymbol{\varepsilon}} & \frac{\partial \mathbf{r}_{\boldsymbol{\varepsilon}}}{\partial \bar{\varepsilon}} & \frac{\partial \mathbf{r}_{\boldsymbol{\varepsilon}}}{\partial \Delta \gamma} \\
\frac{\partial r_{\bar{\varepsilon}}}{\partial \boldsymbol{\varepsilon}} & \frac{\partial r_{\bar{\varepsilon}}}{\partial \bar{\varepsilon}} & \frac{\partial r_{\bar{\varepsilon}}}{\partial \Delta \gamma} \\
\frac{\partial r_{f}}{\partial \boldsymbol{\varepsilon}} & \frac{\partial r_{f}}{\partial \bar{\varepsilon}} & \frac{\partial r_{f}}{\partial \Delta \gamma}
\end{bmatrix}.
\end{equation}
\end{split}\]
初值可选为 \(\mathbf{v}^{(0)}=\left[\boldsymbol{\varepsilon}^{e, \text{trial}},\bar{\varepsilon}_{n}^p,0\right]\). 注意，\(\bar{\varepsilon}_{n}^p\) 使用的是当前时间步上一个外层非线性迭代的结果
单方程简化
Newton-Raphson 方法可以用于求解一般的非线性方程组。然而，由于在每一个高斯积分点都需要对非线性方程组进行求解，尤其在三维情况下，方程的维数达到 8，非线性求解的难度也较大，因此计算代价非常高。在某些特定情况下，可以通过消元，将方程组简化为单个方程，从而降低求解的复杂度
由于塑性流动不可压，因此
\[
\text{tr}(\dot{\boldsymbol{\varepsilon}^{p}})=0\quad \Longrightarrow\quad \text{tr}(\Delta\boldsymbol{\varepsilon}^{p})=0，
\]
故
\[
\mathbb{C}^e:\Delta\boldsymbol{\varepsilon}^{p} = 2G\ \text{dev}(\Delta\boldsymbol{\varepsilon}^{p})+K\text{tr}(\Delta\boldsymbol{\varepsilon}^{p})\mathbf{I}=2G\ \text{dev}(\Delta\boldsymbol{\varepsilon}^{p}),
\]
其中，\(\text{dev}(\cdot)\) 表示偏张量的部分，由于 \(\text{dev}(\Delta\boldsymbol{\varepsilon}^{p})=\Delta\boldsymbol{\varepsilon}^{p}\)，故
\[
\mathbb{C}^e:(\Delta\gamma\mathbf{N}_{n+1})=\mathbb{C}^e:(\Delta\boldsymbol{\varepsilon}^{p})=2G\text{dev}(\Delta\boldsymbol{\varepsilon}^{p})=2G\Delta\boldsymbol{\varepsilon}^{p}=2G\Delta\gamma\mathbf{N}_{n+1},
\]
另一方面
\[
\mathbb{C}^e:(\Delta\boldsymbol{\varepsilon}^{p})=\mathbb{C}^e:(\Delta\boldsymbol{\varepsilon} - \Delta\boldsymbol{\varepsilon}^{e})=\mathbb{C}^e:((\boldsymbol{\varepsilon}_{n}^{e}+\Delta\boldsymbol{\varepsilon}) - (\boldsymbol{\varepsilon}_{n}^{e} + \Delta\boldsymbol{\varepsilon}^{e}))=\boldsymbol{\sigma}_{\text{trial}}-\boldsymbol{\sigma}_{n+1},
\]
于是
\[
\boldsymbol{\sigma}_{n+1} = \boldsymbol{\sigma}^{\text{trial}} - 2G\Delta\gamma\sqrt{\frac{3}{2}} \frac{\mathbf{s}_{n+1}}{\|\mathbf{s}_{n+1}\|},
\]
从而
\[
\mathbf{s}_{n+1} = \mathbf{s}^{\text{trial}} - 2G\Delta\gamma\sqrt{\frac{3}{2}} \frac{\mathbf{s}_{n+1}}{\|\mathbf{s}_{n+1}\|},
\]
这表明 \(\mathbf{s}_{n+1}\) 和 \(\mathbf{s}^{\text{trial}}\) 同向，因此
\[
\frac{\mathbf{s}_{n+1}}{\|\mathbf{s}_{n+1}\|}=\frac{\mathbf{s}^{\text{trial}}}{\|\mathbf{s}^{\text{trial}}\|},
\]
于是
\[
\mathbf{s}_{n+1} = \left(1-\frac{\sqrt{6}G\Delta\gamma}{\|\mathbf{s}^{\text{trial}}\|}\right)\mathbf{s}^{\text{trial}},
\]
将上式和 \(\bar{\varepsilon}_{n+1}^p = \bar{\varepsilon}_n^p + \Delta \gamma\) 代入到屈服面函数，得到关于 \(\Delta\gamma\) 的方程
\[
\sqrt{\frac{3}{2}}\left(\|\mathbf{s}^{\text{trial}}\|-\sqrt{6}G\Delta\gamma\right)-\sigma_{y}(\bar{\varepsilon}_n^p + \Delta \gamma)=0,
\]
即
(104)\[
\sqrt{\frac{3}{2}}\|\mathbf{s}^{\text{trial}}\|-3G\Delta\gamma-\sigma_{y}(\bar{\varepsilon}_n^p + \Delta \gamma)=0.
\]
对于线性硬化，如
\[
\sigma_{y}(\bar{\varepsilon}_n^p + \Delta \gamma) = \sigma_{0} + H(\bar{\varepsilon}_n^p + \Delta \gamma),
\]
方程 (104) 是线性的；对于非线性硬化，方程是非线性的，通常使用 Newton-Rapshon 求解
\(\boldsymbol{\sigma}_{n+1} = \boldsymbol{\sigma}_{n+1}(\Delta\gamma)\)
\[\begin{split}
\begin{equation}
\begin{aligned}
\boldsymbol{\sigma}_{n+1} &= \boldsymbol{\sigma}^{\text{trial}} - 2G\Delta\gamma\sqrt{\frac{3}{2}} \frac{\mathbf{s}_{n+1}}{\|\mathbf{s}_{n+1}\|}\\
&=\mathbb{C}^e : \boldsymbol{\varepsilon}^{e, \text{trial}} - \sqrt{6}G\Delta\gamma\frac{\mathbf{s}^{\text{trial}}}{\|\mathbf{s}^{\text{trial}}\|}
\end{aligned}
\end{equation}
\end{split}\]
由于
\[
\boldsymbol{\sigma}^{\text{trial}} = \mathbb{C}^{e}:\boldsymbol{\varepsilon}^{e,\text{trial}}
\]
因此
\[
\begin{equation}
\begin{aligned}
\text{dev}(\boldsymbol{\sigma}^{\text{trial}}) + \frac{1}{3}\text{tr}(\boldsymbol{\sigma}^{\text{trial}})\mathbf{I} &= \mathbb{C}^{e}:(\text{dev}(\boldsymbol{\varepsilon}^{e,\text{trial}}) + \frac{1}{3}\text{tr}(\boldsymbol{\varepsilon}^{e,\text{trial}})\mathbf{I})
\end{aligned}
\end{equation}
\]
又
\[
\begin{equation}
\begin{aligned}
\mathbb{C}^{e}:\text{tr}(\boldsymbol{\varepsilon}^{e,\text{trial}})\mathbf{I} = \text{tr}(\boldsymbol{\sigma}^{\text{trial}})\mathbf{I}
\end{aligned}
\end{equation}
\]
故
\[\begin{split}
\begin{equation}
\begin{aligned}
\mathbf{s}^{\text{trial}}&=\text{dev}(\boldsymbol{\sigma}^{\text{trial}}) = \mathbb{C}^{e}:\text{dev}(\boldsymbol{\varepsilon}^{e,\text{trial}})\\
&=2G\text{dev}(\boldsymbol{\varepsilon}^{e,\text{trial}}) = 2G\mathbf{I}_{d}:\boldsymbol{\varepsilon}^{e,\text{trial}}
\end{aligned}
\end{equation}
\end{split}\]
其中，\(\mathbf{I}_{d}\) 是因此偏张量投影张量，于是
\[
\begin{equation}
\begin{aligned}
\boldsymbol{\sigma}_{n+1}=\left(\mathbb{C}^e - \frac{2\sqrt{6}G^{2}\Delta\gamma}{\|\mathbf{s}^{\text{trial}}\|}\mathbf{I}_{d}\right) : \boldsymbol{\varepsilon}^{e, \text{trial}}
\end{aligned}
\end{equation}
\]
由于 \(\|\mathbf{s}^{\text{trial}}\| = \sqrt{2J_{2}(\mathbf{s}^{\text{trial}})} = \sqrt{\frac{2}{3}}\cdot\sqrt{3J_{2}(\mathbf{s}^{\text{trial}})}\)，代入上式最终得到
\[
\begin{equation}
\begin{aligned}
\boldsymbol{\sigma}_{n+1}=\left(\mathbb{C}^e - \frac{6G^{2}\Delta\gamma}{\sqrt{3J_{2}(\mathbf{s}^{\text{trial}})}}\mathbf{I}_{d}\right) : \boldsymbol{\varepsilon}^{e, \text{trial}}
\end{aligned}
\end{equation}
\]
其中，\(\sqrt{3J_{2}(\mathbf{s}^{\text{trial}})}\) 是试验等效应力
