Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
ilut.h
Go to the documentation of this file.
1#ifndef _RHEOLEF_ILUT_H
2#define _RHEOLEF_ILUT_H
3//
4// This file is part of Rheolef.
5//
6// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
7//
8// Rheolef is free software; you can redistribute it and/or modify
9// it under the terms of the GNU General Public License as published by
10// the Free Software Foundation; either version 2 of the License, or
11// (at your option) any later version.
12//
13// Rheolef is distributed in the hope that it will be useful,
14// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16// GNU General Public License for more details.
17//
18// You should have received a copy of the GNU General Public License
19// along with Rheolef; if not, write to the Free Software
20// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21//
22// =========================================================================
23// AUTHOR: Pierre.Saramito@imag.fr
24// DATE: 28 february 2018
25
26namespace rheolef {
84} // namespace rheolef
85
86#include "rheolef/solver.h"
87#ifdef _RHEOLEF_HAVE_EIGEN
88#pragma GCC diagnostic push
89#pragma GCC diagnostic ignored "-Weffc++"
90#include <Eigen/Sparse>
91#pragma GCC diagnostic pop
92
93namespace rheolef {
94// =======================================================================
95// rep
96// =======================================================================
97template<class T, class M>
99public:
100// typedef:
101
103 typedef typename base::size_type size_type;
104 typedef typename base::determinant_type determinant_type;
105
106// allocator:
107
108 explicit solver_ilut_rep (const csr<T,M>& a, size_type fill_factor, T drop_tol);
111 bool good() const { return true; }
112 void update_values (const csr<T,M>& a);
113
114// accessors:
115
116 vec<T,M> trans_solve (const vec<T,M>& rhs) const;
117 vec<T,M> solve (const vec<T,M>& rhs) const;
118
119protected:
120// data:
124 Eigen::IncompleteLUT<T> _ilut_a;
125};
126// -----------------------------------------------------------------------------
127// inlined
128// -----------------------------------------------------------------------------
129template<class T, class M>
130inline
133 _a(a),
134 _fill_factor(fill_factor),
135 _drop_tol(drop_tol),
136 _ilut_a()
137{
138 update_values (a);
139}
140template<class T, class M>
141inline
143 : solver_abstract_rep<T,M>(x.option()),
144 _a(x._a),
145 _fill_factor(x._fill_factor),
146 _drop_tol(x._drop_tol),
147 _ilut_a()
148{
149 // Eigen::IncompleteLUT copy cstor is non-copyable, so re-initialize for a copy
151}
152template <class T, class M>
153inline
154solver_abstract_rep<T,M>*
156{
157 typedef solver_ilut_rep<T,M> rep;
158 return new_macro (rep(*this));
159}
160// =======================================================================
161// preconditioner interface
162// =======================================================================
163//<verbatim:
164template <class T, class M>
166 const csr<T,M>& a,
167 size_t fill_factor = 10,
168 T drop_tol = 1e3*std::numeric_limits<T>::epsilon());
169//>verbatim:
170
171template <class T, class M>
172inline
173solver_basic<T,M> ilut (const csr<T,M>& a, size_t fill_factor, T drop_tol)
174{
175 using rep = solver_ilut_rep<T,M>;
177 p.solver_basic<T,M>::base::operator= (new_macro(rep(a,fill_factor,drop_tol)));
178 return p;
179}
180
181} // namespace rheolef
182#endif // _RHEOLEF_HAVE_EIGEN
183#endif // _RHEOLEF_ILUT_H
see the csr page for the full documentation
Definition csr.h:317
csr< T, M >::size_type size_type
Definition solver.h:193
solver_ilut_rep(const csr< T, M > &a, size_type fill_factor, T drop_tol)
Definition ilut.h:131
Eigen::IncompleteLUT< T > _ilut_a
Definition ilut.h:124
size_type _fill_factor
Definition ilut.h:122
void update_values(const csr< T, M > &a)
Definition ilut.cc:34
base::size_type size_type
Definition ilut.h:103
solver_abstract_rep< T, M > base
Definition ilut.h:102
base::determinant_type determinant_type
Definition ilut.h:104
solver_ilut_rep(const solver_ilut_rep &)
vec< T, M > trans_solve(const vec< T, M > &rhs) const
Definition ilut.cc:73
bool good() const
Definition ilut.h:111
solver_abstract_rep< T, M > * clone() const
Definition ilut.h:155
vec< T, M > solve(const vec< T, M > &rhs) const
Definition ilut.cc:61
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.
solver_basic< T, M > ilut(const csr< T, M > &a, size_t fill_factor=10, T drop_tol=1e3 *std::numeric_limits< T >::epsilon())
Definition ilut.h:173
Definition sphere.icc:25
Expr1::memory_type M