Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
mic.h
Go to the documentation of this file.
1#ifndef _RHEOLEF_MIC_H
2#define _RHEOLEF_MIC_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 {
78} // namespace rheolef
79
80#include "rheolef/solver.h"
81#ifdef _RHEOLEF_HAVE_EIGEN
82#pragma GCC diagnostic push
83#pragma GCC diagnostic ignored "-Weffc++"
84#include <Eigen/Sparse>
85#pragma GCC diagnostic pop
86
87namespace rheolef {
88
89// =======================================================================
90// rep
91// =======================================================================
92template<class T, class M>
94public:
95// typedef:
96
98 typedef typename base::size_type size_type;
99 typedef typename base::determinant_type determinant_type;
100
101// allocator:
102
103 explicit solver_mic_rep (const csr<T,M>& a, const T& shift);
106 bool good() const { return true; }
107 void update_values (const csr<T,M>& a);
108
109// accessors:
110
111 vec<T,M> trans_solve (const vec<T,M>& rhs) const;
112 vec<T,M> solve (const vec<T,M>& rhs) const;
113
114protected:
115// data:
118 Eigen::IncompleteCholesky<T, Eigen::Lower, Eigen::AMDOrdering<int> > _mic_a;
119};
120// -----------------------------------------------------------------------------
121// inlined
122// -----------------------------------------------------------------------------
123template<class T, class M>
124inline
127 _a(a),
128 _shift(shift),
129 _mic_a()
130{
131 update_values (a);
132}
133template<class T, class M>
134inline
136 : solver_abstract_rep<T,M>(x.option()),
137 _a(x._a),
138 _shift(x._shift),
139 _mic_a()
140{
141 // Eigen::IncompleteCholesky copy cstor is non-copyable, so re-initialize for a copy
143}
144template <class T, class M>
145inline
146solver_abstract_rep<T,M>*
148{
149 typedef solver_mic_rep<T,M> rep;
150 return new_macro (rep(*this));
151}
152// =======================================================================
153// preconditioner interface
154// =======================================================================
155//<verbatim:
156template <class T, class M>
158 const csr<T,M>& a,
159 const T & shift = 1e-3);
160//>verbatim:
161
162template <class T, class M>
163inline
164solver_basic<T,M> mic (const csr<T,M>& a, const T& shift)
165{
166 using rep = solver_mic_rep<T,M>;
168 p.solver_basic<T,M>::base::operator= (new_macro(rep(a,shift)));
169 return p;
170}
171
172} // namespace rheolef
173#endif // _RHEOLEF_HAVE_EIGEN
174#endif // _RHEOLEF_MIC_H
see the csr page for the full documentation
Definition csr.h:317
csr< T, M > _a
Definition mic.h:116
void update_values(const csr< T, M > &a)
Definition mic.cc:34
solver_mic_rep(const solver_mic_rep &)
base::size_type size_type
Definition mic.h:98
solver_abstract_rep< T, M > base
Definition mic.h:97
solver_mic_rep(const csr< T, M > &a, const T &shift)
Definition mic.h:125
base::determinant_type determinant_type
Definition mic.h:99
vec< T, M > trans_solve(const vec< T, M > &rhs) const
Definition mic.cc:72
bool good() const
Definition mic.h:106
solver_abstract_rep< T, M > * clone() const
Definition mic.h:147
vec< T, M > solve(const vec< T, M > &rhs) const
Definition mic.cc:60
Eigen::IncompleteCholesky< T, Eigen::Lower, Eigen::AMDOrdering< int > > _mic_a
Definition mic.h:118
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 > mic(const csr< T, M > &a, const T &shift=1e-3)
Definition mic.h:164
Definition sphere.icc:25
Expr1::memory_type M