# Finite strain
## Constitutive theory
For finite strain poromechanics, by writing Clausius-Duhem inequality (dissipation inequality) we can derive constitutive equations for solid and fluid phases as {cite}`irwin2024porofinitestrain` (see also {eq}`poromechanics-colemannoll-lambda` and discussion thereof)

$$
  \bm{S}^\rs_{E(\rs)} = 2\frac{\partial(\rho_{0(\rs)}^\rs\psi^\rs)}{ \partial \bm{C}_\rs}, \quad p_\rf = (\rho^{\fR})^2\frac{\partial \psi^\rf}{\partial \rho^{\fR}},
$$ (constitutive-solid-fluid)

where $\bm{S}^\rs_{E(\rs)}$ is effective second Piola-Kirchoff stress tensor, $\psi^\rs$ is strain energy per unit mass for the solid skeleton which is defined in previous sections for neo-Hookean model.

Helmholtz free energy per unit mass for fluid phase $\psi^\rf$ is given by

$$
  \psi^\rf = \frac{1}{2}\phi^\rf \bulk_\rf \left(\log \rho^{\fR}\right)^2,
$$ (energy-fluid-phase)

which using second equation of {eq}`constitutive-solid-fluid`, we arrive at constitutive equation for the real mass density of the pore fluid as

$$
  \rho^{\fR} = \rho^{\fR}_0 \exp\left( \frac{p_\rf - p_{\rf,0}}{\bulk_\rf}\right),
$$ (constitutive-real-density-fluid)

where $\rho^{\fR}_0, p_{\rf,0}$ are the initial real mass density and pore fluid pressure at time $t = 0$.

The total stress for the mixture for neo-Hookean model {eq}`neo-hookean-energy` (dropping now $(\cdot)_\rs$ subscripts where appropriate) is

$$
  \begin{aligned}
    \bm{S} &= J \bm{F}^{-1} \bm{\sigma} \bm{F}^{-T} = \bm{S}^\rs_E - J p_\rf \bm{C}^{-1}, \\
    \bm{S}^\rs_E &= \firstlame_d J V' \bm{C}^{-1} + \secondlame_d \left( \bm{I} - \bm{C}^{-1} \right) = \frac{\firstlame_d}{2} \left( J^2 - 1 \right)\bm{C}^{-1} + \secondlame_d \left( \bm{I} - \bm{C}^{-1} \right),
  \end{aligned}
$$ (constitutive-hyper-poroelastic-initial)

and using {eq}`S-to-tau` we can write symmetric Kirchhoff stress tensor as

$$
  \begin{aligned}
    \bm{\tau} &= \bm{\tau}^\rs_E - J p_\rf \bm{I}, \\
    \bm{\tau}^\rs_E &= \firstlame_d J V' \bm{I} + \secondlame_d \left( \bm{b} - \bm{I} \right) = \frac{\firstlame_d}{2} \left( J^2 - 1 \right)\bm{I} + \secondlame_d \left( \bm{b} - \bm{I} \right).
  \end{aligned}
$$ (constitutive-hyper-poroelastic-current)

## Strong and weak formulations
To derive the strong form for finite strain poromechanics, we use the balance of momentum of the mixture {eq}`poromechanics-momentum-mixture-initial`$_1$ (under the assumption that $\ddot{\bm{u}} = \bm{0}$) and balance of mass of the mixture {eq}`balance-mass-simplified-initial` as

