I'm solving coupled PDEs using finite difference method: Incompressible Navier-Stokes and the divergence-free induction equation (Maxwell's equations) with non-uniform electrical conductivity. The dimensionless equations are:

$1) \partial_t \mathbf{u} + \nabla \bullet (\mathbf{u}\mathbf{u}) = -\nabla p + Re^{-1}\nabla^2 \mathbf{u} + N \mathbf{j} \times \mathbf{B}, \qquad \mathbf{j} = Re_m^{-1} \nabla \times \mathbf{B}$

$2) \partial_t \mathbf{B} + Re_m^{-1} \nabla \times (\bar{\sigma}^{-1} \nabla \times \mathbf{B}) = \nabla \times (\mathbf{u} \times \mathbf{B})$

$3) \nabla \bullet \mathbf{u} = 0, \qquad \nabla \bullet \mathbf{B} = 0$

Where $Re,N,Re_m$ are the Reynolds, Interaction and magnetic Reynolds numbers. Also, $\bar{\sigma}$ is the dimensionless electrical conductivity.

The induction equation becomes stiff / ill-conditioned as $\bar{\sigma}$ approaches zero. In my case, $\bar{\sigma}$ is very small and, as a result, I resort to treating $Re_m^{-1} \nabla \times (\bar{\sigma}^{-1} \nabla \times \mathbf{B})$ implicitly in time using diagonally Preconditioned Conjugate Gradient Method (PCG). I use the limiting stopping criteria among 3 choices:

  • $N_{\text{max-iterations}}$
  • $\frac{||Ax_k-b||}{||Ax_0-b||} \le \varepsilon_{\text{relative}}$
  • $||Ax_k-b|| \le \varepsilon_{\text{absolute}}$

Here, $A=I + \Delta t Re_m^{-1} \nabla \times (\bar{\sigma}^{-1} \nabla \times)$. $b = \Delta t (\mathbf{B}^n + \nabla \times (\mathbf{u}^n \times \mathbf{B}^n))$. $x_k$ is the $k$th estimated solution of $\mathbf{B}^{n+1}$ using PCG. Finally, $N,\varepsilon$ are the maximum allowable PCG iterations and tolerances (relative to initial guess and absolute respectively).

Experimenting with stopping criteria and different time steps ($\Delta t$), it seems that there is a non-trivial trade-off in picking $\Delta t$ and stopping criteria to optimize performance, which I'll define as

$\mathcal P = \frac{\text{convective time per step}}{\text{wall clock time per step}} = \frac{\text{convective time}}{\text{wall clock time}}$.

Ultimately, I'm trying to reach a steady solution, however, the coupled equations are non-linear, stiff and can (depending on $Re,N,Re_m$) have seemingly long time scales. If I choose $\Delta t$ too large, or perform too few iterations of PCG on Eq. 2), the solution diverges.


1) Are there any methods to optimize the choice of stopping criteria or $\Delta t$ in order to optimize $\mathcal P$?

2) Are there any other tricks I can use to reach steady state faster?


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.