Rheolef  7.1
an efficient C++ finite element environment
solver_eigen.h
Go to the documentation of this file.
1 #ifndef _RHEOLEF_SOLVER_EIGEN_H
2 #define _RHEOLEF_SOLVER_EIGEN_H
23 // solver implementation: interface
24 //
25 
26 #include "rheolef/config.h"
27 
28 #ifdef _RHEOLEF_HAVE_EIGEN
29 
30 #include "rheolef/solver.h"
31 #pragma GCC diagnostic push
32 #pragma GCC diagnostic ignored "-Weffc++"
33 #include <Eigen/Sparse>
34 #pragma GCC diagnostic pop
35 
36 namespace rheolef {
37 
38 // =======================================================================
39 // rep
40 // =======================================================================
41 template<class T, class M>
43 public:
44 // typedef:
45 
47  typedef typename base::size_type size_type;
49 
50 // allocator:
51 
52  explicit solver_eigen_rep (const csr<T,M>& a, const solver_option& opt = solver_option());
55  void update_values (const csr<T,M>& a);
56  bool initialized() const { return true; }
57 
58 // accessors:
59 
60  vec<T,M> trans_solve (const vec<T,M>& rhs) const;
61  vec<T,M> solve (const vec<T,M>& rhs) const;
62  determinant_type det() const { return _det; }
63 
64 protected:
65 // data:
67  Eigen::SparseLU<Eigen::SparseMatrix<T,Eigen::ColMajor>, Eigen::COLAMDOrdering<int> >
69  Eigen::SimplicialLDLT<Eigen::SparseMatrix<T,Eigen::ColMajor>, Eigen::Lower, Eigen::AMDOrdering<int> >
72 };
73 // -----------------------------------------------------------------------------
74 // inlined
75 // -----------------------------------------------------------------------------
76 template<class T, class M>
77 inline
79  : solver_abstract_rep<T,M>(opt),
80  _a(a),
81  _superlu_a(),
82  _ldlt_a(),
83  _det()
84 {
85  update_values (a);
86 }
87 template<class T, class M>
88 inline
90  : solver_abstract_rep<T,M>(x.option()),
91  _a(x._a),
92  _superlu_a(),
93  _ldlt_a(),
94  _det()
95 {
96  // Eigen::SuperLU copy cstor is non-copyable, so re-initialize for a copy
97  update_values (_a);
98 }
99 template <class T, class M>
100 inline
103 {
104  typedef solver_eigen_rep<T,M> rep;
105  return new_macro (rep(*this));
106 }
107 
108 } // namespace rheolef
109 #endif // _RHEOLEF_HAVE_EIGEN
110 #endif // _RHEOLEF_SOLVER_EIGEN_H
csr< T, M >::size_type size_type
Definition: solver.h:193
determinant_type det() const
Definition: solver_eigen.h:62
solver_eigen_rep(const csr< T, M > &a, const solver_option &opt=solver_option())
Definition: solver_eigen.h:78
void update_values(const csr< T, M > &a)
Definition: solver_eigen.cc:42
base::size_type size_type
Definition: solver_eigen.h:47
solver_abstract_rep< T, M > base
Definition: solver_eigen.h:46
base::determinant_type determinant_type
Definition: solver_eigen.h:48
Eigen::SparseLU< Eigen::SparseMatrix< T, Eigen::ColMajor >, Eigen::COLAMDOrdering< int > > _superlu_a
Definition: solver_eigen.h:68
vec< T, M > trans_solve(const vec< T, M > &rhs) const
Eigen::SimplicialLDLT< Eigen::SparseMatrix< T, Eigen::ColMajor >, Eigen::Lower, Eigen::AMDOrdering< int > > _ldlt_a
Definition: solver_eigen.h:70
solver_abstract_rep< T, M > * clone() const
Definition: solver_eigen.h:102
vec< T, M > solve(const vec< T, M > &rhs) const
Definition: solver_eigen.cc:85
determinant_type _det
Definition: solver_eigen.h:71
see the solver_option page for the full documentation
Expr1::float_type T
Definition: field_expr.h:261
This file is part of Rheolef.
Expr1::memory_type M
Definition: vec_expr_v2.h:416