$$
  \begin{aligned}
    -\nabla_X \cdot \bm{P} - \rho_0 \bm g &= 0, \qquad \text{in $\Omega_0$} \\
    \frac{\phi^\rf}{\bulk_\rf}\frac{D p_\rf}{D t} + \nabla_x\cdot \dot{\bm{u}} + \nabla_x\cdot \dot{\bm w} + \frac{\dot{\bm w}}{\bulk_\rf}\cdot \nabla_x p_\rf &= 0, \qquad \text{in $\Omega$} \\
    \bm{u} &= \bar{\bm u}, \quad    \text{on $\partial \Omega^{D}_0$} \\
    p_\rf &= \bar{p}_\rf, \quad    \text{on $\partial \Omega^{D}$} \\
    \bm{P} \cdot \bm N &= \bar{\bm t}, \qquad \text{on $\partial \Omega^{N}_0$} \\
    \dot{\bm{w}} \cdot \bm n &= \bar{s}, \qquad \text{on $\partial \Omega^{N}$} \\
    \bm u &= \bm u_0, \qquad \text{in $\Omega_0$} \\
    p_\rf &= p_{\rf,0}, \qquad \text{in $\Omega$}
  \end{aligned}
$$ (poro-hyper-strong-form)

where we have dropped subscript $(\cdot)_\rs$ by convention, $\bm N$ and $\bm n$ are the unit normal vectors on the face in initial and current configurations, $_X$ and $_x$ in $\nabla_X, \nabla_x$ indicate that the gradient are calculated with respect to the initial/current configurations. Multiplying the strong form {eq}`poro-hyper-strong-form` by test functions $(\bm{v}, q)$ and  integrate by parts, we can obtain the weak form for finite strain hyperporoelasticity as: find $(\bm{u}, p_\rf) \in \mathcal{V} \times \mathcal{Q} \subset H^1 \left( \Omega_0 \right) \times H^1 \left( \Omega_0 \right)$ such that

$$
  \begin{aligned}
    r_u(\bm u, p_\rf) &\coloneqq \int_{\Omega_0}{\nabla_x \bm{v} \tcolon \bm{\tau}} \, dV
     - \int_{\Omega_0}{\bm{v} \cdot \rho_0 \bm{g}} \, dV
     - \int_{\partial \Omega_0}{\bm{v} \cdot \bar{\bm t}} \, dA = 0, \quad \forall \bm{v} \in \mathcal{V}, \\
    r_{p_\rf}(\bm u, p_\rf) &\coloneqq \int_{\Omega_0} q \cdot \left( \frac{\phi^\rf}{\bulk_\rf}\dot{p}_\rf + \nabla_x\cdot \dot{\bm{u}} - \frac{\varkappa}{\bulk_\rf}\nabla_x p_\rf \cdot \nabla_x p_\rf \right) J \, dV \\
                  &+\int_{\Omega_0} \nabla_x q \cdot \varkappa \nabla_x p_\rf \, J dV + \int_{\partial \Omega}{ q \, \bar{s}} \, da = 0, \quad \forall q\in \mathcal{Q}, \\
  \end{aligned}
$$ (poro-hyper-weak-form)

where we have used the push forward {eq}`push-forward` to write the first equation with respect to $_x$ and $\bm \tau$, and replaced $\dot{\bm w} = -\varkappa \nabla_x p_\rf$ by assuming $\ddot{\bm u} = \bm{0}$ and $\bm f^\rf = \bm{0}$ in {eq}`generalized-darcy-finite-up`, the latter of which implies that $\bm{g} = \bm{0}$. The definition for $\bm{\tau}$ for the standard neo-Hookean model for the drained mixture (solid skeleton) response is given in {eq}`constitutive-hyper-poroelastic-current`.

## Linearization
As explained in {ref}`Newton-linearization-neo-hookean`, we need the Jacobian form of the {eq}`poro-hyper-weak-form` which is:
find $(\diff \bm{u}, \diff p_\rf) \in \mathcal{V} \times \mathcal{Q}$ such that

$$
  \begin{aligned}
    &\int_{\Omega_0} \nabla_x \bm{v} \tcolon \left( \diff \bm{\tau} -\bm{\tau}(\nabla_x \diff\bm{u})^T \right) dV = -r_u, \quad \forall \bm{v} \in \mathcal{V}, \\
    &\int_{\Omega_0} q \, \left( \diff L J + L \diff J \right) dV + \int_{\Omega_0} \nabla_x q \cdot (\diff K - K\nabla_x\diff\bm{u})\, dV = -r_{p_\rf}, \quad \forall q \in \mathcal{Q}
  \end{aligned}
