dune-fem  2.8-git
istlpreconditioner.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_ISTLPRECONDITIONERWRAPPER_HH
2 #define DUNE_FEM_ISTLPRECONDITIONERWRAPPER_HH
3 
4 #include <memory>
5 // standard diagonal preconditioner
7 
8 
9 #if HAVE_DUNE_ISTL
10 #include <dune/common/version.hh>
11 
12 #include <dune/istl/operators.hh>
13 #include <dune/istl/preconditioners.hh>
14 
15 #include <dune/istl/paamg/amg.hh>
16 #include <dune/istl/paamg/pinfo.hh>
17 #endif
18 
19 #include <dune/fem/solver/parameter.hh> // SolverParameter
21 
22 namespace Dune
23 {
24 
25  namespace Fem
26  {
27 
28  template <class MatrixImp>
30 
31  template <class MatrixImp>
33 
34  template <class MatrixImp>
36 
37  struct ISTLSolverParameter : public LocalParameter< SolverParameter, ISTLSolverParameter >
38  {
40  public:
41  using BaseType :: parameter ;
42  using BaseType :: keyPrefix ;
43 
45  : BaseType( parameter )
46  {}
47 
50  {}
51 
54  {}
55 
56  virtual double overflowFraction () const
57  {
58  return parameter().getValue< double >( keyPrefix() + "matrix.overflowfraction", 1.0 );
59  }
60 
61  virtual bool fastILUStorage () const
62  {
63  return parameter().getValue< bool >( keyPrefix() + "preconditioning.fastilustorage", true );
64  }
65  };
66 
67 #if HAVE_DUNE_ISTL
68  template< class MatrixObject, class X, class Y >
69  class FemDiagonalPreconditioner : public Preconditioner<X,Y>
70  {
71  public:
72  typedef typename MatrixObject :: ColumnDiscreteFunctionType DiscreteFunctionType ;
73  // select preconditioner for assembled operators
74  typedef DiagonalPreconditionerBase< DiscreteFunctionType, MatrixObject, true > PreconditionerType ;
75 
76  protected:
77  PreconditionerType diagonalPrecon_;
78 
79  public:
80  typedef typename X::field_type field_type;
81  FemDiagonalPreconditioner( const MatrixObject& mObj )
82  : diagonalPrecon_( mObj )
83  {}
84 
86  void pre (X& x, Y& b) override {}
87 
89  void apply (X& v, const Y& d) override
90  {
91  diagonalPrecon_.applyToISTLBlockVector( d, v );
92  }
93 
95  void post (X& x) override {}
96 
97 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
99  SolverCategory::Category category () const override { return SolverCategory::sequential; }
100 #endif // #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
101  };
102 
103 
104  template< class X, class Y >
105  class IdentityPreconditionerWrapper : public Preconditioner<X,Y>
106  {
107  template <class XImp, class YImp>
108  struct Apply
109  {
110  inline static void copy(XImp& v, const YImp& d)
111  {
112  }
113  };
114 
115  template <class XImp>
116  struct Apply<XImp,XImp>
117  {
118  inline static void copy(X& v, const Y& d)
119  {
120  v = d;
121  }
122  };
123 
124  public:
126  typedef X domain_type;
128  typedef Y range_type;
130  typedef typename X::field_type field_type;
131 
132 #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
133  enum {
135  category=SolverCategory::sequential };
136 #endif // #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
137 
139  IdentityPreconditionerWrapper(){}
140 
142  void pre (X& x, Y& b) override {}
143 
145  void apply (X& v, const Y& d) override
146  {
147  Apply< X, Y> :: copy( v, d );
148  }
149 
151  void post (X& x) override {}
152 
153 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
154  SolverCategory::Category category () const override { return SolverCategory::sequential; }
155 #endif // #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
156  };
157 
158 
162  template<class MatrixImp>
163  class PreconditionerWrapper
164  : public Preconditioner<typename MatrixImp :: RowBlockVectorType,
165  typename MatrixImp :: ColBlockVectorType>
166  {
167  typedef MatrixImp MatrixType;
168  typedef typename MatrixType :: RowBlockVectorType X;
169  typedef typename MatrixType :: ColBlockVectorType Y;
170 
171  // use BCRSMatrix type because of specializations in dune-istl
172  typedef typename MatrixType :: BaseType ISTLMatrixType ;
173 
174  typedef typename MatrixType :: CollectiveCommunictionType
175  CollectiveCommunictionType;
176 
177  // matrix adapter for AMG
178  //#if HAVE_MPI
179  // typedef Dune::OverlappingSchwarzOperator<
180  // ISTLMatrixType, X, Y, CollectiveCommunictionType> OperatorType ;
181  //#else
182  typedef MatrixAdapter< ISTLMatrixType, X, Y > OperatorType;
183  //#endif
184  mutable std::shared_ptr< OperatorType > op_;
185 
186  // auto pointer to preconditioning object
187  typedef Preconditioner<X,Y> PreconditionerInterfaceType;
188  mutable std::shared_ptr<PreconditionerInterfaceType> preconder_;
189 
190  // flag whether we have preconditioning, and if yes if it is AMG
191  const int preEx_;
192  const bool verbose_;
193 
194  template <class XImp, class YImp>
195  struct Apply
196  {
197  inline static void apply(PreconditionerInterfaceType& preconder,
198  XImp& v, const YImp& d)
199  {
200  }
201 
202  inline static void copy(XImp& v, const YImp& d)
203  {
204  }
205  };
206 
207  template <class XImp>
208  struct Apply<XImp,XImp>
209  {
210  inline static void apply(PreconditionerInterfaceType& preconder,
211  XImp& v, const XImp& d)
212  {
213  preconder.apply( v ,d );
214  }
215 
216  inline static void copy(X& v, const Y& d)
217  {
218  v = d;
219  }
220  };
221  public:
223  typedef X domain_type;
225  typedef Y range_type;
227  typedef typename X::field_type field_type;
228 
229 #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
230  enum {
232  category=SolverCategory::sequential };
233 #endif // #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
234 
236  PreconditionerWrapper (const PreconditionerWrapper& org, bool verbose)
237  : op_( org.op_ )
238  , preconder_(org.preconder_)
239  , preEx_(org.preEx_)
240  , verbose_( verbose )
241  {
242  }
243 
245  PreconditionerWrapper(bool verbose)
246  : op_()
247  , preconder_()
248  , preEx_( 0 )
249  , verbose_( verbose )
250  {}
251 
253  template <class PreconditionerType>
254  PreconditionerWrapper(MatrixType & matrix,
255  int iter,
256  field_type relax,
257  bool verbose,
258  const PreconditionerType* p)
259  : op_()
260  , preconder_( new PreconditionerType( matrix, iter, relax ) )
261  , preEx_( 1 )
262  , verbose_( verbose )
263  {
264  }
265 
266 
269  template <class PreconditionerType>
270  PreconditionerWrapper(MatrixType & matrix, bool verbose,
271  PreconditionerType* p)
272  : op_()
273  , preconder_( p )
274  , preEx_( 1 )
275  , verbose_( verbose )
276  {
277  }
278 
280  template <class PreconditionerType>
281  PreconditionerWrapper(MatrixType & matrix,
282  int iter,
283  field_type relax,
284  bool verbose,
285  const PreconditionerType* p ,
286  const CollectiveCommunictionType& comm)
287  //#if HAVE_MPI
288  // : op_( new OperatorType( matrix, comm ) )
289  //#else
290  : op_( new OperatorType( matrix ) )
291  //#endif
292  , preconder_( createAMGPreconditioner(comm, iter, relax, p) )
293  , preEx_( 2 )
294  , verbose_( verbose )
295  {
296  }
297 
299  void pre (X& x, Y& b) override
300  {
301  // all the sequentiel implemented Preconditioners do nothing in pre and post
302  // apply preconditioner
303  if( preEx_ > 1 )
304  {
305  //X tmp (x);
306  preconder_->pre(x,b);
307  //assert( std::abs( x.two_norm() - tmp.two_norm() ) < 1e-15);
308  }
309  }
310 
312  void apply (X& v, const Y& d) override
313  {
314  if( preEx_ )
315  {
316  // apply preconditioner
317  Apply<X,Y>::apply( *preconder_ , v, d );
318  }
319  else
320  {
321  Apply<X,Y>::copy( v, d );
322  }
323  }
324 
326  void post (X& x) override
327  {
328  // all the implemented Preconditioners do nothing in pre and post
329  // apply preconditioner
330  if( preEx_ > 1 )
331  {
332  preconder_->post(x);
333  }
334  }
335 
337 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
338  SolverCategory::Category category () const override
339  {
340  return (preconder_ ? preconder_->category() : SolverCategory::sequential);
341  }
342 #endif // #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
343 
344  protected:
345  template <class Smoother>
346  PreconditionerInterfaceType*
347  createAMGPreconditioner(const CollectiveCommunictionType& comm,
348  int iter, field_type relax, const Smoother* )
349  {
350  typedef typename Dune::FieldTraits< field_type>::real_type real_type;
351  typedef typename std::conditional< std::is_convertible<field_type,real_type>::value,
352  Dune::Amg::FirstDiagonal, Dune::Amg::RowSum >::type Norm;
353  typedef Dune::Amg::CoarsenCriterion<
354  Dune::Amg::UnSymmetricCriterion<ISTLMatrixType, Norm > > Criterion;
355 
356  typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
357 
358  SmootherArgs smootherArgs;
359 
360  smootherArgs.iterations = iter;
361  smootherArgs.relaxationFactor = relax ;
362 
363  int coarsenTarget=1200;
364  Criterion criterion(15,coarsenTarget);
365  criterion.setDefaultValuesIsotropic(2);
366  criterion.setAlpha(.67);
367  criterion.setBeta(1.0e-8);
368  criterion.setMaxLevel(10);
369  if( verbose_ && Parameter :: verbose() )
370  criterion.setDebugLevel( 1 );
371  else
372  criterion.setDebugLevel( 0 );
373 
374  /*
375  if( comm.size() > 1 )
376  {
377  typedef Dune::OwnerOverlapCopyCommunication<int> ParallelInformation;
378  ParallelInformation pinfo(MPI_COMM_WORLD);
379  typedef Dune::Amg::AMG<OperatorType, X, Smoother, ParallelInformation> AMG;
380  return new AMG(*op_, criterion, smootherArgs, pinfo);
381  }
382  else
383  */
384  {
385  // X == Y is needed for AMG
386  typedef Dune::Amg::AMG<OperatorType, X, Smoother> AMG;
387  return new AMG(*op_, criterion, smootherArgs);
388  //return new AMG(*op_, criterion, smootherArgs, 1, 1, 1, false);
389  }
390  }
391  };
392 
393 
394  // ISTLPreconditionerFactory
395  // -------------------------
396 
397  template < class MatrixObject >
398  class ISTLMatrixAdapterFactory ;
399 
400  template < class DomainSpace, class RangeSpace, class DomainBlock, class RangeBlock,
401  template <class, class, class, class> class MatrixObject >
402  class ISTLMatrixAdapterFactory< MatrixObject< DomainSpace, RangeSpace, DomainBlock, RangeBlock > >
403  {
404  typedef DomainSpace DomainSpaceType ;
405  typedef RangeSpace RangeSpaceType;
406 
407  public:
408  typedef MatrixObject< DomainSpaceType, RangeSpaceType, DomainBlock, RangeBlock > MatrixObjectType;
409  typedef typename MatrixObjectType :: MatrixType MatrixType;
410  typedef ISTLParallelMatrixAdapterInterface< MatrixType > MatrixAdapterInterfaceType;
411 
413  static std::unique_ptr< MatrixAdapterInterfaceType >
414  matrixAdapter( const MatrixObjectType& matrixObj,
415  const ISTLSolverParameter& param)
416  {
417  std::unique_ptr< MatrixAdapterInterfaceType > ptr;
418  if( matrixObj.domainSpace().continuous() )
419  {
420  typedef LagrangeParallelMatrixAdapter< MatrixType > MatrixAdapterImplementation;
421  ptr.reset( matrixAdapterObject( matrixObj, (MatrixAdapterImplementation *) nullptr, param ) );
422  }
423  else
424  {
425  typedef DGParallelMatrixAdapter< MatrixType > MatrixAdapterImplementation;
426  ptr.reset( matrixAdapterObject( matrixObj, (MatrixAdapterImplementation *) nullptr, param ) );
427  }
428  return ptr;
429  }
430 
431  protected:
432  template <class MatrixAdapterType>
433  static MatrixAdapterType*
434  matrixAdapterObject( const MatrixObjectType& matrixObj,
435  const MatrixAdapterType*,
436  const ISTLSolverParameter& param )
437  {
438  typedef typename MatrixAdapterType :: PreconditionAdapterType PreConType;
439  return new MatrixAdapterType( matrixObj.exportMatrix(),
440  matrixObj.domainSpace(), matrixObj.rangeSpace(), PreConType(param.verbose()) );
441  }
442 
443  };
444 
445 #ifndef DISABLE_ISTL_PRECONDITIONING
447  template < class Space, class DomainBlock, class RangeBlock,
448  template <class, class, class, class> class MatrixObject >
449  class ISTLMatrixAdapterFactory< MatrixObject< Space, Space, DomainBlock, RangeBlock > >
450  {
451  public:
452  typedef Space DomainSpaceType ;
453  typedef Space RangeSpaceType;
454 
455  typedef MatrixObject< DomainSpaceType, RangeSpaceType, DomainBlock, RangeBlock > MatrixObjectType;
456  typedef typename MatrixObjectType :: MatrixType MatrixType;
457  typedef typename MatrixType :: BaseType ISTLMatrixType ;
458  typedef Fem::PreconditionerWrapper< MatrixType > PreconditionAdapterType;
459 
460  protected:
461  template <class MatrixAdapterType, class PreconditionerType>
462  static MatrixAdapterType*
463  createMatrixAdapter(const MatrixAdapterType*,
464  const PreconditionerType* preconditioning,
465  MatrixType& matrix,
466  const DomainSpaceType& domainSpace,
467  const RangeSpaceType& rangeSpace,
468  const double relaxFactor,
469  std::size_t numIterations,
470  bool verbose)
471  {
472  typedef typename MatrixAdapterType :: PreconditionAdapterType PreConType;
473  PreConType preconAdapter(matrix, numIterations, relaxFactor, verbose, preconditioning );
474  return new MatrixAdapterType(matrix, domainSpace, rangeSpace, preconAdapter );
475  }
476 
477  template <class MatrixAdapterType, class PreconditionerType>
478  static MatrixAdapterType*
479  createAMGMatrixAdapter(const MatrixAdapterType*,
480  const PreconditionerType* preconditioning,
481  MatrixType& matrix,
482  const DomainSpaceType& domainSpace,
483  const RangeSpaceType& rangeSpace,
484  const double relaxFactor,
485  std::size_t numIterations,
486  bool solverVerbose)
487  {
488  typedef typename MatrixAdapterType :: PreconditionAdapterType PreConType;
489  PreConType preconAdapter(matrix, numIterations, relaxFactor, solverVerbose,
490  preconditioning, domainSpace.gridPart().comm() );
491  return new MatrixAdapterType(matrix, domainSpace, rangeSpace, preconAdapter );
492  }
493 
494  template <class MatrixAdapterType>
495  static MatrixAdapterType*
496  matrixAdapterObject( const MatrixObjectType& matrixObj,
497  const MatrixAdapterType*,
498  const ISTLSolverParameter& param )
499  {
500  int preconditioning = param.preconditionMethod(
501  { SolverParameter::none, // no preconditioner
502  SolverParameter::ssor, // SSOR preconditioner
503  SolverParameter::sor , // SOR preconditioner
504  SolverParameter::ilu , // ILU preconditioner (deprecated)
505  SolverParameter::gauss_seidel,// Gauss-Seidel preconditioner
506  SolverParameter::jacobi, // Jacobi preconditioner
507  SolverParameter::amg_ilu, // AMG with ILU-0 smoother (deprecated)
508  SolverParameter::amg_jacobi, // AMG with Jacobi smoother
509  SolverParameter::ildl // ILDL from istl
510  } );
511  const double relaxFactor = param.relaxation();
512  const size_t numIterations = param.preconditionerIteration();
513 
514  const DomainSpaceType& domainSpace = matrixObj.domainSpace();
515  const RangeSpaceType& rangeSpace = matrixObj.rangeSpace();
516 
517  MatrixType& matrix = matrixObj.exportMatrix();
518  const auto procs = domainSpace.gridPart().comm().size();
519 
520  typedef typename MatrixType :: BaseType ISTLMatrixType;
521  typedef typename MatrixAdapterType :: PreconditionAdapterType PreConType;
522 
523  typedef typename MatrixObjectType :: RowBlockVectorType RowBlockVectorType;
524  typedef typename MatrixObjectType :: ColumnBlockVectorType ColumnBlockVectorType;
525 
526  // no preconditioner
527  if( preconditioning == SolverParameter::none )
528  {
529  return new MatrixAdapterType( matrix, domainSpace, rangeSpace, PreConType(param.verbose()) );
530  }
531  // SSOR
532  else if( preconditioning == SolverParameter::ssor )
533  {
534  if( procs > 1 )
535  DUNE_THROW(InvalidStateException,"ISTL::SeqSSOR not working in parallel computations");
536 
537  typedef SeqSSOR<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType> PreconditionerType;
538  return createMatrixAdapter( (MatrixAdapterType *)nullptr,
539  (PreconditionerType*)nullptr,
540  matrix, domainSpace, rangeSpace, relaxFactor, numIterations, param.verbose() );
541  }
542  // SOR
543  else if(preconditioning == SolverParameter::sor )
544  {
545  if( procs > 1 )
546  DUNE_THROW(InvalidStateException,"ISTL::SeqSOR not working in parallel computations");
547 
548  typedef SeqSOR<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType> PreconditionerType;
549  return createMatrixAdapter( (MatrixAdapterType *)nullptr,
550  (PreconditionerType*)nullptr,
551  matrix, domainSpace, rangeSpace, relaxFactor, numIterations, param.verbose() );
552  }
553  // ILU
554  else if(preconditioning == SolverParameter::ilu)
555  {
556  if( procs > 1 )
557  DUNE_THROW(InvalidStateException,"ISTL::SeqILU not working in parallel computations");
558 
559  typedef SeqILU<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType> PreconditionerType;
560  typedef typename MatrixAdapterType :: PreconditionAdapterType PreConType;
561  // might need to set verbosity here?
562  PreConType preconAdapter( matrix, param.verbose(), new PreconditionerType( matrix, numIterations, relaxFactor, param.fastILUStorage()) );
563  return new MatrixAdapterType( matrix, domainSpace, rangeSpace, preconAdapter );
564  }
565  // Gauss-Seidel
566  else if(preconditioning == SolverParameter::gauss_seidel)
567  {
568  if( procs > 1 )
569  DUNE_THROW(InvalidStateException,"ISTL::SeqGS not working in parallel computations");
570 
571  typedef SeqGS<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType> PreconditionerType;
572  return createMatrixAdapter( (MatrixAdapterType *)nullptr,
573  (PreconditionerType*)nullptr,
574  matrix, domainSpace, rangeSpace, relaxFactor, numIterations, param.verbose() );
575  }
576  // Jacobi
577  else if(preconditioning == SolverParameter::jacobi)
578  {
579  if( numIterations == 1 ) // diagonal preconditioning
580  {
581  typedef FemDiagonalPreconditioner< MatrixObjectType, RowBlockVectorType, ColumnBlockVectorType > PreconditionerType;
582  typedef typename MatrixAdapterType :: PreconditionAdapterType PreConType;
583  PreConType preconAdapter( matrix, param.verbose(), new PreconditionerType( matrixObj ) );
584  return new MatrixAdapterType( matrix, domainSpace, rangeSpace, preconAdapter );
585  }
586  else if ( procs == 1 )
587  {
588  typedef SeqJac<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType> PreconditionerType;
589  return createMatrixAdapter( (MatrixAdapterType *)nullptr,
590  (PreconditionerType*)nullptr,
591  matrix, domainSpace, rangeSpace, relaxFactor, numIterations,
592  param.verbose());
593  }
594  else
595  {
596  DUNE_THROW(InvalidStateException,"ISTL::SeqJac(Jacobi) only working with istl.preconditioning.iterations: 1 in parallel computations");
597  }
598  }
599  // AMG ILU-0
600  else if(preconditioning == SolverParameter::amg_ilu)
601  {
602  // use original SeqILU because of some AMG traits classes.
603  typedef SeqILU<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType> PreconditionerType;
604  return createAMGMatrixAdapter( (MatrixAdapterType *)nullptr,
605  (PreconditionerType*)nullptr,
606  matrix, domainSpace, rangeSpace, relaxFactor, numIterations,
607  param.verbose());
608  }
609  // AMG Jacobi
610  else if(preconditioning == SolverParameter::amg_jacobi)
611  {
612  typedef SeqJac<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType,1> PreconditionerType;
613  return createAMGMatrixAdapter( (MatrixAdapterType *)nullptr,
614  (PreconditionerType*)nullptr,
615  matrix, domainSpace, rangeSpace, relaxFactor, numIterations,
616  param.verbose());
617  }
618  // ILDL
619  else if(preconditioning == SolverParameter::ildl)
620  {
621  if( procs > 1 )
622  DUNE_THROW(InvalidStateException,"ISTL::SeqILDL not working in parallel computations");
623 
624  PreConType preconAdapter( matrix, param.verbose(), new SeqILDL<ISTLMatrixType,RowBlockVectorType,ColumnBlockVectorType>( matrix , relaxFactor ) );
625  return new MatrixAdapterType( matrix, domainSpace, rangeSpace, preconAdapter );
626  }
627  else
628  {
629  preConErrorMsg(preconditioning);
630  }
631  return new MatrixAdapterType(matrix, domainSpace, rangeSpace, PreConType(param.verbose()) );
632  }
633 
634  typedef ISTLParallelMatrixAdapterInterface< MatrixType > MatrixAdapterInterfaceType;
635 
636  static void preConErrorMsg(int preCon)
637  {
638  // TODO: re-write
639  std::cerr << "ERROR: Wrong precoditioning number (p = " << preCon;
640  std::cerr <<") in ISTLMatrixObject! \n";
641  std::cerr <<"Valid values are: \n";
642  std::cerr <<"0 == no \n";
643  std::cerr <<"1 == SSOR \n";
644  std::cerr <<"2 == SOR \n";
645  std::cerr <<"3 == ILU-0 \n";
646  std::cerr <<"4 == ILU-n \n";
647  std::cerr <<"5 == Gauss-Seidel \n";
648  std::cerr <<"6 == Jacobi \n";
649 
650  assert( false );
651  DUNE_THROW(NotImplemented,"Wrong precoditioning selected");
652  }
653 
654  public:
656  static std::unique_ptr< MatrixAdapterInterfaceType >
657  matrixAdapter( const MatrixObjectType& matrixObj,
658  const ISTLSolverParameter& param)
659  {
660  const ISTLSolverParameter* parameter = dynamic_cast< const ISTLSolverParameter* > (&param);
661  std::unique_ptr< ISTLSolverParameter > paramPtr;
662  if( ! parameter )
663  {
664  paramPtr.reset( new ISTLSolverParameter( param ) );
665  parameter = paramPtr.operator->();
666  }
667 
668  std::unique_ptr< MatrixAdapterInterfaceType > ptr;
669  if( matrixObj.domainSpace().continuous() )
670  {
671  typedef LagrangeParallelMatrixAdapter< MatrixType > MatrixAdapterImplementation;
672  ptr.reset( matrixAdapterObject( matrixObj, (MatrixAdapterImplementation *) nullptr, *parameter ) );
673  }
674  else
675  {
676  typedef DGParallelMatrixAdapter< MatrixType > MatrixAdapterImplementation;
677  ptr.reset( matrixAdapterObject( matrixObj, (MatrixAdapterImplementation *) nullptr, *parameter ) );
678  }
679  return ptr;
680  }
681 
682  }; // end specialization of ISTLMatrixAdapterFactor for domainSpace == rangeSpace
683 #endif
684 
685 #endif // end HAVE_DUNE_ISTL
686 
687  } // namespace Fem
688 
689 } // namespace Dune
690 
691 #endif // #ifndef DUNE_FEM_PRECONDITIONERWRAPPER_HH
Definition: bindguard.hh:11
T getValue(const std::string &key) const
get mandatory parameter
Definition: reader.hh:159
static bool verbose()
obtain the cached value for fem.verbose
Definition: io/parameter.hh:445
static ParameterContainer & container()
Definition: io/parameter.hh:193
Definition: io/parameter.hh:551
Definition: istlpreconditioner.hh:29
Definition: istlpreconditioner.hh:35
Definition: istlpreconditioner.hh:38
virtual bool fastILUStorage() const
Definition: istlpreconditioner.hh:61
LocalParameter< SolverParameter, ISTLSolverParameter > BaseType
Definition: istlpreconditioner.hh:39
ISTLSolverParameter(const std::string &keyPrefix, const ParameterReader &parameter=Parameter::container())
Definition: istlpreconditioner.hh:48
virtual double overflowFraction() const
Definition: istlpreconditioner.hh:56
ISTLSolverParameter(const ParameterReader &parameter=Parameter::container())
Definition: istlpreconditioner.hh:44
ISTLSolverParameter(const SolverParameter &sp)
Definition: istlpreconditioner.hh:52
Definition: solver/parameter.hh:15
static const int none
Definition: solver/parameter.hh:40
const std::string & keyPrefix() const
Definition: solver/parameter.hh:66
static const int ssor
Definition: solver/parameter.hh:41
const ParameterReader & parameter() const
Definition: solver/parameter.hh:68
static const int sor
Definition: solver/parameter.hh:42
static const int ildl
Definition: solver/parameter.hh:48
static const int gauss_seidel
Definition: solver/parameter.hh:44
static const int ilu
Definition: solver/parameter.hh:43
static const int jacobi
Definition: solver/parameter.hh:45
static const int amg_ilu
Definition: solver/parameter.hh:46
static const int amg_jacobi
Definition: solver/parameter.hh:47