Rheolef  7.2
an efficient C++ finite element environment
solver_umfpack.h
Go to the documentation of this file.
1 #ifndef _RHEOLEF_SOLVER_UMFPACK_H
2 #define _RHEOLEF_SOLVER_UMFPACK_H
23 // solver implementation: interface
24 //
25 
26 #include "rheolef/config.h"
27 
28 #ifdef _RHEOLEF_HAVE_UMFPACK
29 
30 #include "rheolef/solver.h"
31 extern "C" {
32 #include <umfpack.h>
33 }
34 namespace rheolef {
35 
36 // =======================================================================
37 // rep
38 // =======================================================================
39 template<class T, class M>
41 public:
42 // typedef:
43 
45  typedef typename base::size_type size_type;
47 
48 // allocator:
49 
51  explicit solver_umfpack_rep (const csr<T,M>& a, const solver_option& opt = solver_option());
53  bool initialized() const { return true; }
54  void update_values (const csr<T,M>& a);
56 
57 // accessors:
58 
59  vec<T,M> trans_solve (const vec<T,M>& rhs) const;
60  vec<T,M> solve (const vec<T,M>& rhs) const;
61  determinant_type det() const { return _det; }
62 
63 protected:
64 // data:
65  int _n;
66  void* _numeric;
67  int* _ia_p;
68  int* _ja_p;
69  T* _a_p; // TODO: only used when T=double yet
71  double _control [UMFPACK_INFO];
72  mutable double _info [UMFPACK_INFO];
73 
74 // internal:
75  void _set_control();
76  void _destroy();
77  void _solve (int transpose_flag, const vec<T,M>& b, vec<T,M>& x) const;
78 private:
81 };
82 template <class T, class M>
83 inline
86 {
87  typedef solver_umfpack_rep<T,M> rep;
88  return new_macro (rep(*this));
89 }
90 
91 } // namespace rheolef
92 #endif // _RHEOLEF_HAVE_UMFPACK
93 #endif // _RHEOLEF_SOLVER_UMFPACK_H
csr< T, M >::size_type size_type
Definition: solver.h:193
see the solver_option page for the full documentation
determinant_type det() const
void update_values(const csr< T, M > &a)
solver_abstract_rep< T, M > base
base::determinant_type determinant_type
void _solve(int transpose_flag, const vec< T, M > &b, vec< T, M > &x) const
vec< T, M > trans_solve(const vec< T, M > &rhs) const
double _control[UMFPACK_INFO]
solver_abstract_rep< T, M > * clone() const
double _info[UMFPACK_INFO]
vec< T, M > solve(const vec< T, M > &rhs) const
Expr1::float_type T
Definition: field_expr.h:230
This file is part of Rheolef.