$$ (poro-hyper-jacobian-current)

where

$$
  \begin{aligned}
    \diff \bm{\tau} -\bm{\tau}(\nabla_x \diff\bm{u})^T &= \diff \bm{\tau}_E^\rs -\bm{\tau}_E^\rs(\nabla_x \diff\bm{u})^T - \left(\diff J p_\rf + J \diff p_\rf \right) \bm I + J p_\rf (\nabla_x \diff\bm{u})^T, \\
    &= (\nabla_x \diff\bm{u}) \bm{\tau} + \bm{F} \diff \bm{S}_E^\rs \bm{F}^T - \left(\diff J p_\rf + J \diff p_\rf \right) \bm I + 2 J p_\rf \diff\bm\epsilon, \\
    L &= \frac{\phi^\rf}{\bulk_\rf}\dot{p}_\rf + \nabla_x\cdot \dot{\bm{u}} - \frac{\varkappa}{\bulk_\rf}\nabla_x p \cdot \nabla_x p, \\
    \diff L &= \mathrm{shift_v} \left(\frac{\phi^\rf}{\bulk_\rf}\diff p_\rf + \diff \, \left(\nabla_x\cdot \bm{u} \right) \right) - \frac{\diff \varkappa}{\bulk_\rf}\nabla_x p_\rf \cdot \nabla_x p_\rf - 2 \frac{\varkappa}{\bulk_\rf}\nabla_x p_\rf \cdot \diff \, (\nabla_x p_\rf), \\
    K &= \varkappa \nabla_x p_\rf J\, ,\\
    \diff K &= \diff \varkappa \nabla_x p_\rf \, J + \varkappa \diff \, (\nabla_x p_\rf) \, J + \varkappa \nabla_x p_\rf \diff J,
  \end{aligned}
$$

with given $\diff \bm{\tau}_E^\rs -\bm{\tau}_E^\rs(\nabla_x \diff\bm{u})^T$ in {eq}`cur_simp_Jac` for example with the standard neo-Hookean model and

$$
  \begin{aligned}
  \diff J &= J \, \trace (\nabla_x \diff\bm u) = J \, \trace \diff\bm\epsilon, \\
  \diff \varkappa &= \frac{\varkappa_{\circ}}{\mu_\rf \mathcal{F}(\phi^\rf_0)} \diff \mathcal{F}(\phi^\rf) = \frac{\varkappa_{\circ}}{\mu_\rf \mathcal{F}(\phi^\rf_0)} \frac{\partial \mathcal{F}(\phi^\rf)}{\partial \phi^\rf} \diff \phi^\rf = \varkappa\frac{\partial\mathcal{F}(\phi^\rf)}{\partial\phi^\rf} \mathcal{F}(\phi^\rf)^{-1} \diff \phi^\rf, \\
  \diff \phi^\rf &= \diff \left(1 - \phi^\rs \right) = \diff \left(1 - \frac{\phi^\rs_0}{J} \right) = \frac{\phi^\rs_0}{J^2} \diff J = \phi^\rs \trace \diff\bm\epsilon, \\
  \diff \, (\nabla_x p_\rf) &= \diff \, \left( \bm{F}^{-T} \nabla_X p_\rf \right) = \nabla_x \diff p_\rf - (\nabla_x \diff\bm{u})^T \nabla_x p_\rf, \\
  \diff \, \left(\nabla_x \cdot \bm{u} \right) &= \diff \, \trace \left( \nabla_x \bm{u}\right) = \diff \, \trace \left( \nabla_X \bm{u} \bm{F}^{-1}\right) = \trace \left( \nabla_x \diff \bm{u} - \nabla_x \bm{u} \nabla_x \diff \bm{u}\right), \\
  &= \nabla_x \cdot \diff\bm{u} - \nabla_x \bm{u} \tcolon (\nabla_x \diff\bm{u})^T,
  \end{aligned}
