gko::LinOpFactory#

Abstract base for factories that produce a LinOp from a system matrix. The canonical workflow is

auto factory = gko::solver::Cg<double>::build()
                   .with_criteria(gko::stop::max_iters(100u))
                   .on(exec);   // returns a unique_ptr<Cg<double>::Factory>
auto solver = factory->generate(system_matrix);   // returns a unique_ptr<Cg<double>>
solver->apply(b, x);

with_* setters configure the factory parameters; .on(exec) ties the factory to an executor and produces the concrete factory, which is a LinOpFactory; finally generate(matrix) binds the factory to a concrete system and returns the operator. Solver / preconditioner / factorisation classes all expose this pattern through the GKO_ENABLE_LIN_OP_FACTORY macro, which builds on the EnableDefaultLinOpFactory mixin.

Generating with reuse#

When a sequence of matrices shares one sparsity pattern, a factory can keep what it computed for the first matrix. create_empty_reuse_data() creates an empty ReuseData, and generate_reuse(input, reuse_data) generates the operator. The first call is a full generation that fills reuse_data; later calls reuse it and recompute what depends on the values. The input must keep its size and sparsity pattern, and one ReuseData serves the whole sequence. The generated operator is independent of reuse_data. Factories without reuse support call generate().

auto factory = gko::solver::Multigrid::build()
                   .with_mg_level(gko::multigrid::Pgm<double, int>::build())
                   .with_criteria(gko::stop::max_iters(1u))
                   .on(exec);
auto reuse_data = factory->create_empty_reuse_data();
auto mg = factory->generate_reuse(A, *reuse_data);
// ... change the values of A, keeping its sparsity pattern ...
mg = factory->generate_reuse(A, *reuse_data);

Attention

The reuse interface (ReuseData, create_empty_reuse_data() and generate_reuse()) is experimental and may change.

class LinOpFactory #

Inherits from

  • public gko::AbstractFactory<LinOp, std::shared_ptr<const LinOp>>

A LinOpFactory represents a higher order mapping which transforms one linear operator into another.

In Ginkgo, every linear solver is viewed as a mapping. For example, given an s.p.d linear system \(Ax = b\), the solution \(x = A^{-1}b\) can be computed using the CG method. This algorithm can be represented in terms of linear operators and mappings between them as follows:

  • A Cg::Factory is a higher order mapping which, given an input operator \(A\), returns a new linear operator \(A^{-1}\)

    stored in “CG

    format”

  • Storing the operator \(A^{-1}\) in “CG format” means that the data structure used to store the operator is just a simple pointer to the original matrix \(A\). The application \(x = A^{-1}b\) of such an operator can then be implemented by solving the linear system \(Ax = b\) using the CG method. This is achieved in code by having a special class for each of those “formats” (e.g. the “Cg” class defines such a format for the CG solver).

Another example of a LinOpFactory is a preconditioner. A preconditioner for a linear operator \(A\) is a linear operator \(M^{-1}\), which approximates \(A^{-1}\). In addition, it is stored in a way such that both the data of \(M^{-1}\) is cheap to compute from \(A\), and the operation \(x = M^{-1}b\) can be computed quickly. These operators are useful to accelerate the convergence of Krylov solvers. Thus, a preconditioner also fits into the LinOpFactory framework:

  • The factory maps a linear operator \(A\) into a preconditioner \(M^{-1}\) which is stored in suitable format (e.g. as a product of two factors in case of ILU preconditioners).

  • The resulting linear operator implements the application operation \(x = M^{-1}b\) depending on the format the preconditioner is stored in (e.g. as two triangular solves in case of ILU)

#

Example: using CG in Ginkgo#

// Suppose A is a matrix, b a rhs vector, and x an initial guess
// Create a CG which runs for at most 1000 iterations, and stops after
// reducing the residual norm by 6 orders of magnitude
auto cg_factory = solver::Cg<>::build()
    .with_max_iters(1000)
    .with_rel_residual_goal(1e-6)
    .on(cuda);
// create a linear operator which represents the solver
auto cg = cg_factory->generate(A);
// solve the system
cg->apply(b, x);

Subclassed by

Public Functions

virtual std::unique_ptr<ReuseData> create_empty_reuse_data() const#

Creates an empty ReuseData. The default implementation returns a plain ReuseData.

std::unique_ptr<LinOp> generate_reuse(
std::shared_ptr<const LinOp> input,
ReuseData &reuse_data,
) const#

Generates a LinOp like generate(), reusing what earlier calls recorded in reuse_data.

A ReuseData is empty when it comes from create_empty_reuse_data() and initialized after the first successful generate_reuse(). With an empty reuse_data, this is a full generation that also fills it. Once initialized, the input must have the same size and sparsity pattern. What is reused is up to the factory and may depend on values, so the result can differ from generate(). Factories without reuse support just call generate().

As with generate(), the result is independent of reuse_data and of earlier results: it may outlive reuse_data, and later calls leave it unchanged. One ReuseData serves a whole sequence of generations:

auto reuse_data = factory->create_empty_reuse_data();
auto solver = factory->generate_reuse(matrix, *reuse_data);
// change the values of matrix, keeping its sparsity pattern
solver = factory->generate_reuse(matrix, *reuse_data);

If this throws, an empty reuse_data stays empty and an initialized one stays usable.

Note

This interface is experimental.

Parameters:
  • input – the input, copied to this factory’s executor if needed

  • reuse_data – created by this factory, or by one with the same type, parameters and executor

Throws:
  • NotSupported – if reuse_data has the wrong type for this factory

  • DimensionMismatch – if a dimension of input, like its size or its number of stored elements, differs from the input reuse_data was initialized with

  • InvalidStateError – if input does not fit reuse_data otherwise

Returns:

the generated LinOp, with this factory’s loggers

virtual void check_reuse_consistent(
const LinOp *input,
const ReuseData &reuse_data,
) const#

Throws if input can’t be used with reuse_data. generate_reuse() calls this first, also while reuse_data is still empty, so overrides have to accept any input then. Checks that need a converted input may run during the generation instead. The default implementation, used by factories without reuse support, accepts any input, since they generate from scratch.

Throws:
inline AbstractFactory(std::shared_ptr<const Executor> exec)#

Constructs a new factory on the specified executor.

Parameters:

exec – the executor where the factory should be constructed

class ReuseData#

Opaque base for the data a factory records in generate_reuse() to speed up later generations. Factories supporting reuse derive from it.

Note

This interface is experimental.

Subclassed by