DOKK / manpages / debian 10 / librheolef-dev / solver.2rheolef.en
solver(2rheolef) rheolef-7.0 solver(2rheolef)

solver - direct or interative solver interface

The class implements a matrix factorization: LU factorization for an unsymmetric matrix and Choleski fatorisation for a symmetric one.

Let a be a square invertible matrix in csr format (see csr(2)).


csr<Float> a;
We get the factorization by:


solver sa (a);
Each call to the direct solver for a*x = b writes either:


vec<Float> x = sa.solve(b);
When the matrix is modified in a computation loop but conserves its sparsity pattern, an efficient re-factorization writes:


sa.update_values (new_a);
x = sa.solve(b);

The symmetry of the matrix is tested via the a.is_symmetric() property (see csr(2)) while the choice between direct or iterative solver is switched by default from the a.pattern_dimension() value. When the pattern is 3D, an iterative method is faster and less memory consuming. Otherwise, for 1D or 2D problems, the direct method is preferred.

These default choices can be supersetted by using explicit options:


solver_option opt;
opt.iterative = true;
solver sa (a, opt);
When using an iterative method, the sa.solve(b) call the conjugate gradient when the matrix is symmetric, or the generalized minimum residual algorithm when the matrix is unsymmetric.

When using a direct method, the determinant of the a matrix can be computed as:


solver_option opt;
opt.iterative = false;
solver sa (a, opt);
cout << sa.det().mantissa << "*2^" << sa.det().exponant << endl;
The sa.det() method returns an object of type solver::determinant_type that contains a mantissa and an exponent in base 2. This feature is useful e.g. when tracking a change of sign in the determinant of a matrix.

When using an iterative method, the default is to do no preconditionning. A suitable preconditionner can be supplied via:


solver_option opt;
opt.iterative = true;
solver sa (a, opt);
sa.set_preconditionner (pa);
x = sa.solve(b);
The supplied pa variable is also of type solver. A library of commonly used preconditionners is still in development.

template <class T, class M = rheo_default_memory_model>
class solver_basic : public smart_pointer<solver_rep<T,M> > {
public:
// typedefs:

typedef solver_rep<T,M> rep;
typedef smart_pointer<rep> base;
typedef typename rep::size_type size_type;
typedef typename rep::determinant_type determinant_type; // allocator:
solver_basic ();
explicit solver_basic (const csr<T,M>& a, const solver_option& opt = solver_option());
const solver_option& option() const;
void update_values (const csr<T,M>& a);
void set_preconditionner (const solver_basic<T,M>&); // accessors:
vec<T,M> trans_solve (const vec<T,M>& b) const;
vec<T,M> solve (const vec<T,M>& b) const;
determinant_type det() const; }; // factorizations: template <class T, class M> solver_basic<T,M> ldlt(const csr<T,M>& a, const solver_option& opt = solver_option()); template <class T, class M> solver_basic<T,M> lu (const csr<T,M>& a, const solver_option& opt = solver_option()); // preconditionners: template <class T, class M = rheo_default_memory_model> solver_basic<T,M> eye_basic(); inline solver_basic<Float> eye() { return eye_basic<Float>(); } template <class T, class M> solver_basic<T,M> ic0 (const csr<T,M>& a, const solver_option& opt = solver_option()); template <class T, class M> solver_basic<T,M> ilu0(const csr<T,M>& a, const solver_option& opt = solver_option()); template <class T, class M> solver_basic<T,M> ldlt_seq(const csr<T,M>& a, const solver_option& opt = solver_option()); template <class T, class M> solver_basic<T,M> lu_seq (const csr<T,M>& a, const solver_option& opt = solver_option()); typedef solver_basic<Float> solver;

csr(2), csr(2)

Copyright (C) 2000-2018 Pierre Saramito <Pierre.Saramito@imag.fr> GPLv3+: GNU GPL version 3 or later <http://gnu.org/licenses/gpl.html>. This is free software: you are free to change and redistribute it. There is NO WARRANTY, to the extent permitted by law.

rheolef-7.0 rheolef-7.0