$$

where $\mathcal{F}(\phi^\rf)$ is defined, for example, in {eq}`kozeny-carman`.

## Generic weak formulation
For overview on generic weak form implementation for the poromechanics model, see {ref}`matrix-free-poromechanics`. Here we define

$$
  \begin{aligned}
    \bm f_0 &= -\rho_0 \bm g = \bm{0}, \quad \bm f_1 = \bm{\tau} = \bm{\tau}_E^\rs - J p_\rf \bm{I}, \\
    \bm g_0 &= L \, J = \left( \frac{\phi^\rf}{\bulk_\rf}\dot{p}_\rf + \nabla_x\cdot \dot{\bm{u}} - \frac{\varkappa}{\bulk_\rf}\nabla_x p_\rf \cdot \nabla_x p_\rf \right) J, \quad \bm g_1 = K = \varkappa \nabla_x p_\rf \, J, \\
		\diff \bm f_1 &= \diff \bm{\tau} -\bm{\tau}(\nabla_x \diff\bm{u})^T, \quad \diff \bm g_0 = \diff L J + L \diff J, \\
    \diff \bm g_1 &= \diff K - K\nabla_x\diff \bm{u} = \diff \varkappa \nabla_x p_\rf \, J + \varkappa \diff \, (\nabla_x p_\rf) \, J + \varkappa \nabla_x p_\rf \diff J - \varkappa \nabla_x p_\rf J\nabla_x \diff \bm{u} .
  \end{aligned}
$$

Note that, technically speaking, $\bm f_0$ needs to be linearized under the assmption of $\bm{g} \neq \bm{0}$ (and by extension $\bm{f}^\rf \neq \bm{0}$) in the case of evolving volume fractions and/or compressible pore fluid. However, at this time, we continue to assume the body force acting on the mixture to be zero.

