implicitMultiUpdater (1d, 2d, 3d)

The implicitMultiUpdater provides Jacobian Free Newton-Krylov methods with a range of non-linear and linear solvers and preconditioning strategies for solving elliptic or implicit hyperbolic problems. The implicitMultiUpdater casts these problems in residual form:

\[\notag \begin{align} \mathcal{R} \left( \mathbf{q} \right) = 0 \end{align}\]

This formulation can then be used to solve linear problems, such as Poisson’s equation:

\[\notag \begin{align} \mathcal{R} \left( \mathbf{q} \right) = \nabla^2 \phi - \sum_{\mathrm{species}} Q_{\mathrm{species}} \end{align}\]

or, non linear problems such as a backward Euler discretization of a non-linear hyperbolic equation:

\[\notag \begin{align} \mathcal{R} \left( \mathbf{q} \right) = q^{n+1} - q^{n} + \Delta t \left\{ \nabla\cdot\left[ \mathcal{F} \left( \mathbf{w}^{n+1} \right) \right] - \mathcal{S} \left( \mathbf{w}^{n+1} \right) \right\} \end{align}\]

The operations performed by the implictMultiUpdater are specified using the UpdateStep and UpdateSequence pattern used for the main USim input file. Note that only the loop attribute for the UpdateSequence is used by the implicitMultiUpdater. If the implictMultiUpdater is being used to solve a system that includes a time discretization (e.g. the backward Euler example above), then the attribute “operation = integrate” must be specified in the final UpdateStep in the loop in order for USim to perform the time integration. The UpdateStep attribute “operation = operate” is not compatible with the implicitMultiUpdater.

The time integration scheme used by the multiUpdater is specified by the use of a Time Integrator block.

Preconditioning of the implictMultiUpdater can be specified by the addition of a Preconditioner block. one preconditioner can be specified.

Data

in (string vector, required)

The implicitMultiUpdater accepts exactly one input variable:

Vector of Unknowns (nodalArray, n-components, required) The vector of unknowns, \(\mathbf{q}\) at time \(t^{n}\). Must have the attribute useEpetraVector = true defined in the associated dataStruct block.

out (string vector, required)

The implicitMultiUpdater accepts exactly one output variable:

Vector of Unknowns (nodalArray, n-components, required) The vector of unknowns, \(\mathbf{q}\) at time \(t^{n+1}\). Must have the attribute useEpetraVector = true defined in the associated dataStruct block and have the same number of components as the input vector of unknows.

residual (string vector, optional)
A string vector that specifies exactly one nodalArray to store the residual, \(\mathcal{R} \left( \mathbf{q} \right)\) at the end of the solve. The specified nodalArray must have the attribute useEpetraVector = true defined in the associated dataStruct block.
solverPerf (string vector, optional)

A string vector that specifies exactly one dynVector with 6 componets to store the solver performance. Data is appended to this data structure at each Newton iteration. The components correspond to:

  1. Number of non-linear iterations.
  2. The residual norm.
  3. The number of linear iterations.
  4. The achieved tolerance.
  5. Time to solve the linear system.
  6. Time to construct the preconditioner.

Parameters

The implicitMultiUpdater updater accepts the parameters below, in addition to those required by Updater:

maxNonlinearIterations (integer, required)
Maximum number of outer Newton steps the solver will take before returning an error code.
maxLinearIterations (integer, required)
The maximum number of inner linear (Krylov) iterations to take at each Newton iteration.
numItersToStagnation (integer, optional)
The number of iterations that the stagnation condition (see stagnationThreshold) is allowed to be violated for before the outer Newton solve exits. Defaults to \(\mathrm{maxNonlinearIterations}/4\)
linearTolerance (float, required)
The tolerance required to achieve convergence in the inner linear (Krylov) iterations for each Newton step.
relativeResidual (float, optional)
The outer Newton problem converges when \(\mathrm{currentResidual} \le \mathrm{relativeResidual} \times \mathrm{initialResidual}\). Required if convergenceCriteria = [relativeResidual] is specified.
absoluteResidual (float, optional)
The outer Newton problem converges when \(\mathrm{currentResidual} \le \mathrm{absoluteResidual}\). Required if convergenceCriteria = [absoluteResidual] is specified.
stagnationThreshold (float, optional)
Causes the outer Newton solve to exit if \(\mathrm{currentResidual} > \mathrm{stagnationThreshold} \times \mathrm{previousResidual}\) for numItersToStagnation iterations. Defaults to 1.0.
convergenceTest (string vector, optional)

List of convergence tests to determine when the outer Newton solve is converged. Options include:

relativeResidual The outer Newton problem converges when \(\mathrm{currentResidual} \le \mathrm{relativeResidual} \times \mathrm{initialResidual}\).

absoluteResidual The outer Newton problem converges when \(\mathrm{currentResidual} \le \mathrm{absoluteResidual}\)

normSolutionUpdate The outer Newton problem converges when the norm of the change in the solution vector falls below 1.0e-3.

rmsSolutionUpdate The outer Newton problem converges when the weighted root mean square norm fo the solution update satisfies

