SpECTRE
v2025.01.30
|
A GMRES solver for nonsymmetric linear systems of equations \(Ax=b\). More...
#include <Gmres.hpp>
Public Types | |
using | fields_tag = FieldsTag |
using | options_group = OptionsGroup |
using | source_tag = SourceTag |
using | operand_tag = std::conditional_t< Preconditioned, db::add_tag_prefix< LinearSolver::Tags::Preconditioned, db::add_tag_prefix< LinearSolver::Tags::Operand, fields_tag > >, db::add_tag_prefix< LinearSolver::Tags::Operand, fields_tag > > |
Apply the linear operator to this tag in each iteration. | |
using | preconditioner_source_tag = db::add_tag_prefix< LinearSolver::Tags::Operand, fields_tag > |
Invoke a linear solver on the operand_tag sourced by the preconditioner_source_tag before applying the operator in each step. | |
using | component_list = tmpl::list< detail::ResidualMonitor< Metavariables, FieldsTag, OptionsGroup > > |
The parallel components used by the GMRES linear solver. | |
using | initialize_element = detail::InitializeElement< FieldsTag, OptionsGroup, Preconditioned > |
using | register_element = tmpl::list<> |
using | amr_projectors = initialize_element |
template<typename ApplyOperatorActions , typename ObserveActions = tmpl::list<>, typename Label = OptionsGroup> | |
using | solve = tmpl::list< detail::PrepareSolve< FieldsTag, OptionsGroup, Preconditioned, Label, SourceTag, ArraySectionIdTag >, ObserveActions, detail::NormalizeInitialOperand< FieldsTag, OptionsGroup, Preconditioned, Label, ArraySectionIdTag >, detail::PrepareStep< FieldsTag, OptionsGroup, Preconditioned, Label, ArraySectionIdTag >, ApplyOperatorActions, detail::PerformStep< FieldsTag, OptionsGroup, Preconditioned, Label, ArraySectionIdTag >, detail::OrthogonalizeOperand< FieldsTag, OptionsGroup, Preconditioned, Label, ArraySectionIdTag >, detail::NormalizeOperandAndUpdateField< FieldsTag, OptionsGroup, Preconditioned, Label, ArraySectionIdTag >, ObserveActions, detail::CompleteStep< FieldsTag, OptionsGroup, Preconditioned, Label, ArraySectionIdTag > > |
Static Public Attributes | |
static constexpr bool | preconditioned = Preconditioned |
A GMRES solver for nonsymmetric linear systems of equations \(Ax=b\).
The only operation we need to supply to the algorithm is the result of the operation \(A(p)\) (see Linear Solver). Each step of the algorithm expects that \(A(q)\) is computed and stored in the DataBox as db::add_tag_prefix<LinearSolver::Tags::OperatorAppliedTo, operand_tag>
. To perform a solve, add the solve
action list to an array parallel component. Pass the actions that compute \(A(q)\) as the first template parameter to solve
. The second template parameter allows you to pass any further actions you wish to run in each step of the algorithm (such as observations). If you add the solve
action list multiple times, use the third template parameter to label each solve with a different type.
This linear solver supports preconditioning. Enable preconditioning by setting the Preconditioned
template parameter to true
. If you do, run a preconditioner (e.g. another parallel linear solver) in each step. The preconditioner should approximately solve the linear problem \(A(q)=b\) where \(q\) is the operand_tag
and \(b\) is the preconditioner_source_tag
. Make sure the tag db::add_tag_prefix<LinearSolver::Tags::OperatorAppliedTo, operand_tag>
is updated with the preconditioned result in each step of the algorithm, i.e. that it is \(A(q)\) where \(q\) is the preconditioner's approximate solution to \(A(q)=b\). The preconditioner always begins at an initial guess of zero. It does not need to compute the operator applied to the initial guess, since it's zero as well due to the linearity of the operator.
Note that the operand \(q\) for which \(A(q)\) needs to be computed is not the field \(x\) we are solving for but db::add_tag_prefix<LinearSolver::Tags::Operand, FieldsTag>
. This field is initially set to the residual \(q_0 = b - A(x_0)\) where \(x_0\) is the initial value of the FieldsTag
.
When the algorithm step is performed after the operator action \(A(q)\) has been computed and stored in the DataBox, the GMRES algorithm implemented here will converge the field \(x\) towards the solution and update the operand \(q\) in the process. This requires reductions over all elements that are received by a ResidualMonitor
singleton parallel component, processed, and then broadcast back to all elements. Since the reductions are performed to find a vector that is orthogonal to those used in previous steps, the number of reductions increases linearly with iterations. No restarting mechanism is currently implemented. The actions are implemented in the gmres::detail
namespace and constitute the full algorithm in the following order:
PerformStep
(on elements): Start an Arnoldi orthogonalization by computing the inner product between \(A(q)\) and the first of the previously determined set of orthogonal vectors.StoreOrthogonalization
(on ResidualMonitor
): Keep track of the computed inner product in a Hessenberg matrix, then broadcast.OrthogonalizeOperand
(on elements): Proceed with the Arnoldi orthogonalization by computing inner products and reducing to StoreOrthogonalization
on the ResidualMonitor
until the new orthogonal vector is constructed. Then compute its magnitude and reduce.StoreOrthogonalization
(on ResidualMonitor
): Perform a QR decomposition of the Hessenberg matrix to produce a residual vector. Broadcast to NormalizeOperandAndUpdateField
along with a termination flag if the Convergence::Tags::Criteria
are met.NormalizeOperandAndUpdateField
(on elements): Set the operand \(q\) as the new orthogonal vector and normalize. Use the residual vector and the set of orthogonal vectors to determine the solution \(x\).Parallel::Section
). Set the ArraySectionIdTag
template parameter to restrict the solver to elements in that section. Only a single section must be associated with the ArraySectionIdTag
. The default is void
, which means running over all elements in the array. Note that the actions in the ApplyOperatorActions
list passed to solve
will not be restricted to run only on section elements, so all elements in the array may participate in preconditioning (see LinearSolver::multigrid::Multigrid).