The discretized form of $\delta W_{\text{int}}$ in (3.3-3) may be written as $$\begin{equation} \delta W_{\text{int}}=\sum\limits _{e=1}^{n_{e}}\sum\limits _{k=1}^{n_{\text{int}}^{\left(e\right)}}W_{k}J_{\eta}\sum\limits _{a=1}^{m}\left[\begin{array}{ccc} \delta\mathbf{v}_{a} & \delta\tilde{p}_{a} & \delta\tilde{c}_{a}\end{array}\right]\cdot\left[\begin{array}{c} \mathbf{r}_{a}^{u}\\ r_{a}^{p}\\ r_{a}^{c} \end{array}\right]\,,\label{eq265} \end{equation}$$ where $n_{e}$ is the number of elements in $b$ , $n_{\text{int}}^{\left(e\right)}$ is the number of integration points in the $e-$ th element, $W_{k}$ is the quadrature weight associated with the $k-$ th integration point, and $J_{\eta}$ is the Jacobian of the transformation from the current spatial configuration to the parametric space of the element. In the above expression, $$\begin{equation} \begin{aligned}\mathbf{r}_{a}^{u} & =\boldsymbol{\sigma}\cdot\grad N_{a}\,,\\ r_{a}^{p} & =\mathbf{w}\cdot\grad N_{a}-N_{a}\frac{1}{J}\frac{\partial J}{\partial t}\,,\\ r_{a}^{c} & =\mathbf{j}\cdot\grad N_{a}-N_{a}\frac{1}{J}\frac{\partial}{\partial t}\left(J\varphi^{w}\tilde{\kappa}\tilde{c}\right)\,, \end{aligned} \label{eq266} \end{equation}$$ and it is understood that $J_{\eta}$ , $\mathbf{r}_{a}^{u}$ , $r_{a}^{p}$ and $r_{a}^{c}$ are evaluated at the parametric coordinates of the $k-$ th integration point. Since the parametric space is invariant, time derivatives are evaluated in a material frame. For example, the time derivative $D^{s}J\left(\mathbf{x},t\right)/Dt$ appearing in (3.3-3) becomes $\partial J\left(\eta_{k},t\right)/\partial t$ when evaluated at the parametric coordinates $\eta_{k}=\left(\eta_{k}^{1},\eta_{k}^{2},\eta_{k}^{3}\right)$ of the $k-$ th integration point.

