Logo  0.95.0-final
Finite Element Embedded Library and Language in C++
Feel++ Feel++ on Github Feel++ community
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
solverlinear.hpp
1 // $Id: linear_solver.h,v 1.2 2005/02/22 22:17:34 jwpeterson Exp $
2 
3 // The libMesh Finite Element Library.
4 // Copyright (C) 2002-2005 Benjamin S. Kirk, John W. Peterson
5 
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 3.0 of the License, or (at your option) any later version.
10 
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
15 
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
19 
20 #ifndef __linear_solver_h__
21 #define __linear_solver_h__
22 
23 #include <boost/mpi/communicator.hpp>
24 
25 #include <feel/feelcore/parameter.hpp>
26 #include <feel/feelalg/enums.hpp>
28 #include <feel/feelcore/traits.hpp>
29 
30 namespace Feel
31 {
32 namespace mpi=boost::mpi;
33 template<typename T> class Vector;
34 template<typename T> class MatrixSparse;
35 
44 template <typename T>
46 {
47 
48 public:
49 
50  typedef SolverLinear<T> self_type;
51  typedef boost::shared_ptr<SolverLinear<T> > self_ptrtype;
52 
53  typedef T value_type;
54  typedef typename type_traits<T>::real_type real_type;
55 
56  typedef boost::shared_ptr<Preconditioner<T> > preconditioner_ptrtype;
57 
61  SolverLinear ( WorldComm const& worldComm = Environment::worldComm() );
62 
66  SolverLinear ( po::variables_map const& vm, WorldComm const& worldComm = Environment::worldComm() );
67 
71  virtual ~SolverLinear ();
72 
73  WorldComm const& worldComm() const
74  {
75  return M_worldComm;
76  }
77  void setWorldComm( WorldComm const& worldComm )
78  {
79  M_worldComm=worldComm;
80  }
81 
86  bool initialized () const
87  {
88  return M_is_initialized;
89  }
90 
91 
95  virtual void clear () {}
96 
100  virtual void init () = 0;
101 
105  po::variables_map vm() const
106  {
107  return M_vm;
108  }
109 
113  value_type rTolerance() const
114  {
115  return M_rtolerance;
116  }
117 
121  value_type dTolerance() const
122  {
123  return M_dtolerance;
124  }
125 
129  value_type aTolerance() const
130  {
131  return M_atolerance;
132  }
133 
138  {
139  return M_solver_type;
140  }
141 
146  {
147  return M_maxit;
148  }
149 
153  std::string const& prefix() const{ return M_prefix; }
154 
158  void setPrefix( std::string const& p ) { M_prefix = p; }
159 
164  BOOST_PARAMETER_MEMBER_FUNCTION( ( void ),
165  setTolerances,
166  tag,
167  ( required
168  ( rtolerance,( double ) )
169  )
170  ( optional
171  ( maxit,( size_type ), 1000 )
172  ( atolerance,( double ), 1e-50 )
173  ( dtolerance,( double ), 1e5 )
174  ) )
175  {
176  M_rtolerance = rtolerance;
177  M_dtolerance = dtolerance;
178  M_atolerance = atolerance;
179  M_maxit=maxit;
180  }
181 
185  void setSolverType ( const SolverType st )
186  {
187  M_solver_type = st;
188  }
189 
194  {
195  if ( M_preconditioner )
196  return M_preconditioner->type();
197 
198  return M_preconditioner_type;
199  }
200 
205  {
206  if ( M_preconditioner )
207  M_preconditioner->setType( pct );
208 
209  else
210  M_preconditioner_type = pct;
211  }
212 
216  void attachPreconditioner( preconditioner_ptrtype preconditioner )
217  {
218  if ( this->M_is_initialized )
219  {
220  std::cerr<<"Preconditioner must be attached before the solver is initialized!"<<std::endl;
221  }
222 
223  M_preconditioner_type = SHELL_PRECOND;
224  M_preconditioner = preconditioner;
225  }
226 
227  void setFieldSplitType( const FieldSplitType fst )
228  {
229  M_fieldSplit_type = fst;
230  }
231 
232  FieldSplitType fieldSplitType() const
233  {
234  return M_fieldSplit_type;
235  }
236 
240  void setMatSolverPackageType ( const MatSolverPackageType mspackt )
241  {
242  M_matSolverPackage_type = mspackt;
243  }
247  MatSolverPackageType matSolverPackageType () const
248  {
250  }
251 
257  {
258  return M_prec_matrix_structure;
259  }
260 
265  virtual void setPrecMatrixStructure( MatrixStructure mstruct )
266  {
267  // warning : in boths cases!
268  if ( M_preconditioner )
269  M_preconditioner->setPrecMatrixStructure(mstruct);
270 
271  M_prec_matrix_structure = mstruct;
272  }
273 
277  bool showKSPMonitor() const { return M_showKSPMonitor; }
278  void setShowKSPMonitor( bool b ) { M_showKSPMonitor=b; }
279 
283  bool showKSPConvergedReason() const { return M_showKSPConvergedReason; }
284  void setShowKSPConvergedReason( bool b ) { M_showKSPConvergedReason=b; }
285 
300  virtual /*std::pair<unsigned int, real_type>*/
301  boost::tuple<bool,unsigned int, real_type>
302  solve ( MatrixSparse<T> const& mat,
303  Vector<T>& x,
304  Vector<T> const& b,
305  const double tolerance,
306  const unsigned int maxit,
307  bool transpose
308  ) = 0;
309 
310 
311 
326  virtual
327  /*std::pair<unsigned int, real_type>*/boost::tuple<bool,unsigned int, real_type>
328  solve ( MatrixSparse<T> const& mat,
329  MatrixSparse<T> const& prec,
330  Vector<T>& x,
331  Vector<T> const& b,
332  const double tolerance,
333  const unsigned int maxit,
334  bool transpose
335  ) = 0;
336 
337 
338 protected:
339 
343  void setInitialized( bool init )
344  {
346  }
347 
348 private:
349 
350  //mpi communicator
351  WorldComm M_worldComm;
352 
353 protected:
354 
356  po::variables_map M_vm;
357 
358  std::string M_prefix;
359 
361  double M_rtolerance;
362 
364  double M_dtolerance;
365 
367  double M_atolerance;
368 
371 
376 
381 
385  preconditioner_ptrtype M_preconditioner;
386 
387 
388  FieldSplitType M_fieldSplit_type;
389 
393  MatSolverPackageType M_matSolverPackage_type;
394 
395 
400 
401  MatrixStructure M_prec_matrix_structure;
402 
403  bool M_showKSPMonitor;
404  bool M_showKSPConvergedReason;
405 };
406 
407 
408 
409 
410 /*----------------------- inline functions ----------------------------------*/
411 template <typename T>
412 inline
413 SolverLinear<T>::SolverLinear ( WorldComm const& worldComm ) :
414  M_worldComm( worldComm ),
415  M_solver_type ( GMRES ),
416  M_preconditioner_type ( LU_PRECOND ),
417  M_preconditioner(),
418  M_is_initialized ( false ),
419  M_prec_matrix_structure( SAME_NONZERO_PATTERN ),
420  M_showKSPMonitor( false ),
421  M_showKSPConvergedReason( false )
422 {
423 }
424 
425 template <typename T>
426 inline
427 SolverLinear<T>::SolverLinear ( po::variables_map const& vm, WorldComm const& worldComm ) :
428  M_worldComm( worldComm ),
429  M_vm( vm ),
430  M_solver_type ( GMRES ),
431  M_preconditioner_type ( LU_PRECOND ),
432  M_is_initialized ( false ),
433  M_prec_matrix_structure( SAME_NONZERO_PATTERN ),
434  M_showKSPMonitor( false ),
435  M_showKSPConvergedReason( false )
436 {
437 }
438 
439 
440 
441 template <typename T>
442 inline
444 {
445  this->clear ();
446 }
447 
448 
449 } // Feel
450 
451 #endif // #ifdef __solver_h__
bool M_is_initialized
Definition: solverlinear.hpp:399
size_type maxIterations() const
Definition: solverlinear.hpp:145
Definition: solverlinear.hpp:33
value_type dTolerance() const
Definition: solverlinear.hpp:121
SolverType solverType() const
Definition: solverlinear.hpp:137
PreconditionerType
Definition: feelalg/enums.hpp:114
SolverLinear(WorldComm const &worldComm=Environment::worldComm())
Definition: solverlinear.hpp:413
double M_rtolerance
relative tolerance
Definition: solverlinear.hpp:361
virtual void init()=0
bool initialized() const
Definition: solverlinear.hpp:86
std::string const & prefix() const
Definition: solverlinear.hpp:153
SolverType M_solver_type
Definition: solverlinear.hpp:375
MatSolverPackageType matSolverPackageType() const
Definition: solverlinear.hpp:247
MatSolverPackageType M_matSolverPackage_type
Definition: solverlinear.hpp:393
bool showKSPMonitor() const
Definition: solverlinear.hpp:277
MatrixStructure
Definition: feelalg/enums.hpp:146
void setSolverType(const SolverType st)
Definition: solverlinear.hpp:185
void attachPreconditioner(preconditioner_ptrtype preconditioner)
Definition: solverlinear.hpp:216
po::variables_map vm() const
Definition: solverlinear.hpp:105
SolverType
Definition: feelalg/enums.hpp:89
FieldSplitType
Definition: feelalg/enums.hpp:138
virtual boost::tuple< bool, unsigned int, real_type > solve(MatrixSparse< T > const &mat, Vector< T > &x, Vector< T > const &b, const double tolerance, const unsigned int maxit, bool transpose)=0
void setPreconditionerType(const PreconditionerType pct)
Definition: solverlinear.hpp:204
size_type M_maxit
maximum number of iterations
Definition: solverlinear.hpp:370
virtual void setPrecMatrixStructure(MatrixStructure mstruct)
Definition: solverlinear.hpp:265
void setMatSolverPackageType(const MatSolverPackageType mspackt)
Definition: solverlinear.hpp:240
virtual ~SolverLinear()
Definition: solverlinear.hpp:443
size_t size_type
Indices (starting from 0)
Definition: feelcore/feel.hpp:319
preconditioner_ptrtype M_preconditioner
Definition: solverlinear.hpp:385
Definition: solverlinear.hpp:45
virtual MatrixStructure precMatrixStructure() const
Definition: solverlinear.hpp:256
Definition: matrixsparse.hpp:76
value_type rTolerance() const
Definition: solverlinear.hpp:113
PreconditionerType M_preconditioner_type
Definition: solverlinear.hpp:380
value_type aTolerance() const
Definition: solverlinear.hpp:129
bool showKSPConvergedReason() const
Definition: solverlinear.hpp:283
double M_atolerance
absolute tolerance
Definition: solverlinear.hpp:367
PreconditionerType preconditionerType() const
Definition: solverlinear.hpp:193
void setPrefix(std::string const &p)
Definition: solverlinear.hpp:158
virtual void clear()
Definition: solverlinear.hpp:95
void setInitialized(bool init)
Definition: solverlinear.hpp:343
double M_dtolerance
divergence tolerance
Definition: solverlinear.hpp:364

Generated on Sun Dec 22 2013 13:11:13 for Feel++ by doxygen 1.8.5