Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
solver_cholmod.h
Go to the documentation of this file.
1#ifndef _RHEOLEF_SOLVER_CHOLMOD_H
2#define _RHEOLEF_SOLVER_CHOLMOD_H
23// solver implementation: interface
24//
25#include "rheolef/solver.h"
26
27#if defined(_RHEOLEF_HAVE_CHOLMOD)
28#include <suitesparse/cholmod.h>
29#pragma GCC diagnostic push
30#pragma GCC diagnostic ignored "-Weffc++"
31#include <Eigen/Sparse>
32#include <Eigen/src/CholmodSupport/CholmodSupport.h>
33#pragma GCC diagnostic pop
34namespace rheolef {
35
36// =======================================================================
37// rep
38// =======================================================================
39template<class T, class M>
41public:
42// typedef:
43
45 typedef typename base::size_type size_type;
46 typedef typename base::determinant_type determinant_type;
47
48// allocator:
49
50 explicit solver_cholmod_rep (const csr<T,M>& a, const solver_option& opt = solver_option());
53 void update_values (const csr<T,M>& a);
54 bool initialized() const { return true; }
55
56// accessors:
57
58 vec<T,M> trans_solve (const vec<T,M>& rhs) const;
59 vec<T,M> solve (const vec<T,M>& rhs) const;
60 determinant_type det() const { return _det; }
61
62protected:
63// data:
65 Eigen::CholmodDecomposition<Eigen::SparseMatrix<double> >
68};
69// -----------------------------------------------------------------------------
70// inlined
71// -----------------------------------------------------------------------------
72template<class T, class M>
73inline
75 : solver_abstract_rep<T,M>(opt),
76 _a(a),
77 _llt_a(),
78 _det()
79{
80 update_values (a);
81}
82template<class T, class M>
83inline
85 : solver_abstract_rep<T,M>(x.option()),
86 _a(x._a),
87 _llt_a(),
88 _det()
89{
90 // Eigen::CholmodSupernodalLLT copy cstor is non-copyable, so re-initialize for a copy
92}
93template <class T, class M>
94inline
97{
98 typedef solver_cholmod_rep<T,M> rep;
99 return new_macro (rep(*this));
100}
101
102} // namespace rheolef
103#endif // _RHEOLEF_HAVE_CHOLMOD
104#endif // _RHEOLEF_SOLVER_CHOLMOD_H
see the csr page for the full documentation
Definition csr.h:317
determinant_type det() const
solver_cholmod_rep(const csr< T, M > &a, const solver_option &opt=solver_option())
Eigen::CholmodDecomposition< Eigen::SparseMatrix< double > > _llt_a
void update_values(const csr< T, M > &a)
solver_abstract_rep< T, M > base
base::determinant_type determinant_type
vec< T, M > trans_solve(const vec< T, M > &rhs) const
solver_abstract_rep< T, M > * clone() const
vec< T, M > solve(const vec< T, M > &rhs) const
see the solver_option page for the full documentation
see the vec page for the full documentation
Definition vec.h:79
Expr1::float_type T
Definition field_expr.h:230
This file is part of Rheolef.
Expr1::memory_type M