For [dynamic example](example-dynamic) acceleration is **not zero** and we use the PETSc Timestepper (TS) object to manage timestepping. Specifically, Ratel uses the Generalized-Alpha method for second order systems (`TSALPHA2`) to solve the second order system of equations. More information on the PETSc `TSALPHA2` time stepper can be found in the [PETSc documentation TS](https://petsc.org/release/docs/manualpages/TS/index.html).

We add

$$
  \begin{aligned}
    &\int_{\Omega_0} \bm v \cdot \left(\rho^b_0 \, \ddot{\bm u} \right) \, dV = \int_{\Omega_0} \bm v \cdot \left(J \rho^b \, \ddot{\bm u} \right) \, dV, \\
    &\int_{\Omega_0} q \left(-\varkappa \, \ddot{\bm u} \cdot \nabla_x p_\rf\right) \frac{J \, \rho^{\fR}}{\bulk_\rf} dV + \int_{\Omega_0} \nabla_x q \cdot \left(\varkappa \rho^{\fR} J \, \ddot{\bm u}\right) \, dV,
  \end{aligned}
$$ (poro-hyper-dynamics-residual)

to the first and second equations of {eq}`poro-hyper-weak-form`, such that

$$
  \begin{aligned}
    \bm f_0 &\mathrel{+}= J\rho^b\ddot{\bm{u}}\, ,\\
    \bm g_0 &\mathrel{+}= -\frac{\varkappa}{\bulk_\rf}\nabla_x p_\rf\cdot(\rho^\fR\ddot{\bm{u}})J\,\\
    \bm g_1 &\mathrel{+}= \varkappa\rho^\fR\ddot{\bm{u}} J\, ,
  \end{aligned}
$$

and the jacobian is

$$
  \begin{aligned}
    &\int_{\Omega_0} \bm v \cdot \left(\diff J \rho^b \, \ddot{\bm u} + J \diff \rho^b \, \ddot{\bm u} + J \rho^b \, \mathrm{shift_a} \, \diff \bm u \right) \, dV, \\
    &\int_{\Omega_0} q \left[-\diff \varkappa  \, \ddot{\bm u} \cdot \nabla_x p_\rf - \varkappa \, \mathrm{shift_a} \diff \bm u \cdot \nabla_x p_\rf - \varkappa \, \ddot{\bm u} \cdot \diff \left(\nabla_x p_\rf \right) \right] \frac{J \, \rho^{\fR}}{\bulk_\rf} dV \\
    &\int_{\Omega_0} q \left(-\varkappa \, \ddot{\bm u} \cdot \nabla_x p_\rf\right) \frac{\left(\diff J \rho^{\fR} + J \, \diff \rho^{\fR}\right)}{\bulk_\rf} \, dV + \\
    &\int_{\Omega_0} \nabla_x q \cdot \left[ \left(\diff \varkappa \rho^{\fR} J + \varkappa \diff \rho^{\fR} J + \varkappa \rho^{\fR} \diff J \right) \ddot{\bm u} + \varkappa \rho^{\fR} J \, \mathrm{shift_a}\diff \bm u - \varkappa \rho^{\fR} J (\nabla_x \diff \bm u) \ddot{\bm u} \right] \, dV,
  \end{aligned}
$$ (poro-hyper-dynamics-jacobian)

such that,

$$
  \begin{aligned}
		\diff \bm f_0 &\mathrel{+}= \diff J \rho^b \, \ddot{\bm u} + J \diff \rho^b \, \ddot{\bm u} + J \rho^b \, \mathrm{shift_a} \, \diff \bm u\, ,\\
    \diff \bm g_0 &\mathrel{+}= \left[-\diff \varkappa  \, \ddot{\bm u} \cdot \nabla_x p_\rf - \varkappa \, \mathrm{shift_a} \diff \bm u \cdot \nabla_x p_\rf - \varkappa \, \ddot{\bm u} \cdot \diff \left(\nabla_x p_\rf \right) \right] \frac{J \, \rho^{\fR}}{\bulk_\rf}\\
    &\phantom{\mathrel{+}=} + \left(-\varkappa \, \ddot{\bm u} \cdot \nabla_x p_\rf\right) \frac{\left(\diff J \rho^{\fR} + J \, \diff \rho^{\fR}\right)}{\bulk_\rf}, \\
    \diff \bm g_1 &\mathrel{+}= \left(\diff \varkappa \rho^{\fR} J + \varkappa \diff \rho^{\fR} J + \varkappa \rho^{\fR} \diff J \right) \ddot{\bm u} + \varkappa \rho^{\fR} J \, \mathrm{shift_a}\diff \bm u - \varkappa \rho^{\fR} J (\nabla_x \diff \bm u) \ddot{\bm u}\, ,
  \end{aligned}
$$

where using {eq}`constitutive-real-density-fluid` we have

$$
  \begin{aligned}
    \diff \rho^b &= \left(\rho^{\fR}-\rho^{\sR} \right) \diff \phi^\rf + \phi^\rf \diff \rho^{\fR}, \\
    \diff \rho^{\fR} &= \frac{\rho^{\fR}_0 \, \diff p_\rf}{\bulk_\rf} \exp\left( \frac{p_\rf - p_{\rf,0}}{\bulk_\rf}\right) = \frac{\diff p_\rf}{\bulk_\rf}\rho^{\fR}.
  \end{aligned}
$$

Note that unlike linear case where the porosity and density are constant, in finite strain we have assumed that the real mass density of the pore fluid is a function of pore fluid pressure (real mass density of solid is assumed constant, see {eq}`real-solid-density`). Therefore, we use the mixture density in the current configuration $\rho^b$ in residual and jacobian to account for its variation with respect to pore fluid pressure.