The Accelerated Pseudo-Transient (APT) method consists in augmenting the right-hand-side of the target PDE with a pseudo-time derivative (where $\psi$ is the pseudo-time) of the primary variables. We then solve the resulting system of equations with an iterative method. The pseudo-time derivative is then gradually reduced, until the original PDE is solved and the changes in the primary variables are below a preset tolerance.

Heat diffusion

The APT heat-diffusion equation is:

$$\begin{array}{r}\tilde{\rho }\frac{\partial T}{\partial \psi }+\rho C_p\frac{\partial T}{\partial t}=\mathrm{\nabla }\cdot (\kappa \mathrm{\nabla }T)=-\mathrm{\nabla }q\end{array}$$

We use a second order APT scheme were continuation is also done on the flux, so that:

$$\begin{array}{r}\tilde{\theta }\frac{\partial q}{\partial \psi }+q=-\kappa \mathrm{\nabla }T\end{array}$$

Stokes equations

For example, the APT formulation of the Stokes equations yields:

$$\begin{array}{r}\tilde{\rho }\frac{\partial \mathit{u}}{\partial \psi }=\mathrm{\nabla }\cdot \mathit{\tau }-\mathrm{\nabla }p-\mathit{f}\end{array}$$$$\begin{array}{r}\frac{1}{\tilde{K}}\frac{\partial p}{\partial \psi }+\mathrm{\nabla }\cdot \mathit{v}=-\beta \frac{\partial p}{\partial t}+\alpha \frac{\partial T}{\partial t}\end{array}$$

Constitutive equations

A APT continuation is also done on the constitutive law:

$$\begin{array}{r}\frac{1}{2\tilde{G}}\frac{\partial \mathit{\tau }}{\partial \psi }+\frac{1}{2G}\frac{D\mathit{\tau }}{Dt}+\frac{\mathit{\tau }}{2\eta }=\dot{\mathit{\epsilon }}\end{array}$$

where the wide tile denotes the effective damping coefficients and $\psi$ is the pseudo-time step. These are defined as in Räss et al. (2022):

$$\begin{array}{r}\tilde{\rho }=Re\frac{\eta }{\tilde{V}L},\tilde{G}=\frac{\tilde{\rho }{\tilde{V}}^2}{r+2},\tilde{K}=r\tilde{G}\end{array}$$

and

$$\begin{array}{r}\tilde{V}=\sqrt{\frac{\tilde{K}+2\tilde{G}}{\tilde{\rho }}},r=\frac{\tilde{K}}{\tilde{G}},Re=\frac{\tilde{\rho }\tilde{V}L}{\eta }\end{array}$$

where the P-wave $\tilde{V}=V_p$ is the characteristic velocity scale for Stokes, and $Re$ is the Reynolds number.

Free-surface stabilization

For a vertically stratified 2D model, the free-surface stabilization used by the standard Stokes and DYREL solvers adds

$$\begin{array}{r}{R}_y^{FS}=R_y+\theta \mathrm{\Delta }t{v}_y\frac{\partial (\rho g_y)}{\partial y}\end{array}$$

to the vertical momentum residual. This velocity-dependent term contributes the local operator diagonal

$$\begin{array}{r}{c}_{FS}=-\theta \mathrm{\Delta }t\frac{\partial (\rho g_y)}{\partial y}.\end{array}$$

The standard APT update uses ${\eta }_{d\tau }/({\eta }_{\tau }+{\eta }_{d\tau }{c}_{FS})$ for the stabilized vertical pseudo-time step. DYREL adds the same $c_{FS}$ to its vertical diagonal preconditioner and Gershgorin eigenvalue bound. This keeps the iterative coefficients consistent with the stabilized residual. A non-positive standard APT denominator is rejected through the solver's non-finite residual check because it indicates that the physical timestep is too large for the local density inversion.

<!–

Physical parameters

Symbol	Parameter
$T$	Temperature
$q$	Flux
$\mathit{\tau }$	Deviatoric stress
$\dot{\mathit{\epsilon }}$	Deviatoric strain rate
$\mathit{u}$	Velocity
$\mathit{f}$	External forces
$P$	Pressure
$\eta$	Viscosity
$\rho$	Density
$\beta$	Compressibility
$G$	Shear modulus
$\alpha$	Thermal expansivity
$C_p$	Heat capacity
$\kappa$	Heat conductivity

–>

Pseudo-transient parameters

Symbol	Parameter
$\psi$	Pseudo time step
$\tilde{K}$	Pseudo bulk modulus
$\tilde{G}$	Pseudo shear modulus
$\tilde{V}$	Characteristic velocity scale
$\tilde{\rho }$	Pseudo density
$\tilde{\theta }$	Relaxation time
$Re$	Reynolds number