dune-pdelab  2.0.0
ovlp_amg_dg_backend.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_OVLP_AMG_DG_BACKEND_HH
2 #define DUNE_PDELAB_OVLP_AMG_DG_BACKEND_HH
3 
4 #include <dune/istl/matrixmatrix.hh>
5 
6 #include <dune/grid/common/datahandleif.hh>
7 
16 
17 namespace Dune {
18  namespace PDELab {
19  //***********************************************************
20  // A data handle / function to communicate matrix entries
21  // in the overlap. It turned out that it is actually not
22  // required to do that, but anyway it is there now.
23  //***********************************************************
24 
27  template<class GFS>
29  : public Dune::CommDataHandleIF<LocalGlobalMapDataHandle<GFS>,int>
30  {
31  public:
32  // some types from the grid view
33  typedef typename GFS::Traits::GridView GV;
34  typedef typename GV::IndexSet IndexSet;
35  typedef typename IndexSet::IndexType IndexType;
36  typedef typename GV::Grid Grid;
37  typedef typename Grid::Traits::GlobalIdSet GlobalIdSet;
38  typedef typename GlobalIdSet::IdType IdType;
39 
41  typedef int DataType;
42 
43  // map types
44  typedef std::map<IndexType,IdType> LocalToGlobalMap;
45  typedef std::map<IdType,IndexType> GlobalToLocalMap;
46 
48  bool contains (int dim, int codim) const
49  {
50  return (codim==0);
51  }
52 
54  bool fixedsize (int dim, int codim) const
55  {
56  return true;
57  }
58 
63  template<class EntityType>
64  size_t size (EntityType& e) const
65  {
66  return 1;
67  }
68 
70  template<class MessageBuffer, class EntityType>
71  void gather (MessageBuffer& buff, const EntityType& e) const
72  {
73  // fill map
74  IndexType myindex = indexSet.index(e);
75  IdType myid = globalIdSet.id(e);
76  l2g[myindex] = myid;
77  g2l[myid] = myindex;
78  //std::cout << gv.comm().rank() << ": gather local=" << myindex << " global=" << myid << std::endl;
79 
80  buff.write(0); // this is a dummy, we are not interested in the data
81  }
82 
87  template<class MessageBuffer, class EntityType>
88  void scatter (MessageBuffer& buff, const EntityType& e, size_t n)
89  {
90  DataType x;
91  buff.read(x); // this is a dummy, we are not interested in the data
92 
93  IndexType myindex = indexSet.index(e);
94  IdType myid = globalIdSet.id(e);
95  l2g[myindex] = myid;
96  g2l[myid] = myindex;
97  //std::cout << gv.comm().rank() << ": scatter local=" << myindex << " global=" << myid << std::endl;
98  }
99 
102  : gfs(gfs_), gv(gfs.gridView()), indexSet(gv.indexSet()),
103  grid(gv.grid()), globalIdSet(grid.globalIdSet()),
104  l2g(l2g_), g2l(g2l_)
105  {
106  }
107 
108  private:
109  const GFS& gfs;
110  GV gv;
111  const IndexSet& indexSet;
112  const Grid& grid;
113  const GlobalIdSet& globalIdSet;
114  LocalToGlobalMap& l2g;
115  GlobalToLocalMap& g2l;
116  };
117 
118 
119  // A DataHandle class to exchange rows of a sparse matrix
120  template<class GFS, class M> // mapper type and vector type
122  : public Dune::CommDataHandleIF<MatrixExchangeDataHandle<GFS,M>,
123  std::pair<typename GFS::Traits::GridView::Grid::Traits::GlobalIdSet::IdType,
124  typename M::block_type> >
125  {
126  public:
127  // some types from the grid view
128  typedef typename GFS::Traits::GridView GV;
129  typedef typename GV::IndexSet IndexSet;
130  typedef typename IndexSet::IndexType IndexType;
131  typedef typename GV::Grid Grid;
132  typedef typename Grid::Traits::GlobalIdSet GlobalIdSet;
133  typedef typename GlobalIdSet::IdType IdType;
134 
135  // some types from the matrix
136  typedef typename M::block_type B;
137 
139  typedef std::pair<IdType,B> DataType;
140 
141  // map types
142  typedef std::map<IndexType,IdType> LocalToGlobalMap;
143  typedef std::map<IdType,IndexType> GlobalToLocalMap;
144 
146  bool contains (int dim, int codim) const
147  {
148  return (codim==0);
149  }
150 
152  bool fixedsize (int dim, int codim) const
153  {
154  return false;
155  }
156 
161  template<class EntityType>
162  size_t size (EntityType& e) const
163  {
164  IndexType myindex = indexSet.index(e);
165  typename M::row_type myrow = m[myindex];
166  typename M::row_type::iterator endit=myrow.end();
167  size_t count=0;
168  for (typename M::row_type::iterator it=myrow.begin(); it!=endit; ++it)
169  {
170  typename LocalToGlobalMap::const_iterator find=l2g.find(it.index());
171  if (find!=l2g.end())
172  {
173  count++;
174  }
175  }
176  //std::cout << gv.comm().rank() << ": row=" << myindex << " size=" << count << std::endl;
177  return count;
178  }
179 
181  template<class MessageBuffer, class EntityType>
182  void gather (MessageBuffer& buff, const EntityType& e) const
183  {
184  IndexType myindex = indexSet.index(e);
185  //std::cout << gv.comm().rank() << ": begin gather local=" << myindex << std::endl;
186  typename M::row_type myrow = m[myindex];
187  typename M::row_type::iterator endit=myrow.end();
188  size_t count = 0;
189  for (typename M::row_type::iterator it=myrow.begin(); it!=endit; ++it)
190  {
191  typename LocalToGlobalMap::const_iterator find=l2g.find(it.index());
192  if (find!=l2g.end())
193  {
194  buff.write(std::make_pair(find->second,*it));
195  //std::cout << gv.comm().rank() << ": ==> local=" << find->first << " global=" << find->second << std::endl;
196  count++;
197  }
198  }
199  //std::cout << gv.comm().rank() << ": end gather row=" << myindex << " size=" << count << std::endl;
200  }
201 
206  template<class MessageBuffer, class EntityType>
207  void scatter (MessageBuffer& buff, const EntityType& e, size_t n)
208  {
209  IndexType myindex = indexSet.index(e);
210  std::cout << gv.comm().rank() << ": begin scatter local=" << myindex << " size=" << n << std::endl;
211  DataType x;
212  size_t count = 0;
213  for (size_t i=0; i<n; ++i)
214  {
215  buff.read(x);
216  std::cout << gv.comm().rank() << ": --> received global=" << x.first << std::endl;
217  typename GlobalToLocalMap::const_iterator find=g2l.find(x.first);
218  if (find!=g2l.end())
219  {
220  IndexType nbindex=find->second;
221  if (m.exists(myindex,nbindex))
222  {
223  m[myindex][nbindex] = x.second;
224  B block(x.second);
225  block -= m2[myindex][nbindex];
226  std::cout << gv.comm().rank() << ": compare i=" << myindex << " j=" << nbindex
227  << " norm=" << block.infinity_norm() << std::endl;
228 
229  count++;
230  //std::cout << gv.comm().rank() << ": --> local=" << find->first << " global=" << find->second << std::endl;
231  }
232  }
233  }
234  //std::cout << gv.comm().rank() << ": end scatter row=" << myindex << " size=" << count << std::endl;
235  }
236 
238  MatrixExchangeDataHandle (const GFS& gfs_, M& m_, const LocalToGlobalMap& l2g_, const GlobalToLocalMap& g2l_,M& m2_)
239  : gfs(gfs_), m(m_), gv(gfs.gridView()), indexSet(gv.indexSet()),
240  grid(gv.grid()), globalIdSet(grid.globalIdSet()),
241  l2g(l2g_), g2l(g2l_), m2(m2_)
242  {
243  }
244 
245  private:
246  const GFS& gfs;
247  M& m;
248  GV gv;
249  const IndexSet& indexSet;
250  const Grid& grid;
251  const GlobalIdSet& globalIdSet;
252  const LocalToGlobalMap& l2g;
253  const GlobalToLocalMap& g2l;
254  M& m2;
255  };
256 
257 
260  template<class GFS, class T, class A, int n, int m>
261  void restore_overlap_entries (const GFS& gfs, Dune::BCRSMatrix<Dune::FieldMatrix<T,n,m>,A>& matrix,
262  Dune::BCRSMatrix<Dune::FieldMatrix<T,n,m>,A>& matrix2)
263  {
264  typedef Dune::FieldMatrix<T,n,m> B;
265  typedef Dune::BCRSMatrix<B,A> M; // m is of type M
266  typedef typename LocalGlobalMapDataHandle<GFS>::LocalToGlobalMap LocalToGlobalMap;
267  typedef typename LocalGlobalMapDataHandle<GFS>::GlobalToLocalMap GlobalToLocalMap;
268 
269  // build up two maps to associate local indices and global ids
270  LocalToGlobalMap l2g;
271  GlobalToLocalMap g2l;
272  LocalGlobalMapDataHandle<GFS> lgdh(gfs,l2g,g2l);
273  if (gfs.gridView().comm().size()>1)
274  gfs.gridView().communicate(lgdh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
275 
276  // exchange matrix entries
277  MatrixExchangeDataHandle<GFS,M> mexdh(gfs,matrix,l2g,g2l,matrix2);
278  if (gfs.gridView().comm().size()>1)
279  gfs.gridView().communicate(mexdh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
280  }
281 
282  //***********************************************************
283  // The DG/AMG preconditioner in the overlapping case
284  //***********************************************************
285 
295  template<class DGGFS, class DGMatrix, class DGPrec, class DGCC,
296  class CGGFS, class CGPrec, class CGCC,
297  class P, class DGHelper, class Comm>
299  : public Dune::Preconditioner<typename Dune::PDELab::BackendVectorSelector<DGGFS,typename DGPrec::domain_type::field_type>::Type,
300  typename Dune::PDELab::BackendVectorSelector<DGGFS,typename DGPrec::range_type::field_type>::Type>
301  {
302  const DGGFS& dggfs;
303  DGMatrix& dgmatrix;
304  DGPrec& dgprec;
305  const DGCC& dgcc;
306  const CGGFS& cggfs;
307  CGPrec& cgprec;
308  const CGCC& cgcc;
309  P& p;
310  const DGHelper& dghelper;
311  const Comm& comm;
312  int n1,n2;
313 
314  public:
317  typedef typename V::BaseT X;
318  typedef typename W::BaseT Y;
321 
322  // define the category
323  enum {
325  category=Dune::SolverCategory::overlapping
326  };
327 
335  OvlpDGAMGPrec (const DGGFS& dggfs_, DGMatrix& dgmatrix_, DGPrec& dgprec_, const DGCC& dgcc_,
336  const CGGFS& cggfs_, CGPrec& cgprec_, const CGCC& cgcc_, P& p_,
337  const DGHelper& dghelper_, const Comm& comm_, int n1_, int n2_)
338  : dggfs(dggfs_), dgmatrix(dgmatrix_), dgprec(dgprec_), dgcc(dgcc_),
339  cggfs(cggfs_), cgprec(cgprec_), cgcc(cgcc_), p(p_), dghelper(dghelper_),
340  comm(comm_), n1(n1_), n2(n2_)
341  {
342  }
343 
349  virtual void pre (V& x, W& b)
350  {
352  CGW cgd(cggfs,0.0);
353  CGV cgv(cggfs,0.0);
354  cgprec.pre(Dune::PDELab::istl::raw(cgv),Dune::PDELab::istl::raw(cgd));
355  }
356 
362  virtual void apply (V& x, const W& b)
363  {
364  // need local copies to store defect and solution
365  W d(b);
367  V v(x); // only to get correct size
368 
369  // pre-smoothing on DG matrix
370  for (int i=0; i<n1; i++)
371  {
372  v = 0.0;
375  if (dggfs.gridView().comm().size()>1)
376  dggfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
379  x += v;
380  }
381 
382  // restrict defect to CG subspace
383  dghelper.maskForeignDOFs(d); // DG defect is additive for overlap 1, but in case we use more
384  CGW cgd(cggfs,0.0);
387  if (cggfs.gridView().comm().size()>1)
388  cggfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication); // now we have consistent defect on coarse grid
390  comm.project(Dune::PDELab::istl::raw(cgd));
391  CGV cgv(cggfs,0.0);
392 
393 
394  // call preconditioner
395  cgprec.apply(Dune::PDELab::istl::raw(cgv),Dune::PDELab::istl::raw(cgd));
396 
397  // prolongate correction
401  x += v;
402 
403  // post-smoothing on DG matrix
404  for (int i=0; i<n2; i++)
405  {
406  v = 0.0;
409  if (dggfs.gridView().comm().size()>1)
410  dggfs.gridView().communicate(adddh,Dune::All_All_Interface,Dune::ForwardCommunication);
413  x += v;
414  }
415  }
416 
422  virtual void post (V& x)
423  {
424  dgprec.post(Dune::PDELab::istl::raw(x));
425  CGV cgv(cggfs,0.0);
426  cgprec.post(Dune::PDELab::istl::raw(cgv));
427  }
428  };
429 
430 //***********************************************************
431 // The DG/AMG linear solver backend in the overlapping case
432 //***********************************************************
433 
446 template<class DGGO, class DGCC, class CGGFS, class CGCC, class TransferLOP,
447  template<class,class,class,int> class DGPrec, template<class> class Solver, int s=96>
449  public Dune::PDELab::OVLPScalarProductImplementation<typename DGGO::Traits::TrialGridFunctionSpace>,
451 {
452 public:
453  // DG grid function space
454  typedef typename DGGO::Traits::TrialGridFunctionSpace GFS;
455 
456  // vectors and matrices on DG level
457  typedef typename DGGO::Traits::Jacobian M; // wrapped istl DG matrix
458  typedef typename DGGO::Traits::Domain V; // wrapped istl DG vector
459  typedef typename M::BaseT Matrix; // istl DG matrix
460  typedef typename V::BaseT Vector; // istl DG vector
461  typedef typename Vector::field_type field_type;
462 
463  // vectors and matrices on CG level
464  typedef typename Dune::PDELab::BackendVectorSelector<CGGFS,field_type>::Type CGV; // wrapped istl CG vector
465  typedef typename CGV::BaseT CGVector; // istl CG vector
466 
467  // prolongation matrix
470  typedef TransferLOP CGTODGLOP; // local operator
472  typedef typename PGO::Jacobian PMatrix; // wrapped ISTL prolongation matrix
473  typedef typename PMatrix::BaseT P; // ISTL prolongation matrix
474 
475  // CG subspace matrix
476  typedef typename Dune::TransposedMatMultMatResult<P,Matrix>::type PTADG;
477  typedef typename Dune::MatMultMatResult<PTADG,P>::type CGMatrix; // istl coarse space matrix
478 
479 private:
480 
481  const GFS& gfs;
482  DGGO& dggo;
483  const DGCC& dgcc;
484  CGGFS& cggfs;
485  const CGCC& cgcc;
486  unsigned maxiter;
487  int verbose;
488  bool usesuperlu;
489 
490  CGTODGLOP cgtodglop; // local operator to assemble prolongation matrix
491  PGO pgo; // grid operator to assemble prolongation matrix
492  PMatrix pmatrix; // wrapped prolongation matrix
493 
496  class EmptyLop : public Dune::PDELab::NumericalJacobianApplyVolume<EmptyLop>,
500  {
501  };
502 
503 public:
504 
505  // access to prolongation matrix
507  {
508  return pmatrix;
509  }
510 
513  ISTLBackend_OVLP_AMG_4_DG(DGGO& dggo_, const DGCC& dgcc_, CGGFS& cggfs_, const CGCC& cgcc_,
514  unsigned maxiter_=5000, int verbose_=1, bool usesuperlu_=true) :
515  Dune::PDELab::OVLPScalarProductImplementation<typename DGGO::Traits::TrialGridFunctionSpace>(dggo_.trialGridFunctionSpace()),
516  gfs(dggo_.trialGridFunctionSpace()), dggo(dggo_), dgcc(dgcc_), cggfs(cggfs_), cgcc(cgcc_), maxiter(maxiter_), verbose(verbose_), usesuperlu(usesuperlu_),
517  cgtodglop(), pgo(cggfs,dggo.trialGridFunctionSpace(),cgtodglop), pmatrix(pgo)
518  {
519 #if !HAVE_SUPERLU
520  if (usesuperlu == true)
521  {
522  if (gfs.gridView().comm().rank()==0)
523  std::cout << "WARNING: You are using AMG without SuperLU!"
524  << " Please consider installing SuperLU,"
525  << " or set the usesuperlu flag to false"
526  << " to suppress this warning." << std::endl;
527  }
528 #endif
529 
530  // assemble prolongation matrix; this will not change from one apply to the next
531  pmatrix = 0.0;
532  if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout << "allocated prolongation matrix of size " << pmatrix.N() << " x " << pmatrix.M() << std::endl;
533  CGV cgx(cggfs,0.0); // need vector to call jacobian
534  pgo.jacobian(cgx,pmatrix);
535  }
536 
544  void apply (M& A, V& z, V& r, typename V::ElementType reduction)
545  {
546  // make operator and scalar product for overlapping solver
548  POP pop(dgcc,A);
550  PSP psp(*this);
551 
552  // compute CG matrix
553  // make grid operator with empty local operator => matrix data type and constraints assembly
554  EmptyLop emptylop;
556  CGGO cggo(cggfs,cgcc,cggfs,cgcc,emptylop);
557  typedef typename CGGO::Jacobian CGM;
558 
559  // do triple matrix product ACG = P^T ADG P; this is purely local
560  Dune::Timer watch;
561  watch.reset();
562  tags::attached_container attached_container;
563  CGM acg(attached_container);
564  {
565  PTADG ptadg;
566  Dune::transposeMatMultMat(ptadg,Dune::PDELab::istl::raw(pmatrix),Dune::PDELab::istl::raw(A)); // 1a
567  //Dune::transposeMatMultMat(ptadg,Dune::PDELab::istl::raw(pmatrix),Dune::PDELab::istl::raw(A2)); // 1b
568  Dune::matMultMat(Dune::PDELab::istl::raw(acg),ptadg,Dune::PDELab::istl::raw(pmatrix));
569  }
570  double triple_product_time = watch.elapsed();
571  if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout << "=== triple matrix product " << triple_product_time << " s" << std::endl;
572  //Dune::printmatrix(std::cout,Dune::PDELab::istl::raw(acg),"triple product matrix","row",10,2);
573  CGV cgx(cggfs,0.0); // need vector to call jacobian
574  cggo.jacobian(cgx,acg); // insert trivial rows at processor boundaries
575  //std::cout << "CG constraints: " << cgcc.size() << " out of " << cggfs.globalSize() << std::endl;
576 
577  // NOW we need to insert the processor boundary conditions in DG matrix
579  DGGOEmpty dggoempty(gfs,dgcc,gfs,dgcc,emptylop);
580  dggoempty.jacobian(z,A);
581 
582  // and in the residual
584 
585  // now set up parallel AMG solver for the CG subspace
587  Comm oocc(gfs.gridView().comm());
589  CGHELPER cghelper(cggfs,2);
590  cghelper.createIndexSetAndProjectForAMG(acg,oocc);
591  typedef Dune::OverlappingSchwarzOperator<CGMatrix,CGVector,CGVector,Comm> ParCGOperator;
592  ParCGOperator paroop(Dune::PDELab::istl::raw(acg),oocc);
593  Dune::OverlappingSchwarzScalarProduct<CGVector,Comm> sp(oocc);
594  typedef Dune::Amg::Parameters Parameters; // AMG parameters (might be nice to change from outside)
595  Parameters params(15,2000);
596  params.setDefaultValuesIsotropic(CGGFS::Traits::GridViewType::Traits::Grid::dimension);
597  params.setDebugLevel(verbose);
598  params.setCoarsenTarget(2000);
599  params.setMaxLevel(20);
600  params.setProlongationDampingFactor(1.6);
601  params.setNoPreSmoothSteps(3);
602  params.setNoPostSmoothSteps(3);
603  params.setGamma(1);
604  params.setAdditive(false);
605  //params.setAccumulate(Dune::Amg::AccumulationMode::noAccu); // atOnceAccu results in deadlock
606  typedef Dune::SeqSSOR<CGMatrix,CGVector,CGVector,1> Smoother;
607  typedef Dune::BlockPreconditioner<CGVector,CGVector,Comm,Smoother> ParSmoother;
608  typedef typename Dune::Amg::SmootherTraits<ParSmoother>::Arguments SmootherArgs;
609  SmootherArgs smootherArgs;
610  smootherArgs.iterations = 2;
611  smootherArgs.relaxationFactor = 0.92;
612  typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<CGMatrix,Dune::Amg::FirstDiagonal> > Criterion;
613  Criterion criterion(params);
614  typedef Dune::Amg::AMG<ParCGOperator,CGVector,ParSmoother,Comm> AMG;
615  watch.reset();
616  AMG amg(paroop,criterion,smootherArgs,oocc);
617  double amg_setup_time = watch.elapsed();
618  if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout << "=== AMG setup " <<amg_setup_time << " s" << std::endl;
619 
620  // set up hybrid DG/CG preconditioner
621  typedef DGPrec<Matrix,Vector,Vector,1> DGPrecType;
622  DGPrecType dgprec(Dune::PDELab::istl::raw(A),1,0.92);
623  //DGPrecType dgprec(Dune::PDELab::istl::raw(A),0.92);
626  HybridPrec hybridprec(gfs,Dune::PDELab::istl::raw(A),dgprec,dgcc,cggfs,amg,cgcc,Dune::PDELab::istl::raw(pmatrix),
627  this->parallelHelper(),oocc,3,3);
628 
629  // /********/
630  // /* Test */
631  // /********/
632  // Solver<CGVector> testsolver(paroop,sp,amg,1e-8,100,2);
633  // CGV cgxx(cggfs,0.0);
634  // CGV cgdd(cggfs,1.0);
635  // Dune::InverseOperatorResult statstat;
636  // testsolver.apply(Dune::PDELab::istl::raw(cgxx),Dune::PDELab::istl::raw(cgdd),statstat);
637  // /********/
638 
639  // set up solver
640  int verb=verbose;
641  if (gfs.gridView().comm().rank()>0) verb=0;
642  Solver<V> solver(pop,psp,hybridprec,reduction,maxiter,verb);
643 
644  // solve
645  Dune::InverseOperatorResult stat;
646  watch.reset();
647  solver.apply(z,r,stat);
648  double amg_solve_time = watch.elapsed();
649  if (verbose>0 && gfs.gridView().comm().rank()==0) std::cout << "=== Hybrid total solve time " << amg_solve_time+amg_setup_time+triple_product_time << " s" << std::endl;
650  res.converged = stat.converged;
651  res.iterations = stat.iterations;
652  res.elapsed = amg_solve_time+amg_setup_time+triple_product_time;
653  res.reduction = stat.reduction;
654  res.conv_rate = stat.conv_rate;
655  }
656 
657 };
658 }
659 }
660 #endif
void scatter(MessageBuffer &buff, const EntityType &e, size_t n)
Definition: ovlp_amg_dg_backend.hh:88
std::pair< IdType, B > DataType
export type of data for message buffer
Definition: ovlp_amg_dg_backend.hh:139
GFS::Traits::GridView GV
Definition: ovlp_amg_dg_backend.hh:33
GlobalIdSet::IdType IdType
Definition: ovlp_amg_dg_backend.hh:133
The category the preconditioner is part of.
Definition: ovlp_amg_dg_backend.hh:325
PMatrix::BaseT P
Definition: ovlp_amg_dg_backend.hh:473
Definition: ovlpistlsolverbackend.hh:41
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:172
Definition: ovlp_amg_dg_backend.hh:28
PGO::Jacobian PMatrix
Definition: ovlp_amg_dg_backend.hh:472
Dune::PDELab::GridOperator< CGGFS, GFS, CGTODGLOP, MBE, field_type, field_type, field_type, CC, CC > PGO
Definition: ovlp_amg_dg_backend.hh:471
Definition: ovlp_amg_dg_backend.hh:448
W::BaseT Y
Definition: ovlp_amg_dg_backend.hh:318
void apply(M &A, V &z, V &r, typename V::ElementType reduction)
solve the given linear system
Definition: ovlp_amg_dg_backend.hh:544
bool converged
Definition: solver.hh:32
Grid::Traits::GlobalIdSet GlobalIdSet
Definition: ovlp_amg_dg_backend.hh:132
Tag for requesting a vector or matrix container with a pre-attached underlying object.
Definition: backend/tags.hh:27
void gather(MessageBuffer &buff, const EntityType &e) const
pack data from user to message buffer
Definition: ovlp_amg_dg_backend.hh:182
Default flags for all local operators.
Definition: flags.hh:18
bool fixedsize(int dim, int codim) const
returns true if size per entity of given dim and codim is a constant
Definition: ovlp_amg_dg_backend.hh:152
size_t size(EntityType &e) const
Definition: ovlp_amg_dg_backend.hh:64
Dune::PDELab::BackendVectorSelector< CGGFS, typename CGPrec::range_type::field_type >::Type CGW
Definition: ovlp_amg_dg_backend.hh:320
double elapsed
Definition: solver.hh:34
std::map< IdType, IndexType > GlobalToLocalMap
Definition: ovlp_amg_dg_backend.hh:45
Dune::TransposedMatMultMatResult< P, Matrix >::type PTADG
Definition: ovlp_amg_dg_backend.hh:476
Dune::PDELab::BackendVectorSelector< CGGFS, field_type >::Type CGV
Definition: ovlp_amg_dg_backend.hh:464
Dune::PDELab::EmptyTransformation CC
Definition: ovlp_amg_dg_backend.hh:469
Definition: common/constraintstransformation.hh:127
void set_constrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
construct constraints from given boundary condition function
Definition: common/constraints.hh:772
Definition: ovlp_amg_dg_backend.hh:121
Backend using ISTL matrices.
Definition: istl/descriptors.hh:69
TransferLOP CGTODGLOP
Definition: ovlp_amg_dg_backend.hh:470
Definition: ovlp_amg_dg_backend.hh:298
Dune::PDELab::ISTLMatrixBackend MBE
Definition: ovlp_amg_dg_backend.hh:468
PMatrix & prolongation_matrix()
Definition: ovlp_amg_dg_backend.hh:506
static const int dim
Definition: adaptivity.hh:82
GV::Grid Grid
Definition: ovlp_amg_dg_backend.hh:131
V::BaseT X
Definition: ovlp_amg_dg_backend.hh:317
void gather(MessageBuffer &buff, const EntityType &e) const
pack data from user to message buffer
Definition: ovlp_amg_dg_backend.hh:71
bool contains(int dim, int codim) const
returns true if data for this codim should be communicated
Definition: ovlp_amg_dg_backend.hh:146
Dune::PDELab::BackendVectorSelector< DGGFS, typename DGPrec::domain_type::field_type >::Type V
Definition: ovlp_amg_dg_backend.hh:315
Grid::Traits::GlobalIdSet GlobalIdSet
Definition: ovlp_amg_dg_backend.hh:37
IndexSet::IndexType IndexType
Definition: ovlp_amg_dg_backend.hh:35
V & raw(V &v)
Returns the raw ISTL object associated with v, or v itself it is already an ISTL object.
Definition: backend/istl/utility.hh:26
Dune::PDELab::LinearSolverResult< double > res
Definition: solver.hh:52
const istl::ParallelHelper< DGGO::Traits::TrialGridFunctionSpace > & parallelHelper() const
Definition: ovlpistlsolverbackend.hh:352
std::map< IdType, IndexType > GlobalToLocalMap
Definition: ovlp_amg_dg_backend.hh:143
OvlpDGAMGPrec(const DGGFS &dggfs_, DGMatrix &dgmatrix_, DGPrec &dgprec_, const DGCC &dgcc_, const CGGFS &cggfs_, CGPrec &cgprec_, const CGCC &cgcc_, P &p_, const DGHelper &dghelper_, const Comm &comm_, int n1_, int n2_)
Constructor.
Definition: ovlp_amg_dg_backend.hh:335
GV::IndexSet IndexSet
Definition: ovlp_amg_dg_backend.hh:34
virtual void pre(V &x, W &b)
Prepare the preconditioner.
Definition: ovlp_amg_dg_backend.hh:349
M::BaseT Matrix
Definition: ovlp_amg_dg_backend.hh:459
Dune::MatMultMatResult< PTADG, P >::type CGMatrix
Definition: ovlp_amg_dg_backend.hh:477
CGV::BaseT CGVector
Definition: ovlp_amg_dg_backend.hh:465
Definition: parallelhelper.hh:45
std::map< IndexType, IdType > LocalToGlobalMap
Definition: ovlp_amg_dg_backend.hh:142
Definition: genericdatahandle.hh:623
V::BaseT Vector
Definition: ovlp_amg_dg_backend.hh:460
virtual void post(V &x)
Clean up.
Definition: ovlp_amg_dg_backend.hh:422
LocalGlobalMapDataHandle(const GFS &gfs_, LocalToGlobalMap &l2g_, GlobalToLocalMap &g2l_)
constructor
Definition: ovlp_amg_dg_backend.hh:101
GV::Grid Grid
Definition: ovlp_amg_dg_backend.hh:36
RFType reduction
Definition: solver.hh:35
Dune::PDELab::BackendVectorSelector< CGGFS, typename CGPrec::domain_type::field_type >::Type CGV
Definition: ovlp_amg_dg_backend.hh:319
M::block_type B
Definition: ovlp_amg_dg_backend.hh:136
BackendVectorSelectorHelper< Backend, GridFunctionSpace, FieldType >::Type Type
Definition: backendselector.hh:14
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
Definition: solver.hh:42
void restore_overlap_entries(const GFS &gfs, Dune::BCRSMatrix< Dune::FieldMatrix< T, n, m >, A > &matrix, Dune::BCRSMatrix< Dune::FieldMatrix< T, n, m >, A > &matrix2)
Definition: ovlp_amg_dg_backend.hh:261
DGGO::Traits::Jacobian M
Definition: ovlp_amg_dg_backend.hh:457
Definition: ovlpistlsolverbackend.hh:322
bool contains(int dim, int codim) const
returns true if data for this codim should be communicated
Definition: ovlp_amg_dg_backend.hh:48
sparsity pattern generator
Definition: pattern.hh:14
Implement jacobian_apply_volume() based on alpha_volume()
Definition: defaultimp.hh:321
bool fixedsize(int dim, int codim) const
returns true if size per entity of given dim and codim is a constant
Definition: ovlp_amg_dg_backend.hh:54
virtual void apply(V &x, const W &b)
Apply the precondioner.
Definition: ovlp_amg_dg_backend.hh:362
GV::IndexSet IndexSet
Definition: ovlp_amg_dg_backend.hh:129
size_t size(EntityType &e) const
Definition: ovlp_amg_dg_backend.hh:162
int DataType
export type of data for message buffer
Definition: ovlp_amg_dg_backend.hh:41
Dune::PDELab::BackendVectorSelector< DGGFS, typename DGPrec::range_type::field_type >::Type W
Definition: ovlp_amg_dg_backend.hh:316
Dune::Amg::SequentialInformation type
Definition: parallelhelper.hh:417
std::map< IndexType, IdType > LocalToGlobalMap
Definition: ovlp_amg_dg_backend.hh:44
GlobalIdSet::IdType IdType
Definition: ovlp_amg_dg_backend.hh:38
DGGO::Traits::Domain V
Definition: ovlp_amg_dg_backend.hh:458
Vector::field_type field_type
Definition: ovlp_amg_dg_backend.hh:461
void scatter(MessageBuffer &buff, const EntityType &e, size_t n)
Definition: ovlp_amg_dg_backend.hh:207
const E & e
Definition: interpolate.hh:172
DGGO::Traits::TrialGridFunctionSpace GFS
Definition: ovlp_amg_dg_backend.hh:454
MatrixExchangeDataHandle(const GFS &gfs_, M &m_, const LocalToGlobalMap &l2g_, const GlobalToLocalMap &g2l_, M &m2_)
constructor
Definition: ovlp_amg_dg_backend.hh:238
IndexSet::IndexType IndexType
Definition: ovlp_amg_dg_backend.hh:130
GFS::Traits::GridView GV
Definition: ovlp_amg_dg_backend.hh:128
RFType conv_rate
Definition: solver.hh:36
Definition: ovlpistlsolverbackend.hh:370
Dune::PDELab::BackendMatrixSelector< MBE, Domain, Range, field_type >::Type Jacobian
The type of the jacobian.
Definition: gridoperator.hh:46
const std::string s
Definition: function.hh:1103
ISTLBackend_OVLP_AMG_4_DG(DGGO &dggo_, const DGCC &dgcc_, CGGFS &cggfs_, const CGCC &cgcc_, unsigned maxiter_=5000, int verbose_=1, bool usesuperlu_=true)
Definition: ovlp_amg_dg_backend.hh:513
unsigned int iterations
Definition: solver.hh:33