\[\notag \begin{align}\sqrt{ \frac{1}{N} \sum_{i=1}^N \left( \frac {(q^k_i-q^{k-1}_i)}{ 10^{-2} |q^{k-1}_i| + 10^{-8}} \right) ^2 } \le 1.0 \end{align}\]

finiteResidual Causes the Newton solve to exit if NaN is encountered.

residualStagnation Causes the Newton solve to exit if \(\mathrm{currentResidual} > \mathrm{stagnationThreshold} \times \mathrm{previousResidual}\) for numItersToStagnation iterations.

Sub-Blocks

TimeIntegrator
Time Integrator Currently only one time integration scheme can be specified.
Preconditioner
Preconditioner Currently, only one preconditioner can be specified.
UpdateSequence
UpdateSequence This block is used to set the sequence of update steps
UpdateStep
UpdateStep The steps used in the update

Examples

The code block below demonstrates the implicitMultiUpdater for solving a three-dimensional Poisson problem using an algebraic multigrid preconditioner:

<Updater computeValues>
   kind = implicitMultiUpdater3d
   onGrid = domain
   in = [phi]
   out = [phiNew]
   residual = [f]
   solverPerf = [solverPerformance]
   maxNonlinearIterations = 4
   nonlinearTolerance = 1.e-12
   maxLinearIterations = 128
   linearTolerance = 1.e-14
   computePreconditioningMatrix = 1
   preconditioner = ML
   newtonForceMethod = Constant
   stencilUpdater = [computeNablaPhi]
   writePreconditioningMatrixToFile = 0

   <TimeIntegrator implicitUpdater>
     kind = implicit3d
     ongrid = domain
     scheme = none
   </TimeIntegrator>

   <Preconditioner myPreconditioner>
     # None/ML/AztecOO/Ifpack/New Ifpack
     preconditioner=ML
     # if 0, use a FD preconditioner
     computePreconditioningMatrix = 1
     # write out the preconditioning matrix at startup
     writePreconditioningMatrixToFile = 0
     # maximum age of preconditioner in outer Newton steps (needs a better name)
     linearMaxPrecAge = 10
     # rebuild, reuse or recompute preconditioner (needs a better name)
     linearReusePolicy = Reuse
     mlStrategy=classicSA # SA/DD
     stencilUpdater = [computeNablaPhi]
     testPreconditioner = false
     #mlDampingFactor=0.0675
     #mlSmoother = BSGS-A

     # Block to allow arbitray parameters to be passed to the preconditioner XML list
     #<ParameterList ml>
     #  increasing or decreasing = "increasing"
     #</ParameterList>
   </Preconditioner>

   <UpdateSequence sequence>
     loop = [boundaries,hyper]
   </UpdateSequence>

   <UpdateStep hyper>
      updaters = [computeNablaPhi,poisson,copyBc]
      syncVars = [phiNew]
   </UpdateStep>

   <UpdateStep boundaries>
      updaters = [dirchletBc]
      syncVars = [phi]
   </UpdateStep>

 </Updater>

The code block demonstrates the implicitMultiUpdater for solving a 2 or 3 dimensional compressible flow problem:

<Updater computeValues>
   kind = implicitMultiUpdater$NDIM$d
   onGrid = domain
   inpIndices = [0]
   in = [q]
   out = [qNew]
   residual = [f]
   solverPerf = [solverPerformance]
   maxNonlinearIterations = 100
   nonlinearTolerance = 1.e-6
   maxLinearIterations = 128
   linearTolerance = 1.e-4
   preconditioner = ML
   #aztecPreconditioner = ilut
   stencilUpdater = [hyper]
   newtonForceMethod = Constant
   computePreconditioningMatrix = 0

   <TimeIntegrator implicitUpdater>
     kind = implicit$NDIM$d
     ongrid = domain
     scheme = theta
     theta = 0.5 #Crank-Nicholson
     noInitialGuess = true
   </TimeIntegrator>

   <Preconditioner myPreconditioner>
     kind=autoPreconditioner$NDIM$d
     # None/ML/AztecOO/Ifpack/New Ifpack
     preconditioner=ML
     # if 0, use a FD preconditioner
     computePreconditioningMatrix = 1
     # write out the preconditioning matrix at startup
     writePreconditioningMatrixToFile = 0
     testPreconditioner = false
     # maximum age of preconditioner in outer Newton steps (needs a better name)
     linearMaxPrecAge = 10
     # rebuild, reuse or recompute preconditioner (needs a better name)
     linearReusePolicy = Reuse
     mlStrategy=DD # SA/DD
     #mlSmoother="symmetric Gauss Seidel"
     stencilUpdater = [hyper]
     #mlDampingFactor=0.125
     mlSmoother = BSGS-E

     # Block to allow arbitray parameters to be passed into the preconditioner XML list
     #<ParameterList ml>
     #  increasing or decreasing = "increasing"
     #</ParameterList>
   </Preconditioner>

   <UpdateSequence sequence>
     loop = [boundaries,hyper]
   </UpdateSequence>

   <UpdateStep hyper>
      operation = "integrate"
      updaters = [hyper,periodicFlux,copyQnewToFlux]
      syncVars = [qNew]
   </UpdateStep>

   <UpdateStep boundaries>
      updaters = [periodicQ]
      syncVars = [q]
   </UpdateStep>

 </Updater>