dune-fem  2.8-git
cginverseoperator.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_CGINVERSEOPERATOR_HH
2 #define DUNE_FEM_CGINVERSEOPERATOR_HH
3 
4 #include <limits>
5 #include <memory>
6 #include <type_traits>
7 
13 
14 namespace Dune
15 {
16 
17  namespace Fem
18  {
19 
20  // ConjugateGradientSolver
21  // -----------------------
22 
29  template< class Operator>
31  {
33 
34  public:
37 
42  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
43 
48 
51 
52 
53  private:
54  static_assert( (std::is_same< DomainFunctionType, RangeFunctionType >::value),
55  "DomainFunctionType must equal RangeFunctionType." );
56 
57  public:
66  unsigned int maxIterations,
67  int errorMeasure,
68  bool verbose )
69  : epsilon_( epsilon ),
70  maxIterations_( maxIterations ),
71  errorMeasure_( errorMeasure ),
72  verbose_( verbose && Fem::Parameter::verbose() ),
73  averageCommTime_( 0.0 ),
74  realCount_( 0 )
75  {}
77  unsigned int maxIterations,
78  bool verbose,
79  const ParameterReader &parameter = Parameter::container() )
80  : epsilon_( epsilon ),
81  maxIterations_( maxIterations ),
82  errorMeasure_( 1 ),
83  verbose_( verbose && Fem::Parameter::verbose() ),
84  averageCommTime_( 0.0 ),
85  realCount_( 0 )
86  {}
87 
94  unsigned int maxIterations,
95  const ParameterReader &parameter = Parameter::container() )
96  : epsilon_( epsilon ),
97  maxIterations_( maxIterations ),
98  errorMeasure_( 1 ),
99  verbose_( parameter.getValue< bool >( "fem.solver.verbose", false ) ),
100  averageCommTime_( 0.0 ),
101  realCount_( 0 )
102  {}
103 
104  // ConjugateGradientSolver ( const ThisType & )=delete;
105 
117  void solve ( const OperatorType &op, const RangeFunctionType &b, DomainFunctionType &x ) const;
118 
131  void solve ( const OperatorType &op, const PreconditionerType &p,
132  const RangeFunctionType &b, DomainFunctionType &x ) const;
133 
135  unsigned int iterations () const
136  {
137  return realCount_;
138  }
139 
140  void setMaxIterations ( unsigned int maxIterations ) { maxIterations_ = maxIterations; }
141 
142 
144  double averageCommTime() const
145  {
146  return averageCommTime_;
147  }
148 
149  protected:
151  unsigned int maxIterations_;
153  const bool verbose_;
154  mutable double averageCommTime_;
155  mutable unsigned int realCount_;
156  };
157 
158 
159  namespace Solver
160  {
161  // CGInverseOperator
162  // -----------------
163 
169  template< class DiscreteFunction >
171  : public Fem::Operator< DiscreteFunction, DiscreteFunction >
172  {
175 
176  public:
179 
182 
184  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
185 
186  // non-const version
187  using BaseType::finalize;
188 
189  private:
191 
192  public:
200  CGInverseOperator ( RealType redEps, RealType absLimit,
201  unsigned int maxIter, bool verbose,
202  const ParameterReader &parameter = Parameter::container() )
203  : solver_( absLimit, maxIter, verbose, parameter ),
204  parameter_( parameter )
205  {}
206 
208  : solver_( param.tolerance(),
209  param.maxIterations(),
210  param.errorMeasure(),
211  param.verbose() ),
212  parameter_( param )
213  {
214  }
215 
216 
223  CGInverseOperator ( RealType redEps, RealType absLimit,
224  unsigned int maxIter,
225  const ParameterReader &parameter = Parameter::container() )
226  : solver_( absLimit, maxIter, parameter ),
227  parameter_( parameter )
228  {}
229 
230  CGInverseOperator ( RealType redEps, RealType absLimit,
231  const ParameterReader &parameter = Parameter::container() )
232  : CGInverseOperator( redEps, absLimit, std::numeric_limits< unsigned int >::max(), parameter )
233  {}
234 
244  RealType redEps, RealType absLimit,
245  unsigned int maxIter, bool verbose,
246  const ParameterReader &parameter = Parameter::container() )
247  : CGInverseOperator( redEps, absLimit, maxIter, verbose, parameter )
248  {
249  bind( op );
250  }
251 
252 
261  RealType redEps, RealType absLimit,
262  unsigned int maxIter,
263  const ParameterReader &parameter = Parameter::container() )
264  : CGInverseOperator( redEps, absLimit, maxIter, parameter )
265  {
266  bind( op );
267  }
268 
270  RealType redEps, RealType absLimit,
271  const ParameterReader &parameter = Parameter::container() )
272  : CGInverseOperator( redEps, absLimit, parameter )
273  {
274  bind( op );
275  }
276 
287  const PreconditionerType &precond,
288  RealType redEps, RealType absLimit,
289  unsigned int maxIter, bool verbose,
290  const ParameterReader &parameter = Parameter::container() )
291  : CGInverseOperator( redEps, absLimit, maxIter, verbose )
292  {
293  bind( op, precond );
294  }
295 
305  const PreconditionerType &precond,
306  RealType redEps, RealType absLimit,
307  const ParameterReader &parameter = Parameter::container() )
308  : CGInverseOperator( redEps, absLimit, parameter )
309  {
310  bind( op, precond );
311  }
312 
314  const PreconditionerType &precond,
315  RealType redEps, RealType absLimit,
316  unsigned int maxIter,
317  const ParameterReader &parameter = Parameter::container() )
318  : CGInverseOperator( redEps, absLimit, maxIter, parameter )
319  {
320  bind( op, precond );
321  }
322 
323  void bind ( const OperatorType &op ) { operator_ = &op; preconditioner_ = nullptr; }
324  void bind ( const OperatorType &op, const PreconditionerType &precond )
325  {
326  operator_ = &op;
327  preconditioner_ = &precond;
328  }
329  void unbind () { operator_ = nullptr; preconditioner_ = nullptr; }
330 
339  virtual void operator()( const DomainFunctionType &arg, RangeFunctionType &dest ) const
340  {
341  prepare();
342  apply(arg,dest);
343  const_cast< ThisType& > (*this).finalize();
344  }
345 
346  template<typename... A>
347  inline void prepare(A... ) const
348  {}
349 
358  virtual void apply( const DomainFunctionType &arg, RangeFunctionType &dest ) const
359  {
360  assert(operator_);
361  if(preconditioner_)
362  solver_.solve( *operator_, *preconditioner_, arg, dest );
363  else
364  solver_.solve( *operator_, arg, dest );
365  }
366 
368  unsigned int iterations () const
369  {
370  return solver_.iterations();
371  }
372 
373  void setMaxIterations ( unsigned int maxIter ) { solver_.setMaxIterations( maxIter ); }
374 
376  double averageCommTime() const
377  {
378  return solver_.averageCommTime();
379  }
380 
381  protected:
382  const OperatorType *operator_ = nullptr;
386  };
387  }
388 
389  // CGInverseOperator
390  // -----------------
391 
399  template< class DiscreteFunction,
402  : public Fem::Solver::CGInverseOperator< DiscreteFunction >
403  {
405 
406  public:
410 
412 
414  typedef Op OperatorType;
415 
417  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
418 
419  // Preconditioner is to approximate op^-1 !
421 
423  : BaseType( param )
424  {}
425 
426 
434  CGInverseOperator ( RealType redEps, RealType absLimit,
435  unsigned int maxIter, bool verbose,
436  const ParameterReader &parameter = Parameter::container() )
437  : BaseType( redEps, absLimit, maxIter, verbose, parameter )
438  {}
439 
446  CGInverseOperator ( RealType redEps, RealType absLimit,
447  unsigned int maxIter,
448  const ParameterReader &parameter = Parameter::container() )
449  : BaseType( redEps, absLimit, maxIter, parameter )
450  {}
451 
452  CGInverseOperator ( RealType redEps, RealType absLimit,
453  const ParameterReader &parameter = Parameter::container() )
454  : BaseType( redEps, absLimit, parameter )
455  {}
456 
465  template< class LinearOperator, std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0 >
467  RealType redEps, RealType absLimit,
468  unsigned int maxIter, bool verbose,
469  const ParameterReader &parameter = Parameter::container() )
470  : BaseType( redEps, absLimit, maxIter, verbose, parameter )
471  {
472  bind( op );
473  }
474 
482  template< class LinearOperator, std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0 >
484  RealType redEps, RealType absLimit,
485  unsigned int maxIter,
486  const ParameterReader &parameter = Parameter::container() )
487  : BaseType( redEps, absLimit, maxIter, parameter )
488  {
489  bind( op );
490  }
491 
492  template< class LinearOperator, std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0 >
494  RealType redEps, RealType absLimit,
495  const ParameterReader &parameter = Parameter::container() )
496  : BaseType( redEps, absLimit, parameter )
497  {
498  bind( op );
499  }
500 
510  const PreconditioningType &precond,
511  RealType redEps, RealType absLimit,
512  unsigned int maxIter, bool verbose,
513  const ParameterReader &parameter = Parameter::container() )
514  : BaseType( op, precond, redEps, absLimit, maxIter, verbose, parameter )
515  {}
516 
526  const PreconditioningType &precond,
527  RealType redEps, RealType absLimit,
528  unsigned int maxIter,
529  const ParameterReader &parameter = Parameter::container() )
530  : BaseType( op, precond, redEps, absLimit, maxIter, parameter )
531  {}
532 
534  const PreconditioningType &precond,
535  RealType redEps, RealType absLimit,
536  const ParameterReader &parameter = Parameter::container() )
537  : BaseType( op, precond, redEps, absLimit, parameter )
538  {}
539 
540 
541 
542  template< class LinearOperator, std::enable_if_t< std::is_base_of< OperatorType, LinearOperator >::value, int > = 0 >
543  void bind ( const LinearOperator &op )
544  {
545  BaseType::bind( op );
546  checkPreconditioning( op );
547  }
548  void unbind ()
549  {
551  precondObj_.reset();
552  }
553 
554  protected:
555  template< class LinearOperator >
556  void checkPreconditioning( const LinearOperator &linearOp )
557  {
558  bool preconditioning = false;
559  if (!parameter_.parameter().exists(parameter_.keyPrefix()+"preconditioning.method") &&
560  parameter_.parameter().exists("fem.preconditioning"))
561  {
562  preconditioning = parameter_.parameter().template getValue< bool >( "fem.preconditioning", false );
563  std::cout << "WARNING: using deprecated parameter `fem.preconditioning` use "
564  << parameter_.keyPrefix() << "preconditioning.method instead\n";
565  }
566  else
567  preconditioning = parameter_.preconditionMethod(
568  {
569  SolverParameter::none, // no preconditioning
570  SolverParameter::jacobi // Jacobi preconditioning
571  });
572  if( preconditioning && std::is_base_of< AssembledOperator< DomainFunctionType, DomainFunctionType > ,LinearOperator > :: value )
573  {
574  // create diagonal preconditioner
577  }
578  }
579 
581  using BaseType::parameter_;
582  std::unique_ptr< PreconditioningType > precondObj_;
583  };
584 
585  // Implementation of ConjugateGradientSolver
586  // -----------------------------------------
587 
588  template< class Operator >
591  {
592 
593  RealType tolerance = (epsilon_ * epsilon_);
594  if (errorMeasure_ == 1)
595  tolerance *= b.normSquaredDofs( );
596 
597  averageCommTime_ = 0.0;
598 
599  RangeFunctionType h( b );
600  op( x, h );
601 
602  RangeFunctionType r( h );
603  r -= b;
604 
605  RangeFunctionType p( b );
606  p -= h;
607 
608  RealType prevResiduum = 0;
609  RealType residuum = r.normSquaredDofs( );
610 
611  for( realCount_ = 0; (residuum > tolerance) && (realCount_ < maxIterations_); ++realCount_ )
612  {
613  if( realCount_ > 0 )
614  {
615  assert( residuum/prevResiduum == residuum/prevResiduum );
616  p *= (residuum / prevResiduum);
617  p -= r;
618  }
619 
620  op( p, h );
621 
622  RangeFieldType pdoth = p.scalarProductDofs( h );
623  const RangeFieldType alpha = residuum / pdoth;
624  assert( alpha == alpha );
625  x.axpy( alpha, p );
626  r.axpy( alpha, h );
627 
628  prevResiduum = residuum;
629  residuum = r.normSquaredDofs( );
630 
631  double exchangeTime = h.space().communicator().exchangeTime();
632  if( verbose_ )
633  {
634  std::cout << "CG-Iteration: " << realCount_ << ", Residuum: " << std::sqrt(residuum)
635  << ", Tolerance: " << std::sqrt(tolerance) << std::endl;
636  // only for parallel apps
637  if( b.space().gridPart().comm().size() > 1 )
638  std::cout << "Communication needed: " << exchangeTime << " s" << std::endl;
639  }
640 
641  averageCommTime_ += exchangeTime;
642  }
643  }
644 
645 
646  template< class Operator >
648  ::solve ( const OperatorType &op, const PreconditionerType &precond, const RangeFunctionType &b, DomainFunctionType &x ) const
649  {
650  RealType tolerance = (epsilon_ * epsilon_);
651  if (errorMeasure_ == 1)
652  tolerance *= b.normSquaredDofs( );
653 
654  averageCommTime_ = 0.0;
655 
656  RangeFunctionType h( b );
657  //h=Ax
658  op( x, h );
659 
660  //r=Ax-b
661  RangeFunctionType r( h );
662  r -= b;
663 
664  //p=b-A*x <= r_0 Deufelhard
665  RangeFunctionType p( b );
666  p -= h;
667 
668  //q=B*p <=q Deuf
669  RangeFunctionType q ( b );
670  precond(p,q);
671 
672  RangeFunctionType s (q);
673 
674  RangeFieldType prevResiduum = 0; // note that these will be real_type but require scalar product evaluation
675  RangeFieldType residuum = p.scalarProductDofs( q );//<p,Bp>
676 
677  for( realCount_ = 0; (std::real(residuum) > tolerance) && (realCount_ < maxIterations_); ++realCount_ )
678  {
679  if( realCount_ > 0 )
680  {
681  assert( residuum/prevResiduum == residuum/prevResiduum );
682  const RangeFieldType beta=residuum/prevResiduum;
683  q*=beta;
684  q+=(s);
685  }
686 
687  op( q, h );
688 
689  RangeFieldType qdoth = q.scalarProductDofs( h );
690  const RangeFieldType alpha = residuum / qdoth;//<p,Bp>/<q,Aq>
691  assert( alpha == alpha );
692  x.axpy( alpha, q );
693 
694  p.axpy( -alpha, h );//r_k
695 
696  precond(p,s); //B*r_k
697 
698  prevResiduum = residuum;//<rk-1,B*rk-1>
699 
700  residuum = p.scalarProductDofs( s );//<rk,B*rk>
701 
702  double exchangeTime = h.space().communicator().exchangeTime();
703  if( verbose_ )
704  {
705  std::cout << "CG-Iteration: " << realCount_ << ", Residuum: " << std::sqrt(residuum)
706  << ", Tolerance: " << std::sqrt(tolerance) << std::endl;
707  // only for parallel apps
708  if( b.space().gridPart().comm().size() > 1 )
709  std::cout << "Communication needed: " << exchangeTime << " s" << std::endl;
710  }
711 
712  averageCommTime_ += exchangeTime;
713  }
714  }
715 
716  } // namespace Fem
717 
718 } // namespace Dune
719 
720 #endif // #ifndef DUNE_FEM_CGINVERSEOPERATOR_HH
Definition: bindguard.hh:11
double real(const std::complex< Double > &x)
Definition: double.hh:906
static double sqrt(const Double &v)
Definition: double.hh:886
static double max(const Double &v, const double p)
Definition: double.hh:398
bool exists(const std::string &key) const
check, whether a parameter is defined
Definition: reader.hh:44
Container for User Specified Parameters.
Definition: io/parameter.hh:191
static ParameterContainer & container()
Definition: io/parameter.hh:193
abstract operator
Definition: operator.hh:34
virtual void finalize()
finalization of operator
Definition: operator.hh:61
DomainFunction DomainFunctionType
type of discrete function in the operator's domain
Definition: operator.hh:36
RangeFunction::RangeFieldType RangeFieldType
field type of the operator's range
Definition: operator.hh:43
DomainFunction::RangeFieldType DomainFieldType
field type of the operator's domain
Definition: operator.hh:41
RangeFunction RangeFunctionType
type of discrete function in the operator's range
Definition: operator.hh:38
abstract affine-linear operator
Definition: operator.hh:87
abstract matrix operator
Definition: operator.hh:124
linear solver using the CG algorithm
Definition: cginverseoperator.hh:31
ConjugateGradientSolver(const RealType &epsilon, unsigned int maxIterations, bool verbose, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:76
OperatorType::RangeFunctionType RangeFunctionType
type of the operator's range vectors
Definition: cginverseoperator.hh:47
ConjugateGradientSolver(RealType epsilon, unsigned int maxIterations, const ParameterReader &parameter=Parameter::container())
constructor
Definition: cginverseoperator.hh:93
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: cginverseoperator.hh:42
unsigned int iterations() const
number of iterations needed for last solve
Definition: cginverseoperator.hh:135
void solve(const OperatorType &op, const RangeFunctionType &b, DomainFunctionType &x) const
solve
Definition: cginverseoperator.hh:590
double averageCommTime() const
return average communication time during last solve
Definition: cginverseoperator.hh:144
void setMaxIterations(unsigned int maxIterations)
Definition: cginverseoperator.hh:140
void solve(const OperatorType &op, const PreconditionerType &p, const RangeFunctionType &b, DomainFunctionType &x) const
solve
Definition: cginverseoperator.hh:648
Fem::Operator< RangeFunctionType, DomainFunctionType > PreconditionerType
type of the preconditioner, maps from the range of the operator (the dual space) in it's domain
Definition: cginverseoperator.hh:50
const bool verbose_
Definition: cginverseoperator.hh:153
OperatorType::RangeFieldType RangeFieldType
field type of the operator's range vectors
Definition: cginverseoperator.hh:41
OperatorType::DomainFieldType DomainFieldType
field type of the operator's domain vectors
Definition: cginverseoperator.hh:39
double averageCommTime_
Definition: cginverseoperator.hh:154
OperatorType::DomainFunctionType DomainFunctionType
type of the operator's domain vectors
Definition: cginverseoperator.hh:45
const RealType epsilon_
Definition: cginverseoperator.hh:150
Operator OperatorType
type of the operators to invert
Definition: cginverseoperator.hh:36
ConjugateGradientSolver(const RealType &epsilon, unsigned int maxIterations, int errorMeasure, bool verbose)
constructor
Definition: cginverseoperator.hh:65
unsigned int realCount_
Definition: cginverseoperator.hh:155
int errorMeasure_
Definition: cginverseoperator.hh:152
unsigned int maxIterations_
Definition: cginverseoperator.hh:151
Inverse operator base on CG method. This is the base class for the cg solver and does not imvolve any...
Definition: cginverseoperator.hh:172
CGInverseOperator(const OperatorType &op, RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:260
const PreconditionerType * preconditioner_
Definition: cginverseoperator.hh:383
void prepare(A...) const
Definition: cginverseoperator.hh:347
void unbind()
Definition: cginverseoperator.hh:329
CGInverseOperator(RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:200
virtual void operator()(const DomainFunctionType &arg, RangeFunctionType &dest) const
application operator
Definition: cginverseoperator.hh:339
const OperatorType * operator_
Definition: cginverseoperator.hh:382
CGInverseOperator(const OperatorType &op, const PreconditionerType &precond, RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:304
CGInverseOperator(const OperatorType &op, RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:269
CGInverseOperator(RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:223
Fem::Operator< DomainFunctionType, RangeFunctionType > OperatorType
Definition: cginverseoperator.hh:180
unsigned int iterations() const
number of iterations needed for last solve
Definition: cginverseoperator.hh:368
SolverParameter parameter_
Definition: cginverseoperator.hh:385
virtual void apply(const DomainFunctionType &arg, RangeFunctionType &dest) const
application operator
Definition: cginverseoperator.hh:358
void bind(const OperatorType &op, const PreconditionerType &precond)
Definition: cginverseoperator.hh:324
CGInverseOperator(const OperatorType &op, const PreconditionerType &precond, RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:313
CGInverseOperator(RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:230
OperatorType::RangeFieldType RangeFieldType
Definition: cginverseoperator.hh:183
BaseType::DomainFunctionType DomainFunctionType
Definition: cginverseoperator.hh:177
CGInverseOperator(const OperatorType &op, const PreconditionerType &precond, RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:286
double averageCommTime() const
return average communication time during last solve
Definition: cginverseoperator.hh:376
Fem::Operator< RangeFunctionType, DomainFunctionType > PreconditionerType
Definition: cginverseoperator.hh:181
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: cginverseoperator.hh:184
BaseType::RangeFunctionType RangeFunctionType
Definition: cginverseoperator.hh:178
CGInverseOperator(const SolverParameter &param=SolverParameter(Parameter::container()))
Definition: cginverseoperator.hh:207
SolverType solver_
Definition: cginverseoperator.hh:384
void setMaxIterations(unsigned int maxIter)
Definition: cginverseoperator.hh:373
void bind(const OperatorType &op)
Definition: cginverseoperator.hh:323
CGInverseOperator(const OperatorType &op, RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:243
Inverse operator base on CG method. Uses a runtime parameter fem.preconditioning which enables diagon...
Definition: cginverseoperator.hh:403
CGInverseOperator(RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:434
CGInverseOperator(const OperatorType &op, const PreconditioningType &precond, RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:533
Op OperatorType
type of operator
Definition: cginverseoperator.hh:414
CGInverseOperator(RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:452
CGInverseOperator(const LinearOperator &op, RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:483
Dune::FieldTraits< RangeFieldType >::real_type RealType
Definition: cginverseoperator.hh:417
SolverParameter SolverParameterType
Definition: cginverseoperator.hh:407
void bind(const LinearOperator &op)
Definition: cginverseoperator.hh:543
CGInverseOperator(RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:446
CGInverseOperator(const LinearOperator &op, RealType redEps, RealType absLimit, const ParameterReader &parameter=Parameter::container())
Definition: cginverseoperator.hh:493
OperatorType::RangeFieldType RangeFieldType
Definition: cginverseoperator.hh:416
CGInverseOperator(const OperatorType &op, const PreconditioningType &precond, RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:509
Fem::Operator< RangeFunctionType, DomainFunctionType > PreconditioningType
Definition: cginverseoperator.hh:420
CGInverseOperator(const LinearOperator &op, RealType redEps, RealType absLimit, unsigned int maxIter, bool verbose, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:466
CGInverseOperator(const SolverParameter &param=SolverParameter(Parameter::container()))
Definition: cginverseoperator.hh:422
BaseType::RangeFunctionType RangeFunctionType
Definition: cginverseoperator.hh:409
CGInverseOperator(const OperatorType &op, const PreconditioningType &precond, RealType redEps, RealType absLimit, unsigned int maxIter, const ParameterReader &parameter=Parameter::container())
constructor of CGInverseOperator
Definition: cginverseoperator.hh:525
DomainFunctionType DestinationType
Definition: cginverseoperator.hh:411
void checkPreconditioning(const LinearOperator &linearOp)
Definition: cginverseoperator.hh:556
BaseType::DomainFunctionType DomainFunctionType
Definition: cginverseoperator.hh:408
std::unique_ptr< PreconditioningType > precondObj_
Definition: cginverseoperator.hh:582
void unbind()
Definition: cginverseoperator.hh:548
Precondtioner, multiplies with inverse of the diagonal works with.
Definition: diagonalpreconditioner.hh:152
Definition: solver/parameter.hh:15
static const int none
Definition: solver/parameter.hh:40
const std::string & keyPrefix() const
Definition: solver/parameter.hh:66
virtual int preconditionMethod(const std::vector< int > standardMethods, const std::vector< std::string > &additionalMethods={}, int defaultMethod=0) const
Definition: solver/parameter.hh:185
const ParameterReader & parameter() const
Definition: solver/parameter.hh:68
static const int jacobi
Definition: solver/parameter.hh:45