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
solvernonlinear.hpp
Go to the documentation of this file.
1 /* -*- mode: c++; coding: utf-8; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; show-trailing-whitespace: t -*- vim:fenc=utf-8:ft=tcl:et:sw=4:ts=4:sts=4
2 
3  This file is part of the Feel library
4 
5  Author(s): Christophe Prud'homme <christophe.prudhomme@feelpp.org>
6  Date: 2007-07-02
7 
8  Copyright (C) 2007-2011 Universite Joseph Fourier (Grenoble I)
9 
10  This library is free software; you can redistribute it and/or
11  modify it under the terms of the GNU Lesser General Public
12  License as published by the Free Software Foundation; either
13  version 3.0 of the License, or (at your option) any later version.
14 
15  This library is distributed in the hope that it will be useful,
16  but WITHOUT ANY WARRANTY; without even the implied warranty of
17  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18  Lesser General Public License for more details.
19 
20  You should have received a copy of the GNU Lesser General Public
21  License along with this library; if not, write to the Free Software
22  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
23 */
29 #ifndef __SolverNonLinear_H
30 #define __SolverNonLinear_H 1
31 
32 #include <boost/function.hpp>
33 #include <boost/bind.hpp>
34 
35 #include <Eigen/Core>
36 #include <Eigen/LU>
37 #include <Eigen/Dense>
38 
39 #include <feel/feelalg/enums.hpp>
40 #include <feel/feelalg/glas.hpp>
41 #include <feel/feelcore/traits.hpp>
43 
44 
45 namespace Feel
46 {
47 
48 template<typename T> class Vector;
49 template<typename T> class MatrixSparse;
50 
61 template <typename T>
63 {
64 public:
65 
66 
70 
72  typedef boost::shared_ptr<SolverNonLinear<T> > self_ptrtype;
74  typedef boost::shared_ptr<self_type> solvernonlinear_ptrtype;
75 
76  typedef T value_type;
77  typedef typename type_traits<T>::real_type real_type;
78 
79  typedef boost::shared_ptr<Preconditioner<T> > preconditioner_ptrtype;
80 
81  typedef boost::shared_ptr<Vector<value_type> > vector_ptrtype;
82  typedef boost::shared_ptr<MatrixSparse<value_type> > sparse_matrix_ptrtype;
83 
84  typedef ublas::matrix<value_type> dense_matrix_type;
85  typedef ublas::vector<value_type> dense_vector_type;
86 
87  typedef Eigen::Map< Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> > map_dense_matrix_type;
88  typedef Eigen::Map< Eigen::Matrix<double, Eigen::Dynamic, 1> > map_dense_vector_type;
89 
90  typedef boost::function<void ( const vector_ptrtype& X,
91  vector_ptrtype& R )> residual_function_type;
92  typedef boost::function<void ( const vector_ptrtype& X,
93  sparse_matrix_ptrtype& J )> jacobian_function_type;
94  typedef boost::function<void ( const vector_ptrtype& X,
95  vector_ptrtype& R,
96  sparse_matrix_ptrtype& J )> matvec_function_type;
97 
98  typedef boost::function<void ( dense_vector_type const& X,
99  dense_vector_type & R )> dense_residual_function_type;
100  typedef boost::function<void ( dense_vector_type const& X,
101  dense_matrix_type& J )> dense_jacobian_function_type;
102  typedef boost::function<void ( dense_vector_type const& X,
103  dense_vector_type& R,
104  dense_matrix_type& J )> dense_matvec_function_type;
105 
106  //eigen
107  typedef boost::function<void ( map_dense_vector_type const& X,
108  map_dense_vector_type & R )> map_dense_residual_function_type;
109  typedef boost::function<void ( map_dense_vector_type const& X,
110  map_dense_matrix_type& J )> map_dense_jacobian_function_type;
111  typedef boost::function<void ( map_dense_vector_type const& X,
112  map_dense_vector_type& R,
113  map_dense_matrix_type& J )> map_dense_matvec_function_type;
114 
115 
117 
121 
125  SolverNonLinear(WorldComm const& worldComm = Environment::worldComm() );
126 
130  SolverNonLinear( SolverNonLinear const & );
131 
135  virtual ~SolverNonLinear();
136 
137 
142  static solvernonlinear_ptrtype build( po::variables_map const& vm, std::string const& prefix = "",
143  WorldComm const& worldComm = Environment::worldComm() );
144 
149  static solvernonlinear_ptrtype build( SolverPackage solver_package, WorldComm const& worldComm = Environment::worldComm() );
150 
154  virtual void init () = 0;
155 
156 
158 
162 
163 
165 
169 
173  WorldComm const& comm() const { return M_worldComm; }
174  WorldComm const& worldComm() const { return M_worldComm; }
175 
180  bool initialized () const
181  {
182  return M_is_initialized;
183  }
184 
188  virtual void clear () {}
189 
195  {
196  return M_prec_matrix_structure;
197  }
198 
199  SolverNonLinearType getType() const
200  {
201  return M_snl_type;
202  }
203 
204  double getAbsoluteResidualTol() const
205  {
206  return M_absoluteResidualTol;
207  }
208  double getRelativeResidualTol() const
209  {
210  return M_relativeResidualTol;
211  }
212  double getAbsoluteSolutionTol() const
213  {
214  return M_absoluteSolutionTol;
215  }
216 
217  uint getNbItMax() const
218  {
219  return M_nbItMax;
220  }
221 
222  int reuseJacobian() const
223  {
224  return M_reuse_jac;
225  }
226  int reusePreconditioner() const
227  {
228  return M_reuse_prec;
229  }
230 
234  std::string const& prefix() const{ return M_prefix; }
235 
236 
238 
242 
246  void setPrefix( std::string const& p ) { M_prefix = p; }
247 
252  virtual void setPrecMatrixStructure( MatrixStructure mstruct )
253  {
254  // warning : in boths cases!
255  if ( M_preconditioner )
256  M_preconditioner->setPrecMatrixStructure(mstruct);
257 
258  M_prec_matrix_structure = mstruct;
259  }
260 
261 
265  void setType( SolverNonLinearType snl_type )
266  {
267  M_snl_type=snl_type;
268  }
273  {
274  return M_snl_type;
275  }
276 
280  void setKspSolverType ( const SolverType st )
281  {
282  M_kspSolver_type = st;
283  }
284 
289  {
290  return M_kspSolver_type;
291  }
292 
297  {
298  if ( M_preconditioner )
299  return M_preconditioner->type();
300 
301  return M_preconditioner_type;
302  }
303 
308  {
309  if ( M_preconditioner )
310  M_preconditioner->setType( pct );
311 
312  else
313  M_preconditioner_type = pct;
314  }
315 
319  void attachPreconditioner( preconditioner_ptrtype preconditioner )
320  {
321  if ( this->M_is_initialized )
322  {
323  std::cerr<<"Preconditioner must be attached before the solver is initialized!"<<std::endl;
324  }
325 
326  M_preconditioner_type = SHELL_PRECOND;
327  M_preconditioner = preconditioner;
328  }
329 
333  void setMatSolverPackageType ( const MatSolverPackageType mspackt )
334  {
335  M_matSolverPackage_type = mspackt;
336  }
340  MatSolverPackageType matSolverPackageType () const
341  {
343  }
344 
345 
346 
356  virtual void setReuse( int jac=1, int prec=1 )
357  {
358  M_reuse_jac=jac;
359  M_reuse_prec=prec;
360  }
361 
365  void setRelativeResidualTol( double tol )
366  {
367  M_relativeResidualTol = tol;
368  }
369  void setAbsoluteResidualTol( double tol )
370  {
371  M_absoluteResidualTol = tol;
372  }
373  void setAbsoluteSolutionTol( double tol )
374  {
375  M_absoluteSolutionTol = tol;
376  }
377 
381  void setNbItMax( uint n )
382  {
383  M_nbItMax=n;
384  }
386 
387 
391  bool showSNESMonitor() const { return M_showSNESMonitor; }
392  void setShowSNESMonitor(bool b) { M_showSNESMonitor=b; }
393 
397  bool showKSPMonitor() const { return M_showKSPMonitor; }
398  void setShowKSPMonitor(bool b) { M_showKSPMonitor=b; }
399 
403  bool showSNESConvergedReason() const { return M_showSNESConvergedReason; }
404  void setShowSNESConvergedReason( bool b ) { M_showSNESConvergedReason=b; }
405 
409  bool showKSPConvergedReason() const { return M_showKSPConvergedReason; }
410  void setShowKSPConvergedReason( bool b ) { M_showKSPConvergedReason=b; }
411 
412 
416  double rtoleranceKSP() const { return M_rtoleranceKSP; }
417  void setRtoleranceKSP( double tol ) { M_rtoleranceKSP=tol; }
421  double dtoleranceKSP() const { return M_dtoleranceKSP; }
422  void setDtoleranceKSP( double tol ) { M_dtoleranceKSP=tol; }
426  double atoleranceKSP() const { return M_atoleranceKSP; }
427  void setAtoleranceKSP( double tol ) { M_atoleranceKSP=tol; }
431  size_type maxitKSP() const { return M_maxitKSP; }
432  void setMaxitKSP( size_type n) { M_maxitKSP=n; }
433 
434 
438 
442  virtual std::pair<int, real_type> solve ( sparse_matrix_ptrtype&, // System Jacobian Matrix
443  vector_ptrtype&, // Solution vector
444  vector_ptrtype&, // Residual vector
445  const double, // Stopping tolerance
446  const unsigned int ) = 0; // N. Iterations
447 
451  virtual std::pair<unsigned int, real_type> solve ( dense_matrix_type&, // System Jacobian Matrix
452  dense_vector_type&, // Solution vector
453  dense_vector_type&, // Residual vector
454  const double, // Stopping tolerance
455  const unsigned int ) = 0; // N. Iterations
456 
460  virtual std::pair<unsigned int, real_type> solve ( map_dense_matrix_type&, // System Jacobian Matrix
461  map_dense_vector_type&, // Solution vector
462  map_dense_vector_type&, // Residual vector
463  const double, // Stopping tolerance
464  const unsigned int ) = 0; // N. Iterations
465 
470  residual_function_type residual;
471 
476  jacobian_function_type jacobian;
477 
484  matvec_function_type matvec;
485 
490  dense_residual_function_type dense_residual;
491 
496  map_dense_residual_function_type map_dense_residual;
497 
502  dense_jacobian_function_type dense_jacobian;
503 
508  map_dense_jacobian_function_type map_dense_jacobian;
509 
516  dense_matvec_function_type dense_matvec;
517 
524  map_dense_matvec_function_type map_dense_matvec;
525 
527 
528 
529 
530 protected:
531 
532  WorldComm M_worldComm;
533 
538 
539  MatrixStructure M_prec_matrix_structure;
540 
541  std::string M_prefix;
542 
547 
552 
557 
561  preconditioner_ptrtype M_preconditioner;
562 
566  MatSolverPackageType M_matSolverPackage_type;
567 
572  double M_absoluteResidualTol;
573 
578 
582  uint M_nbItMax;
583 
588 
593 
594  bool M_showKSPMonitor, M_showSNESMonitor;
595  bool M_showKSPConvergedReason, M_showSNESConvergedReason;
596 
613 
614 
615 };
616 
617 
618 } // Feel
619 #endif /* __SolverNonLinear_H */

Generated on Fri Oct 25 2013 14:24:24 for Feel++ by doxygen 1.8.4