Solver parameters#
Subsection: Solver parameters#
Parameter name: Composition solver tolerance#
Default value: 1e-12
Pattern: [Double 0…1 (inclusive)]
Documentation: The relative tolerance up to which the linear system for the composition system gets solved. See ‘Stokes solver parameters/Linear solver tolerance’ for more details.
Parameter name: Temperature solver tolerance#
Default value: 1e-12
Pattern: [Double 0…1 (inclusive)]
Documentation: The relative tolerance up to which the linear system for the temperature system gets solved. See ‘Stokes solver parameters/Linear solver tolerance’ for more details.
Subsection: Solver parameters / AMG parameters#
Parameter name: AMG aggregation threshold#
Default value: 0.001
Pattern: [Double 0…1 (inclusive)]
Documentation: This threshold tells the AMG setup how the coarsening should be performed. In the AMG used by ML, all points that strongly couple with the tentative coarse-level point form one aggregate. The term strong coupling is controlled by the variable aggregation_threshold, meaning that all elements that are not smaller than aggregation_threshold times the diagonal element do couple strongly. The default is strongly recommended. There are indications that for the Newton solver a different value might be better. For extensive benchmarking of various settings of the AMG parameters in this section for the Stokes problem and others, see https://github.com/geodynamics/aspect/pull/234.
Parameter name: AMG output details#
Default value: false
Pattern: [Bool]
Documentation: Turns on extra information on the AMG solver. Note that this will generate much more output.
Parameter name: AMG smoother sweeps#
Default value: 2
Pattern: [Integer range 0…2147483647 (inclusive)]
Documentation: Determines how many sweeps of the smoother should be performed. When the flag elliptic is set to true, (which is true for ASPECT), the polynomial degree of the Chebyshev smoother is set to this value. The term sweeps refers to the number of matrix-vector products performed in the Chebyshev case. In the non-elliptic case, this parameter sets the number of SSOR relaxation sweeps for post-smoothing to be performed. The default is strongly recommended. There are indications that for the Newton solver a different value might be better. For extensive benchmarking of various settings of the AMG parameters in this section for the Stokes problem and others, see https://github.com/geodynamics/aspect/pull/234.
Parameter name: AMG smoother type#
Default value: Chebyshev
Pattern: [Selection Chebyshev|symmetric Gauss-Seidel ]
Documentation: This parameter sets the type of smoother for the AMG. The default is strongly recommended for any normal runs with ASPECT. There are some indications that the symmetric Gauss-Seidel might be better and more stable for the Newton solver. For extensive benchmarking of various settings of the AMG parameters in this section for the Stokes problem and others, see https://github.com/geodynamics/aspect/pull/234.
Subsection: Solver parameters / Advection solver parameters#
Parameter name: GMRES solver restart length#
Default value: 50
Pattern: [Integer range 1…2147483647 (inclusive)]
Documentation: This is the number of iterations that define the GMRES solver restart length. Increasing this parameter makes the solver more robust and decreases the number of iterations. Be aware that increasing this number increases the memory usage of the advection solver, and makes individual iterations more expensive.
Subsection: Solver parameters / Diffusion solver parameters#
Parameter name: Diffusion length scale#
Default value: 1.e4
Pattern: [Double 0…MAX_DOUBLE (inclusive)]
Documentation: Set a length scale for the diffusion of advection fields if the “prescribed field with diffusion” method is selected for a field. More precisely, this length scale represents the square root of the product of diffusivity and time in the diffusion equation, and controls the distance over which features are diffused. Units: \si{\meter}.
Subsection: Solver parameters / Matrix Free#
Parameter name: Execute solver timings#
Default value: false
Pattern: [Bool]
Documentation: Executes different parts of the Stokes solver repeatedly and print timing information. This is for internal benchmarking purposes: It is useful if you want to see how the solver performs. Otherwise, you don’t want to enable this, since it adds additional computational cost to get the timing information.
Parameter name: Output details#
Default value: false
Pattern: [Bool]
Documentation: Turns on extra information for the matrix free GMG solver to be printed.
Subsection: Solver parameters / Newton solver parameters#
Parameter name: Max Newton line search iterations#
Default value: 5
Pattern: [Integer range 0…2147483647 (inclusive)]
Documentation: The maximum number of line search iterations allowed. If the criterion is not reached after this number of iterations, we apply the scaled increment even though it does not satisfy the necessary criteria and simply continue with the next Newton iteration.
Parameter name: Max pre-Newton nonlinear iterations#
Default value: 10
Pattern: [Integer range 0…2147483647 (inclusive)]
Documentation: If the ’Nonlinear Newton solver switch tolerance’ is reached before the maximal number of Picard iterations, then the solver switches to Newton solves anyway.
Parameter name: Maximum linear Stokes solver tolerance#
Default value: 1e-2
Pattern: [Double 0…1 (inclusive)]
Documentation: The linear Stokes solver tolerance is dynamically chosen for the Newton solver, based on the Eisenstat Walker (1994) paper (https://doi.org/10.1137/0917003), equation 2.2. Because this value can become larger than one, we limit this value by this parameter.
Parameter name: Nonlinear Newton solver switch tolerance#
Default value: 1e-5
Pattern: [Double 0…1 (inclusive)]
Documentation: A relative tolerance with respect to the residual of the first iteration, up to which the nonlinear Picard solver will iterate, before changing to the Newton solver.
Parameter name: SPD safety factor#
Default value: 0.9
Pattern: [Double 0…1 (inclusive)]
Documentation: When stabilizing the Newton matrix, we can encounter situations where the coefficient inside the elliptic (top-left) block becomes negative or zero. This coefficient has the form \(1+x\) where \(x\) can sometimes be smaller than \(-1\). In this case, the top-left block of the matrix is no longer positive definite, and both preconditioners and iterative solvers may fail. To prevent this, the stabilization computes an \(\alpha\) so that \(1+\alpha x\) is never negative. This \(\alpha\) is chosen as \(1\) if \(x\ge -1\), and \(\alpha=-\frac 1x\) otherwise. (Note that this always leads to \(0\le \alpha \le 1\).) On the other hand, we also want to stay away from \(1+\alpha x=0\), and so modify the choice of \(\alpha\) to be \(1\) if \(x\ge -c\), and \(\alpha=-\frac cx\) with a \(c\) between zero and one. This way, if \(c<1\), we are assured that \(1-\alpha x>c\), i.e., bounded away from zero.
Parameter name: Stabilization preconditioner#
Default value: SPD
Pattern: [Selection SPD|PD|symmetric|none ]
Documentation: This parameters allows for the stabilization of the preconditioner. If one derives the Newton method without any modifications, the matrix created for the preconditioning is not necessarily Symmetric Positive Definite. This is problematic (see [Fraters et al., 2019]). When ‘none’ is chosen, the preconditioner is not stabilized. The ‘symmetric’ parameters symmetrizes the matrix, and ‘PD’ makes the matrix Positive Definite. ‘SPD’ is the full stabilization, where the matrix is guaranteed Symmetric Positive Definite.
Parameter name: Stabilization velocity block#
Default value: SPD
Pattern: [Selection SPD|PD|symmetric|none ]
Documentation: This parameters allows for the stabilization of the velocity block. If one derives the Newton method without any modifications, the matrix created for the velocity block is not necessarily Symmetric Positive Definite. This is problematic (see [Fraters et al., 2019]). When ‘none’ is chosen, the velocity block is not stabilized. The ‘symmetric’ parameters symmetrizes the matrix, and ‘PD’ makes the matrix Positive Definite. ‘SPD’ is the full stabilization, where the matrix is guaranteed Symmetric Positive Definite.
Parameter name: Use Eisenstat Walker method for Picard iterations#
Default value: false
Pattern: [Bool]
Documentation: If set to true, the Picard iteration uses the Eisenstat Walker method to determine how accurately linear systems need to be solved. The Picard iteration is used, for example, in the first few iterations of the Newton method before the matrix is built including derivatives of the model, since the Picard iteration generally converges even from points where Newton’s method does not.
Once derivatives are used in a Newton method, ASPECT always uses the Eisenstat Walker method.
Parameter name: Use Newton failsafe#
Default value: false
Pattern: [Bool]
Documentation: When this parameter is true and the linear solver fails, we try again, but now with SPD stabilization for both the preconditioner and the velocity block. The SPD stabilization will remain active until the next timestep, when the default values are restored.
Parameter name: Use Newton residual scaling method#
Default value: false
Pattern: [Bool]
Documentation: This method allows to slowly introduce the derivatives based on the improvement of the residual. If set to false, the scaling factor for the Newton derivatives is set to one immediately when switching on the Newton solver. When this is set to true, the derivatives are slowly introduced by the following equation: \(\max(0.0, (1.0-(residual/switch\_initial\_residual)))\), where switch_initial_residual is the residual at the time when the Newton solver is switched on.
Subsection: Solver parameters / Operator splitting parameters#
Parameter name: Reaction time step#
Default value: 1000.0
Pattern: [Double 0…MAX_DOUBLE (inclusive)]
Documentation: Set a time step size for computing reactions of compositional fields and the temperature field in case operator splitting is used. This is only used when the parameter “Use operator splitting” is set to true. The reaction time step must be greater than 0. If you want to prescribe the reaction time step only as a relative value compared to the advection time step as opposed to as an absolute value, you should use the parameter “Reaction time steps per advection step” and set this parameter to the same (or larger) value as the “Maximum time step” (which is 5.69e+300 by default). Units: Years or seconds, depending on the “Use years in output instead of seconds” parameter.
Parameter name: Reaction time steps per advection step#
Default value: 0
Pattern: [Integer range 0…2147483647 (inclusive)]
Documentation: The number of reaction time steps done within one advection time step in case operator splitting is used. This is only used if the parameter “Use operator splitting” is set to true. If set to zero, this parameter is ignored. Otherwise, the reaction time step size is chosen according to this criterion and the “Reaction time step”, whichever yields the smaller time step. Units: none.
Subsection: Solver parameters / Stokes solver parameters#
Parameter name: GMRES solver restart length#
Default value: 50
Pattern: [Integer range 1…2147483647 (inclusive)]
Documentation: This is the number of iterations that define the GMRES solver restart length. Increasing this parameter helps with convergence issues arising from high localized viscosity jumps in the domain. Be aware that increasing this number increases the memory usage of the Stokes solver, and makes individual Stokes iterations more expensive.
Parameter name: IDR(s) parameter#
Default value: 2
Pattern: [Integer range 1…2147483647 (inclusive)]
Documentation: This is the sole parameter for the IDR(s) Krylov solver and will dictate the number of matrix-vector products and preconditioner applications per iteration (s+1) and the total number of temporary vectors required (5+3*s). For s=1, this method is analogous to BiCGStab. As s is increased this method is expected to converge to GMRES in terms of matrix-vector/preconditioner applications to solution.
Parameter name: Krylov method for cheap solver steps#
Default value: GMRES
Pattern: [Selection GMRES|IDR(s) ]
Documentation: This is the Krylov method used to solve the Stokes system. Both options, GMRES and IDR(s), solve non-symmetric, indefinite systems. GMRES guarantees the residual will be reduced in each iteration while IDR(s) has no such property. On the other hand, the vector storage requirement for GMRES is dependent on the restart length and can be quite restrictive (since, for the matrix-free GMG solver, memory is dominated by these vectors) whereas IDR(s) has a short term recurrence. Note that the IDR(s) Krylov method is not available for the AMG solver since it is not a flexible method, i.e., it cannot handle a preconditioner which may change in each iteration (the AMG-based preconditioner contains a CG solve in the pressure space which may have different number of iterations each step).
Parameter name: Linear solver A block tolerance#
Default value: 1e-2
Pattern: [Double 0…1 (inclusive)]
Documentation: A relative tolerance up to which the approximate inverse of the \(A\) block of the Stokes system is computed. This approximate \(A\) is used in the preconditioning used in the GMRES solver. The exact definition of this block preconditioner for the Stokes equation can be found in [Kronbichler et al., 2012].
Parameter name: Linear solver S block tolerance#
Default value: 1e-6
Pattern: [Double 0…1 (inclusive)]
Documentation: A relative tolerance up to which the approximate inverse of the \(S\) block (i.e., the Schur complement matrix \(S = BA^{-1}B^{T}\)) of the Stokes system is computed. This approximate inverse of the \(S\) block is used in the preconditioning used in the GMRES solver. The exact definition of this block preconditioner for the Stokes equation can be found in [Kronbichler et al., 2012].
Parameter name: Linear solver tolerance#
Default value: 1e-7
Pattern: [Double 0…1 (inclusive)]
Documentation: A relative tolerance up to which the linear Stokes systems in each time or nonlinear step should be solved. The absolute tolerance will then be \(\| M x_0 - F \| \cdot \text{tol}\), where \(x_0 = (0,p_0)\) is the initial guess of the pressure, \(M\) is the system matrix, \(F\) is the right-hand side, and tol is the parameter specified here. We include the initial guess of the pressure to remove the dependency of the tolerance on the static pressure. A given tolerance value of 1 would mean that a zero solution vector is an acceptable solution since in that case the norm of the residual of the linear system equals the norm of the right hand side. A given tolerance of 0 would mean that the linear system has to be solved exactly, since this is the only way to obtain a zero residual.
In practice, you should choose the value of this parameter to be so that if you make it smaller the results of your simulation do not change any more (qualitatively) whereas if you make it larger, they do. For most cases, the default value should be sufficient. In fact, a tolerance of 1e-4 might be accurate enough.
Parameter name: Maximum number of expensive Stokes solver steps#
Default value: 1000
Pattern: [Integer range 0…2147483647 (inclusive)]
Documentation: This sets the maximum number of iterations used in the expensive Stokes solver. If this value is set too low for the size of the problem, the Stokes solver will not converge and return an error message pointing out that the user didn’t allow a sufficiently large number of iterations for the iterative solver to converge.
Parameter name: Number of cheap Stokes solver steps#
Default value: 200
Pattern: [Integer range 0…2147483647 (inclusive)]
Documentation: As explained in the paper that describes ASPECT (Kronbichler, Heister, and Bangerth, 2012, see [Kronbichler et al., 2012]) we first try to solve the Stokes system in every time step using a GMRES iteration with a poor but cheap preconditioner. By default, we try whether we can converge the GMRES solver in 200 such iterations before deciding that we need a better preconditioner. This is sufficient for simple problems with variable viscosity and we never need the second phase with the more expensive preconditioner. On the other hand, for more complex problems, and in particular for problems with strongly nonlinear viscosity, the 200 cheap iterations don’t actually do very much good and one might skip this part right away. In that case, this parameter can be set to zero, i.e., we immediately start with the better but more expensive preconditioner.
Parameter name: Stokes solver type#
Default value: block AMG
Pattern: [Selection block AMG|direct solver|block GMG ]
Documentation: This is the type of solver used on the Stokes system. The block geometric multigrid solver currently has a limited implementation and therefore may trigger Asserts in the code when used. If this is the case, please switch to ’block AMG’. Additionally, the block GMG solver requires using material model averaging.
Parameter name: Use direct solver for Stokes system#
Default value: false
Pattern: [Bool]
Documentation: If set to true the linear system for the Stokes equation will be solved using Trilinos klu, otherwise an iterative Schur complement solver is used. The direct solver is only efficient for small problems.
Parameter name: Use full A block as preconditioner#
Default value: false
Pattern: [Bool]
Documentation: This parameter determines whether we use an simplified approximation of the \(A\) block as preconditioner for the Stokes solver, or the full \(A\) block. The simplified approximation only contains the terms that describe the coupling of identical components (plus boundary conditions) as described in [Kronbichler et al., 2012]. The full block is closer to the description in [Rudi et al., 2017].
There is no clear way to determine which preconditioner performs better. The default value (simplified approximation) requires more outer GMRES iterations, but is faster to apply in each iteration. The full block needs less assembly time (because the block is available anyway), converges in less GMRES iterations, but requires more time per iteration. There are also differences in the amount of memory consumption between the two approaches.
The default value should be good for relatively simple models, but in particular for very strong viscosity contrasts the full \(A\) block can be advantageous.