Similarly, the discretized form of $D\delta W_{\text{int}}=D\delta W_{\text{int}}\left[\Delta\mathbf{u}\right]+D\delta W_{\text{int}}\left[\Delta\tilde{p}\right]+D\delta W_{\text{int}}\left[\Delta\tilde{c}\right]$ may be written as $$\begin{equation} D\delta W_{\text{int}}=\sum\limits _{e=1}^{n_{e}}\sum\limits _{k=1}^{n_{\text{int}}^{\left(e\right)}}W_{k}J_{\eta}\sum\limits _{a=1}^{m}\sum\limits _{b=1}^{m}\left[\begin{array}{ccc} \delta\mathbf{v}_{a} & \delta\tilde{p}_{a} & \delta\tilde{c}_{a}\end{array}\right]\cdot\left[\begin{array}{ccc} \mathbf{K}_{ab}^{uu} & \mathbf{k}_{ab}^{up} & \mathbf{k}_{ab}^{uc}\\ \mathbf{k}_{ab}^{pu} & k_{ab}^{pp} & k_{ab}^{pc}\\ \mathbf{k}_{ab}^{cu} & k_{ab}^{cp} & k_{ab}^{cc} \end{array}\right]\cdot\left[\begin{array}{c} \Delta\mathbf{u}_{b}\\ \Delta\tilde{p}_{b}\\ \Delta\tilde{c}_{b} \end{array}\right]\,,\label{eq267} \end{equation}$$ where the terms in the first column are the discretized form of the linearization along $\Delta\mathbf{u}$ : $$\begin{equation} \mathbf{K}_{ab}^{uu}=\grad N_{a}\cdot\boldsymbol{\mathcal{C}}\cdot\grad N_{b}+\left(\grad N_{a}\cdot\boldsymbol{\sigma}\cdot\grad N_{b}\right)\mathbf{I}\,,\label{eq268} \end{equation}$$ $$\begin{equation} \mathbf{k}_{ab}^{pu}=\left(\mathbf{w}_{b}^{u}\right)^{T}\cdot\grad N_{a}+N_{a}\mathbf{q}_{b}^{pu}\,,\label{eq269} \end{equation}$$ $$\begin{equation} \mathbf{k}_{ab}^{cu}=\left(\mathbf{j}_{b}^{u}\right)^{T}\cdot\grad N_{a}+N_{a}\mathbf{q}_{b}^{cu}\,,\label{eq270} \end{equation}$$ where $$\begin{equation} \begin{aligned}\mathbf{j}_{b}^{u} & =J\frac{\partial\tilde{\kappa}}{\partial J}\left[\mathbf{d}\cdot\left(-\varphi^{w}\grad\tilde{c}+\frac{\tilde{c}}{d_{0}}\mathbf{w}\right)\right]\otimes\grad N_{b}+\tilde{\kappa}\left({-\varphi^{w}\grad\tilde{c}+\frac{\tilde{c}}{d_{0}}\mathbf{w}}\right)\cdot\mathbf{d}\cdot\grad N_{b}\\ & +\tilde{\kappa}\left({-\varphi^{s}\left({\mathbf{d}\cdot\grad\tilde{c}}\right)\otimes\grad N_{b}+\frac{\tilde{c}}{d_{0}}\left[{2\left({\grad N_{b}\cdot\mathbf{w}}\right)\mathbf{d}-\left({\mathbf{d}\cdot\mathbf{w}}\right)\otimes\grad N_{b}}\right]}\right)+\tilde{\kappa}\frac{\tilde{c}}{d_{0}}\mathbf{d}\cdot\mathbf{w}_{b}^{u}\;, \end{aligned} \label{eq272} \end{equation}$$ $$\begin{equation} \mathbf{q}_{b}^{pu}=-\frac{1}{\Delta t}\grad N_{b},\label{eq273} \end{equation}$$ $$\begin{equation} \mathbf{q}_{b}^{cu}=\tilde{c}\frac{\partial\left({J\phi^{w}\tilde{\kappa}}\right)}{\partial J}\mathbf{q}_{b}^{pu}.\label{eq274} \end{equation}$$ The terms in the second column of the stiffness matrix in (3.3.3-4) are the discretized form of the linearization along $\Delta\tilde{p}$ : $$\begin{equation} \mathbf{k}_{ab}^{up}=-N_{b}\grad N_{a},\label{eq275} \end{equation}$$ $$\begin{equation} k_{ab}^{pp}=-\grad N_{a}\cdot\tilde{\mathbf{k}}\cdot\grad N_{b},\label{eq276} \end{equation}$$ $$\begin{equation} k_{ab}^{cp}=-\frac{\tilde{\kappa}\tilde{c}}{d_{0}}\grad N_{a}\cdot\mathbf{d}\cdot\tilde{\mathbf{k}}\cdot\grad N_{b}.\label{eq277} \end{equation}$$ The terms in the third column of the stiffness matrix in (3.3.3-4) are the discretized form of the linearization along $\Delta\tilde{c}$ : $$\begin{equation} \mathbf{k}_{ab}^{uc}=N_{b}\left(\boldsymbol{\sigma}_{c}^{\prime}\cdot\grad N_{a}-R\theta\frac{\partial\left(\Phi\tilde{\kappa}\tilde{c}\right)}{\partial\tilde{c}}\grad N_{a}\right)\,,\label{eq278} \end{equation}$$ $$\begin{equation} k_{ab}^{pc}=\grad N_{a}\cdot\mathbf{w}_{b}^{c}\,,\label{eq279} \end{equation}$$ $$\begin{equation} k_{ab}^{cc}=\grad N_{a}\cdot\mathbf{j}_{b}^{c}+N_{a}q_{b}^{c}\,,\label{eq280} \end{equation}$$ where $$\begin{equation} \begin{aligned}\mathbf{w}_{b}^{c} & =-N_{b}\tilde{\mathbf{k}}_{c}^{\prime}\cdot\left(\grad\tilde{p}+R\theta\frac{\tilde{\kappa}}{d_{0}}\mathbf{d}\cdot\grad\tilde{c}\right)\\ & -R\theta\tilde{\mathbf{k}}\cdot\left[N_{b}\left(\frac{\partial}{\partial\tilde{c}}\left(\frac{\tilde{\kappa}}{d_{0}}\right)\mathbf{d}+\frac{\tilde{\kappa}}{d_{0}}\mathbf{d}_{c}^{\prime}\right)\cdot\grad\tilde{c}+\frac{\tilde{\kappa}}{d_{0}}\mathbf{d}\cdot\grad N_{b}\right]\,, \end{aligned} \label{eq281} \end{equation}$$ $$\begin{equation} \mathbf{j}_{b}^{c}=N_{b}\left(\frac{\partial\tilde{\kappa}}{\partial\tilde{c}}\mathbf{d}+\tilde{\kappa}\mathbf{d}_{c}^{\prime}\right)\cdot\left(-\varphi^{w}\grad\tilde{c}+\frac{\tilde{c}}{d_{0}}\mathbf{w}\right)+\tilde{\kappa}\mathbf{d}\cdot\left(-\varphi^{w}\grad N_{b}+\frac{\tilde{c}}{d_{0}}\mathbf{w}_{b}^{c}\right)\,,\label{eq282} \end{equation}$$ $$\begin{equation} q_{b}^{c}=-N_{b}\frac{\phi^{w}}{\Delta t}\frac{\partial\left(\tilde{\kappa}\tilde{c}\right)}{\partial\tilde{c}}\,.\label{eq283} \end{equation}$$ The discretization of $\delta W_{\text{ext}}$ in (3.3-3) has the form $$\begin{equation} \delta W_{\text{ext}}=\sum\limits _{e=1}^{n_{e}}\sum\limits _{k=1}^{n_{\mbox{int}}^{\left(e\right)}}W_{k}J_{\eta}\sum\limits _{a=1}^{m}\left[\begin{array}{ccc} \delta\mathbf{v}_{a} & \delta\tilde{p}_{a} & \delta\tilde{c}_{a}\end{array}\right]\cdot\left[\begin{array}{c} N_{a}t_{n}\mathbf{n}\\ N_{a}w_{n}\\ N_{a}j_{n} \end{array}\right]\,,\label{eq284} \end{equation}$$ where $J_{\eta}=\left|\mathbf{g}_{1}\times\mathbf{g}_{2}\right|$ . The summation is performed over all surface elements on which these boundary conditions are prescribed. The discretization of $-D\delta W_{\text{ext}}$ has the form $$\begin{equation} -D\delta W_{\text{ext}}=\sum\limits _{e=1}^{n_{e}}\sum\limits _{k=1}^{n_{\text{int}}^{\left(e\right)}}W_{k}J_{\eta}\sum\limits _{a=1}^{m}\sum\limits _{b=1}^{m}\left[\begin{array}{ccc} \delta\mathbf{v}_{a} & \delta\tilde{p}_{a} & \delta\tilde{c}_{a}\end{array}\right]\cdot\left[\begin{array}{ccc} \mathbf{K}_{ab}^{uu} & \mathbf{0} & \mathbf{0}\\ \mathbf{k}_{ab}^{pu} & 0 & 0\\ \mathbf{k}_{ab}^{cu} & 0 & 0 \end{array}\right]\cdot\left[\begin{array}{c} \Delta\mathbf{u}_{b}\\ \Delta\tilde{p}_{b}\\ \Delta\tilde{c}_{b} \end{array}\right]\,,\label{eq285} \end{equation}$$ where $$\begin{equation} \begin{aligned}\mathbf{K}_{ab}^{uu} & =t_{n}N_{a}\boldsymbol{\mathcal{A}}\left\{ \frac{\partial N_{b}}{\partial\eta^{1}}\mathbf{g}_{2}-\frac{\partial N_{b}}{\partial\eta^{2}}\mathbf{g}_{1}\right\} \,,\\ \mathbf{k}_{ab}^{pu} & =-w_{n}N_{a}\left(\frac{\partial N_{b}}{\partial\eta^{1}}\mathbf{g}_{2}-\frac{\partial N_{b}}{\partial\eta^{2}}\mathbf{g}_{1}\right)\times\mathbf{n}\,,\\ \mathbf{k}_{ab}^{cu} & =-j_{n}N_{a}\left(\frac{\partial N_{b}}{\partial\eta^{1}}\mathbf{g}_{2}-\frac{\partial N_{b}}{\partial\eta^{2}}\mathbf{g}_{1}\right)\times\mathbf{n}\,. \end{aligned} \label{eq286} \end{equation}$$ In this expression, $\boldsymbol{\mathcal{A}}\left\{ \mathbf{v}\right\}$ is the antisymmetric tensor whose dual vector is $\mathbf{v}$ (such that $\boldsymbol{\mathcal{A}}\left\{ \mathbf{v}\right\} \cdot\mathbf{q}=\mathbf{v}\times\mathbf{q}$ for any vector $\mathbf{q})$ .