dune-pdelab  2.5-dev
pdelab.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 #ifndef DUNE_PDELAB_BOILERPLATE_PDELAB_HH
4 #define DUNE_PDELAB_BOILERPLATE_PDELAB_HH
5 
16 // first of all we include a lot of dune grids and pdelab files
17 #include <iostream>
18 #include <memory>
19 
20 #include <dune/common/parallel/mpihelper.hh> // include mpi helper class
21 #include <dune/common/parametertreeparser.hh>
22 #include <dune/common/classname.hh>
23 #include <dune/common/exceptions.hh>
24 #include <dune/common/fvector.hh>
25 
26 #include <dune/geometry/type.hh>
27 #include <dune/geometry/quadraturerules.hh>
28 
29 #include <dune/grid/onedgrid.hh>
30 #include <dune/grid/io/file/vtk.hh>
31 #include <dune/grid/yaspgrid.hh>
32 #if HAVE_UG
33 #include <dune/grid/uggrid.hh>
34 #endif
35 #if HAVE_ALBERTA
36 #include<dune/grid/albertagrid.hh>
37 #include <dune/grid/albertagrid/dgfparser.hh>
38 #endif
39 #if HAVE_UG
40 #include<dune/grid/uggrid.hh>
41 #endif
42 #if HAVE_DUNE_ALUGRID
43 #include<dune/alugrid/grid.hh>
44 #include <dune/alugrid/dgf.hh>
45 #endif
46 #include <dune/grid/utility/structuredgridfactory.hh>
47 #include <dune/grid/io/file/gmshreader.hh>
48 
49 #include <dune/istl/bvector.hh>
50 #include <dune/istl/operators.hh>
51 #include <dune/istl/solvers.hh>
52 #include <dune/istl/solvercategory.hh>
53 #include <dune/istl/preconditioners.hh>
54 #include <dune/istl/io.hh>
55 
56 #include <dune/istl/paamg/amg.hh>
82 
83 namespace Dune {
84  namespace PDELab {
85 
86  // make grids
87  template<typename T>
89  {
90  public:
91  // export types
92  typedef T Grid;
93  typedef typename T::ctype ctype;
94  static const int dim = T::dimension;
95  static const int dimworld = T::dimensionworld;
96 
97  // constructors
98  StructuredGrid (Dune::GeometryType::BasicType meshtype, unsigned int cells)
99  {
100  FieldVector<ctype,dimworld> lowerLeft(0.0);
101  FieldVector<ctype,dimworld> upperRight(1.0);
102  array<unsigned int,dim> elements; elements.fill(cells);
103 
104  StructuredGridFactory<T> factory;
105 
106  if (meshtype==Dune::GeometryType::cube)
107  gridp = factory.createCubeGrid(lowerLeft,upperRight,elements);
108  else if (meshtype==Dune::GeometryType::simplex)
109  gridp = factory.createSimplexGrid(lowerLeft,upperRight,elements);
110  else
111  {
112  DUNE_THROW(GridError, className<StructuredGrid>()
113  << "::StructuredGrid(): grid type must be simplex or cube ");
114  }
115  }
116 
117 
118  StructuredGrid (Dune::GeometryType::BasicType meshtype,
119  array<double,dimworld> lower_left, array<double,dimworld> upper_right,
120  array<unsigned int,dim> cells)
121  {
122  FieldVector<ctype,dimworld> lowerLeft;
123  FieldVector<ctype,dimworld> upperRight;
124  array<unsigned int,dim> elements;
125 
126  // copy data to correct types for StructuredGridFactory
127  for (size_t i=0; i<dimworld; i++)
128  {
129  lowerLeft[i] = lower_left[i];
130  upperRight[i] = upper_right[i];
131  }
132  for (size_t i=0; i<dim; i++)
133  {
134  elements[i] = cells[i];
135  }
136 
137  StructuredGridFactory<T> factory;
138 
139  if (meshtype==Dune::GeometryType::cube)
140  gridp = factory.createCubeGrid(lowerLeft,upperRight,elements);
141  else if (meshtype==Dune::GeometryType::simplex)
142  gridp = factory.createSimplexGrid(lowerLeft,upperRight,elements);
143  else
144  {
145  DUNE_THROW(GridError, className<StructuredGrid>()
146  << "::StructuredGrid(): grid type must be simplex or cube ");
147  }
148  }
149 
150  // return shared pointer
151  std::shared_ptr<T> getSharedPtr ()
152  {
153  return gridp;
154  }
155 
156  // return grid reference
157  T& getGrid ()
158  {
159  return *gridp;
160  }
161 
162  // return grid reference const version
163  const T& getGrid () const
164  {
165  return *gridp;
166  }
167 
169  {
170  return *gridp;
171  }
172 
174  {
175  return gridp.operator->();
176  }
177 
178  const T& operator*() const
179  {
180  return *gridp;
181  }
182 
183  const T* operator->() const
184  {
185  return gridp.operator->();
186  }
187 
188 
189  private:
190  std::shared_ptr<T> gridp; // hold a shared pointer to a grid
191  };
192 
193  // specialization for yaspgrid; treats paralle case right
194  template<int dim>
195  class StructuredGrid<YaspGrid<dim> >
196  {
197  public:
198 
199  // export types
200  typedef YaspGrid<dim> Grid;
201  typedef typename Grid::ctype ctype;
202  static const int dimworld = Grid::dimensionworld;
203 
204  // simple constructor for the unit cube
205  StructuredGrid (Dune::GeometryType::BasicType meshtype, unsigned int cells, int overlap=1)
206  {
207  // check element type
208  if (meshtype!=Dune::GeometryType::cube)
209  std::cout << "StructuredGrid(): element type " << meshtype << " is ignored" << std::endl;
210 
211  // copy data to correct types for YaspGrid
212  Dune::FieldVector<double,dimworld> L(1.0);
213  std::array<int,dimworld> N(Dune::fill_array<int,dimworld>(cells));
214  std::bitset<dimworld> B(false);
215 
216  // instantiate the grid
217  gridp = std::shared_ptr<Grid>(new Grid(L,N,B,overlap,Dune::MPIHelper::getCollectiveCommunication()));
218  }
219 
220  // constructor with sizes given
221  StructuredGrid (Dune::GeometryType::BasicType meshtype,
222  array<double,dimworld> lower_left, array<double,dimworld> upper_right,
223  array<unsigned int,dim> cells, int overlap=1)
224  {
225  // check that lower right corner is the origin
226  for(int d = 0; d < dimworld; ++d)
227  if(std::abs(lower_left[d]) > std::abs(upper_right[d])*1e-10)
228  DUNE_THROW(GridError, className<StructuredGrid>()
229  << "::createCubeGrid(): The lower coordinates "
230  "must be at the origin for YaspGrid.");
231 
232  // check element type
233  if (meshtype!=Dune::GeometryType::cube)
234  std::cout << "StructuredGrid(): element type " << meshtype << " is ignored" << std::endl;
235 
236  // copy data to correct types for YaspGrid
237  Dune::FieldVector<double,dimworld> L;
238  std::array<int,dimworld> N;
239  std::bitset<dimworld> B(false);
240  for (size_t i=0; i<dimworld; i++)
241  {
242  L[i] = upper_right[i];
243  N[i] = cells[i];
244  }
245 
246  // instantiate the grid
247  gridp = std::shared_ptr<Grid>(new Grid(L,N,B,overlap,Dune::MPIHelper::getCollectiveCommunication()));
248  }
249 
250  // constructor with periodicity argument
251  StructuredGrid (Dune::GeometryType::BasicType meshtype,
252  array<double,dimworld> lower_left, array<double,dimworld> upper_right,
253  array<unsigned int,dim> cells, array<bool,dim> periodic, int overlap=1)
254  {
255  // check that lower right corner is the origin
256  for(int d = 0; d < dimworld; ++d)
257  if(std::abs(lower_left[d]) > std::abs(upper_right[d])*1e-10)
258  DUNE_THROW(GridError, className<StructuredGrid>()
259  << "::createCubeGrid(): The lower coordinates "
260  "must be at the origin for YaspGrid.");
261 
262  // check element type
263  if (meshtype!=Dune::GeometryType::cube)
264  std::cout << "StructuredGrid(): element type " << meshtype << " is ignored" << std::endl;
265 
266  // copy data to correct types for YaspGrid
267  Dune::FieldVector<double,dimworld> L;
268  std::array<int,dimworld> N;
269  std::bitset<dimworld> B(false);
270  for (size_t i=0; i<dimworld; i++)
271  {
272  L[i] = upper_right[i];
273  N[i] = cells[i];
274  B[i] = periodic[i];
275  }
276 
277  // instantiate the grid
278  gridp = std::shared_ptr<Grid>(new Grid(L,N,B,overlap,Dune::MPIHelper::getCollectiveCommunication()));
279  }
280 
281  // return shared pointer
282  std::shared_ptr<Grid> getSharedPtr ()
283  {
284  return gridp;
285  }
286 
287  // return grid reference
288  Grid& getGrid ()
289  {
290  return *gridp;
291  }
292 
293  // return grid reference const version
294  const Grid& getGrid () const
295  {
296  return *gridp;
297  }
298 
299  Grid& operator*()
300  {
301  return *gridp;
302  }
303 
304  Grid* operator->()
305  {
306  return gridp.operator->();
307  }
308 
309  const Grid& operator*() const
310  {
311  return *gridp;
312  }
313 
314  const Grid* operator->() const
315  {
316  return gridp.operator->();
317  }
318 
319  private:
320  std::shared_ptr<Grid> gridp; // hold a shared pointer to a grid
321  };
322 
323  // unstructured grid read from gmsh file
324  template<typename T>
326  {
327  public:
328  // export types
329  typedef T Grid;
330  typedef typename T::ctype ctype;
331  static const int dim = T::dimension;
332  static const int dimworld = T::dimensionworld;
333 
334  // constructors
335  UnstructuredGrid (std::string filename, bool verbose = true, bool insert_boundary_segments=true)
336  {
337  Dune::GridFactory<T> factory;
338  Dune::GmshReader<T>::read(factory,filename,verbose,insert_boundary_segments);
339  gridp = std::shared_ptr<T>(factory.createGrid());
340  }
341 
342  // return shared pointer
343  std::shared_ptr<T> getSharedPtr ()
344  {
345  return gridp;
346  }
347 
348  // return grid reference
349  T& getGrid ()
350  {
351  return *gridp;
352  }
353 
354  // return grid reference const version
355  const T& getGrid () const
356  {
357  return *gridp;
358  }
359 
361  {
362  return *gridp;
363  }
364 
366  {
367  return gridp.operator->();
368  }
369 
370  const T& operator*() const
371  {
372  return *gridp;
373  }
374 
375  const T* operator->() const
376  {
377  return gridp.operator->();
378  }
379 
380  private:
381  std::shared_ptr<T> gridp; // hold a shared pointer to a grid
382  };
383 
384 
385  //============================================================================
386  // Continuous Lagrange Finite Element Space
387  //============================================================================
388 
389  // finite element map base template
390  template<typename GV, typename C, typename R, unsigned int degree, unsigned int dim, Dune::GeometryType::BasicType gt>
391  class CGFEMBase
392  {};
393 
394  template<typename GV, typename C, typename R, unsigned int degree, unsigned int dim>
395  class CGFEMBase<GV,C,R,degree,dim,Dune::GeometryType::simplex>
396  {
397  public:
399 
400  CGFEMBase (const GV& gridview)
401  {
402  femp = std::shared_ptr<FEM>(new FEM(gridview));
403  }
404 
405  FEM& getFEM() {return *femp;}
406  const FEM& getFEM() const {return *femp;}
407 
408  private:
409  std::shared_ptr<FEM> femp;
410  };
411 
412  template<typename GV, typename C, typename R, unsigned int degree, unsigned int dim>
413  class CGFEMBase<GV,C,R,degree,dim,Dune::GeometryType::cube>
414  {
415  public:
417 
418  CGFEMBase (const GV& gridview)
419  {
420  femp = std::shared_ptr<FEM>(new FEM(gridview));
421  }
422 
423  FEM& getFEM() {return *femp;}
424  const FEM& getFEM() const {return *femp;}
425 
426  private:
427  std::shared_ptr<FEM> femp;
428  };
429 
430  //============================================================================
431 
432  // define enumeration type that differentiate conforming and nonconforming meshes
433  enum MeshType {
436  };
437 
438  // constraints base template
439  template<typename Grid, unsigned int degree, Dune::GeometryType::BasicType gt, MeshType mt, SolverCategory::Category st, typename BCType, typename GV = typename Grid::LeafGridView>
440  class CGCONBase
441  {};
442 
443  template<typename Grid, typename BCType, typename GV>
444  class CGCONBase<Grid,1,Dune::GeometryType::simplex,MeshType::nonconforming,SolverCategory::sequential,BCType,GV>
445  {
446  public:
448 
449  CGCONBase (Grid& grid, const BCType& bctype, const GV& gv)
450  {
451  conp = std::shared_ptr<CON>(new CON(grid,true,bctype));
452  }
453 
454  CGCONBase (Grid& grid, const BCType& bctype)
455  {
456  conp = std::shared_ptr<CON>(new CON(grid,true,bctype));
457  }
458 
459  template<typename GFS>
460  void postGFSHook (const GFS& gfs) {}
461  CON& getCON() {return *conp;}
462  const CON& getCON() const {return *conp;}
463  template<typename GFS, typename DOF>
464  void make_consistent (const GFS& gfs, DOF& x) const {}
465  private:
466  std::shared_ptr<CON> conp;
467  };
468 
469  template<typename Grid, typename BCType, typename GV>
470  class CGCONBase<Grid,1,Dune::GeometryType::cube,MeshType::nonconforming,SolverCategory::sequential,BCType,GV>
471  {
472  public:
474 
475  CGCONBase (Grid& grid, const BCType& bctype, const GV& gv)
476  {
477  conp = std::shared_ptr<CON>(new CON(grid,true,bctype));
478  }
479 
480  CGCONBase (Grid& grid, const BCType& bctype)
481  {
482  conp = std::shared_ptr<CON>(new CON(grid,true,bctype));
483  }
484 
485  template<typename GFS>
486  void postGFSHook (const GFS& gfs) {}
487  CON& getCON() {return *conp;}
488  const CON& getCON() const {return *conp;}
489  template<typename GFS, typename DOF>
490  void make_consistent (const GFS& gfs, DOF& x) const {}
491  private:
492  std::shared_ptr<CON> conp;
493  };
494 
495  template<typename Grid, unsigned int degree, Dune::GeometryType::BasicType gt,typename BCType, typename GV>
496  class CGCONBase<Grid,degree,gt,MeshType::conforming,SolverCategory::sequential,BCType,GV>
497  {
498  public:
500 
501  CGCONBase (Grid& grid, const BCType& bctype, const GV& gv)
502  {
503  conp = std::shared_ptr<CON>(new CON());
504  }
505 
506  CGCONBase (Grid& grid, const BCType& bctype)
507  {
508  conp = std::shared_ptr<CON>(new CON());
509  }
510 
511  template<typename GFS>
512  void postGFSHook (const GFS& gfs) {}
513  CON& getCON() {return *conp;}
514  const CON& getCON() const {return *conp;}
515  template<typename GFS, typename DOF>
516  void make_consistent (const GFS& gfs, DOF& x) const {}
517  private:
518  std::shared_ptr<CON> conp;
519  };
520 
521  template<typename Grid, unsigned int degree, Dune::GeometryType::BasicType gt,typename BCType, typename GV>
522  class CGCONBase<Grid,degree,gt,MeshType::conforming,SolverCategory::overlapping,BCType,GV>
523  {
524  public:
526 
527  CGCONBase (Grid& grid, const BCType& bctype, const GV& gv)
528  {
529  conp = std::shared_ptr<CON>(new CON());
530  }
531 
532  CGCONBase (Grid& grid, const BCType& bctype)
533  {
534  conp = std::shared_ptr<CON>(new CON());
535  }
536 
537  template<typename GFS>
538  void postGFSHook (const GFS& gfs) {}
539  CON& getCON() {return *conp;}
540  const CON& getCON() const {return *conp;}
541  template<typename GFS, typename DOF>
542  void make_consistent (const GFS& gfs, DOF& x) const
543  {
544  // make vector consistent; this is needed for all overlapping solvers
545  ISTL::ParallelHelper<GFS> helper(gfs);
546  helper.maskForeignDOFs(Backend::native(x));
548  if (gfs.gridView().comm().size()>1)
549  gfs.gridView().communicate(adddh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
550  }
551  private:
552  std::shared_ptr<CON> conp;
553  };
554 
555  template<typename Grid, unsigned int degree, Dune::GeometryType::BasicType gt,typename BCType, typename GV>
556  class CGCONBase<Grid,degree,gt,MeshType::conforming,SolverCategory::nonoverlapping,BCType,GV>
557  {
558  public:
560  CGCONBase (Grid& grid, const BCType& bctype)
561  {
562  conp = std::shared_ptr<CON>(new CON);
563  }
564 
565  template<typename GFS>
566  CON& getCON() {return *conp;}
567  const CON& getCON() const {return *conp;}
568  template<typename GFS, typename DOF>
569  void make_consistent (const GFS& gfs, DOF& x) const {}
570  private:
571  std::shared_ptr<CON> conp;
572  };
573 
574 
575  // continuous Lagrange finite elements
576  template<typename T, typename N, unsigned int degree, typename BCType,
577  Dune::GeometryType::BasicType gt, MeshType mt, SolverCategory::Category st = SolverCategory::sequential,
578  typename VBET=ISTL::VectorBackend<> >
579  class CGSpace {
580  public:
581 
582  // export types
583  typedef T Grid;
584  typedef typename T::LeafGridView GV;
585  typedef typename T::ctype ctype;
586  static const int dim = T::dimension;
587  static const int dimworld = T::dimensionworld;
588 
591 
592  typedef typename FEMB::FEM FEM;
593  typedef typename CONB::CON CON;
594 
595  typedef VBET VBE;
597 
598  typedef N NT;
601  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
603 
604  // constructor making the grid function space an all that is needed
605  CGSpace (Grid& grid, const BCType& bctype)
606  : gv(grid.leafGridView()), femb(gv), conb(grid,bctype)
607  {
608  gfsp = std::shared_ptr<GFS>(new GFS(gv,femb.getFEM(),conb.getCON()));
609  gfsp->name("cgspace");
610  // initialize ordering
611  gfsp->update();
612  conb.postGFSHook(*gfsp);
613  ccp = std::shared_ptr<CC>(new CC());
614  }
615 
616  FEM& getFEM()
617  {
618  return femb.getFEM();
619  }
620 
621  const FEM& getFEM() const
622  {
623  return femb.getFEM();
624  }
625 
626  // return gfs reference
627  GFS& getGFS ()
628  {
629  return *gfsp;
630  }
631 
632  // return gfs reference const version
633  const GFS& getGFS () const
634  {
635  return *gfsp;
636  }
637 
638  // return gfs reference
639  CC& getCC ()
640  {
641  return *ccp;
642  }
643 
644  // return gfs reference const version
645  const CC& getCC () const
646  {
647  return *ccp;
648  }
649 
650  void assembleConstraints (const BCType& bctype)
651  {
652  ccp->clear();
653  constraints(bctype,*gfsp,*ccp);
654  }
655 
657  {
658  ccp->clear();
659  }
660 
661  void setConstrainedDOFS (DOF& x, NT nt) const
662  {
663  set_constrained_dofs(*ccp,nt,x);
664  conb.make_consistent(*gfsp,x);
665  }
666 
667  void setNonConstrainedDOFS (DOF& x, NT nt) const
668  {
669  set_nonconstrained_dofs(*ccp,nt,x);
670  conb.make_consistent(*gfsp,x);
671  }
672 
673  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
674  {
675  copy_constrained_dofs(*ccp,xin,xout);
676  conb.make_consistent(*gfsp,xout);
677  }
678 
679  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
680  {
681  copy_nonconstrained_dofs(*ccp,xin,xout);
682  conb.make_consistent(*gfsp,xout);
683  }
684 
685  private:
686  GV gv; // need this object here because FEM and GFS store a const reference !!
687  FEMB femb;
688  CONB conb;
689  std::shared_ptr<GFS> gfsp;
690  std::shared_ptr<CC> ccp;
691  };
692 
693  // template specialization for nonoverlapping case
694  template<typename T, typename N, unsigned int degree, typename BCType,
695  Dune::GeometryType::BasicType gt, MeshType mt,
696  typename VBET>
697  class CGSpace<T, N, degree, BCType, gt, mt, SolverCategory::nonoverlapping, VBET> {
698  public:
699 
700  // export types
701  typedef T Grid;
702  typedef typename T::LeafGridView GV;
703  typedef typename T::ctype ctype;
705  static const int dim = T::dimension;
706  static const int dimworld = T::dimensionworld;
707 
710 
711  typedef typename FEMB::FEM FEM;
712  typedef typename CONB::CON CON;
713 
714  typedef VBET VBE;
716 
717  typedef N NT;
720  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
722 
723  // constructor making the grid function space an all that is needed
724  CGSpace (Grid& grid, const BCType& bctype)
725  : gv(grid.leafGridView()), es(gv), femb(es), conb(grid,bctype)
726  {
727  gfsp = std::shared_ptr<GFS>(new GFS(es,femb.getFEM(),conb.getCON()));
728  gfsp->name("cgspace");
729  // initialize ordering
730  gfsp->update();
731  // conb.postGFSHook(*gfsp);
732  ccp = std::shared_ptr<CC>(new CC());
733  }
734 
735  FEM& getFEM()
736  {
737  return femb.getFEM();
738  }
739 
740  const FEM& getFEM() const
741  {
742  return femb.getFEM();
743  }
744 
745  // return gfs reference
746  GFS& getGFS ()
747  {
748  return *gfsp;
749  }
750 
751  // return gfs reference const version
752  const GFS& getGFS () const
753  {
754  return *gfsp;
755  }
756 
757  // return gfs reference
758  CC& getCC ()
759  {
760  return *ccp;
761  }
762 
763  // return gfs reference const version
764  const CC& getCC () const
765  {
766  return *ccp;
767  }
768 
769  void assembleConstraints (const BCType& bctype)
770  {
771  ccp->clear();
772  constraints(bctype,*gfsp,*ccp);
773  }
774 
776  {
777  ccp->clear();
778  }
779 
780  void setConstrainedDOFS (DOF& x, NT nt) const
781  {
782  set_constrained_dofs(*ccp,nt,x);
783  conb.make_consistent(*gfsp,x);
784  }
785 
786  void setNonConstrainedDOFS (DOF& x, NT nt) const
787  {
788  set_nonconstrained_dofs(*ccp,nt,x);
789  conb.make_consistent(*gfsp,x);
790  }
791 
792  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
793  {
794  copy_constrained_dofs(*ccp,xin,xout);
795  conb.make_consistent(*gfsp,xout);
796  }
797 
798  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
799  {
800  copy_nonconstrained_dofs(*ccp,xin,xout);
801  conb.make_consistent(*gfsp,xout);
802  }
803 
804  private:
805  GV gv; // need this object here because FEM and GFS store a const reference !!
806  ES es;
807  FEMB femb;
808  CONB conb;
809  std::shared_ptr<GFS> gfsp;
810  std::shared_ptr<CC> ccp;
811  };
812 
813 
814 
815  //============================================================================
816  // Discontinuous Finite Element Space
817  //============================================================================
818 
819  // constraints base template
820  template<SolverCategory::Category st>
821  class DGCONBase
822  {};
823 
824  template<>
825  class DGCONBase<SolverCategory::sequential>
826  {
827  public:
830  {
831  conp = std::shared_ptr<CON>(new CON());
832  }
833  CON& getCON() {return *conp;}
834  const CON& getCON() const {return *conp;}
835  template<typename GFS, typename DOF>
836  void make_consistent (const GFS& gfs, DOF& x) const {}
837  private:
838  std::shared_ptr<CON> conp;
839  };
840 
841  template<>
842  class DGCONBase<SolverCategory::nonoverlapping>
843  {
844  public:
847  {
848  conp = std::shared_ptr<CON>(new CON());
849  }
850  CON& getCON() {return *conp;}
851  const CON& getCON() const {return *conp;}
852  template<typename GFS, typename DOF>
853  void make_consistent (const GFS& gfs, DOF& x) const {}
854  private:
855  std::shared_ptr<CON> conp;
856  };
857 
858  template<>
859  class DGCONBase<SolverCategory::overlapping>
860  {
861  public:
864  {
865  conp = std::shared_ptr<CON>(new CON());
866  }
867  CON& getCON() {return *conp;}
868  const CON& getCON() const {return *conp;}
869  template<typename GFS, typename DOF>
870  void make_consistent (const GFS& gfs, DOF& x) const
871  {
872  // make vector consistent; this is needed for all overlapping solvers
873  ISTL::ParallelHelper<GFS> helper(gfs);
874  helper.maskForeignDOFs(Backend::native(x));
876  if (gfs.gridView().comm().size()>1)
877  gfs.gridView().communicate(adddh,Dune::InteriorBorder_All_Interface,Dune::ForwardCommunication);
878  }
879  private:
880  std::shared_ptr<CON> conp;
881  };
882 
883  // Discontinuous space
884  // default implementation, use only specializations below
885  template<typename T, typename N, unsigned int degree,
886  Dune::GeometryType::BasicType gt, SolverCategory::Category st = SolverCategory::sequential,
888  class DGPkSpace
889  {
890  public:
891 
892  // export types
893  typedef T Grid;
894  typedef typename T::LeafGridView GV;
895  typedef typename T::ctype ctype;
896  static const int dim = T::dimension;
897  static const int dimworld = T::dimensionworld;
898  typedef N NT;
899 #if HAVE_GMP
901 #else
903 #endif
905  typedef typename CONB::CON CON;
906  typedef VBET VBE;
910  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
912 
913  // constructor making the grid function space an all that is needed
914  DGPkSpace (const GV& gridview) : gv(gridview), conb()
915  {
916  femp = std::shared_ptr<FEM>(new FEM());
917  gfsp = std::shared_ptr<GFS>(new GFS(gv,*femp));
918  // initialize ordering
919  gfsp->update();
920  ccp = std::shared_ptr<CC>(new CC());
921  }
922 
923  FEM& getFEM() { return *femp; }
924  const FEM& getFEM() const { return *femp; }
925 
926  // return gfs reference
927  GFS& getGFS () { return *gfsp; }
928 
929  // return gfs reference const version
930  const GFS& getGFS () const {return *gfsp;}
931 
932  // return gfs reference
933  CC& getCC () { return *ccp;}
934 
935  // return gfs reference const version
936  const CC& getCC () const { return *ccp;}
937 
938  template<class BCTYPE>
939  void assembleConstraints (const BCTYPE& bctype)
940  {
941  ccp->clear();
942  constraints(bctype,*gfsp,*ccp);
943  }
944 
946  {
947  ccp->clear();
948  }
949 
950  void setConstrainedDOFS (DOF& x, NT nt) const
951  {
952  set_constrained_dofs(*ccp,nt,x);
953  conb.make_consistent(*gfsp,x);
954  }
955 
956  void setNonConstrainedDOFS (DOF& x, NT nt) const
957  {
958  set_nonconstrained_dofs(*ccp,nt,x);
959  conb.make_consistent(*gfsp,x);
960  }
961 
962  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
963  {
964  copy_constrained_dofs(*ccp,xin,xout);
965  conb.make_consistent(*gfsp,xout);
966  }
967 
968  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
969  {
970  copy_nonconstrained_dofs(*ccp,xin,xout);
971  conb.make_consistent(*gfsp,xout);
972  }
973 
974  private:
975  GV gv; // need this object here because FEM and GFS store a const reference !!
976  CONB conb;
977  std::shared_ptr<FEM> femp;
978  std::shared_ptr<GFS> gfsp;
979  std::shared_ptr<CC> ccp;
980  };
981 
982  // Discontinuous space
983  // default implementation, use only specializations below
984  template<typename T, typename N, unsigned int degree,
985  Dune::GeometryType::BasicType gt, SolverCategory::Category st = SolverCategory::sequential,
986  //typename VBET=ISTL::VectorBackend<ISTL::Blocking::fixed,Dune::PB::PkSize<degree,T::dimension>::value> >
987  typename VBET=ISTL::VectorBackend<> >
989  {
990  public:
991 
992  // export types
993  typedef T Grid;
994  typedef typename T::LeafGridView GV;
995  typedef typename T::ctype ctype;
996  static const int dim = T::dimension;
997  static const int dimworld = T::dimensionworld;
998  typedef N NT;
999 #if HAVE_GMP
1001 #else
1003 #endif
1005  typedef typename CONB::CON CON;
1006  typedef VBET VBE;
1010  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
1012 
1013  // constructor making the grid function space an all that is needed
1014  DGQkOPBSpace (const GV& gridview) : gv(gridview), conb()
1015  {
1016  femp = std::shared_ptr<FEM>(new FEM());
1017  gfsp = std::shared_ptr<GFS>(new GFS(gv,*femp));
1018  // initialize ordering
1019  gfsp->update();
1020  ccp = std::shared_ptr<CC>(new CC());
1021  }
1022 
1023  FEM& getFEM() { return *femp; }
1024  const FEM& getFEM() const { return *femp; }
1025 
1026  // return gfs reference
1027  GFS& getGFS () { return *gfsp; }
1028 
1029  // return gfs reference const version
1030  const GFS& getGFS () const {return *gfsp;}
1031 
1032  // return gfs reference
1033  CC& getCC () { return *ccp;}
1034 
1035  // return gfs reference const version
1036  const CC& getCC () const { return *ccp;}
1037 
1038  template<class BCTYPE>
1039  void assembleConstraints (const BCTYPE& bctype)
1040  {
1041  ccp->clear();
1042  constraints(bctype,*gfsp,*ccp);
1043  }
1044 
1046  {
1047  ccp->clear();
1048  }
1049 
1050  void setConstrainedDOFS (DOF& x, NT nt) const
1051  {
1052  set_constrained_dofs(*ccp,nt,x);
1053  conb.make_consistent(*gfsp,x);
1054  }
1055 
1056  void setNonConstrainedDOFS (DOF& x, NT nt) const
1057  {
1058  set_nonconstrained_dofs(*ccp,nt,x);
1059  conb.make_consistent(*gfsp,x);
1060  }
1061 
1062  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
1063  {
1064  copy_constrained_dofs(*ccp,xin,xout);
1065  conb.make_consistent(*gfsp,xout);
1066  }
1067 
1068  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
1069  {
1070  copy_nonconstrained_dofs(*ccp,xin,xout);
1071  conb.make_consistent(*gfsp,xout);
1072  }
1073 
1074  private:
1075  GV gv; // need this object here because FEM and GFS store a const reference !!
1076  CONB conb;
1077  std::shared_ptr<FEM> femp;
1078  std::shared_ptr<GFS> gfsp;
1079  std::shared_ptr<CC> ccp;
1080  };
1081 
1082  // Discontinuous space
1083  // default implementation, use only specializations below
1084  template<typename T, typename N, unsigned int degree,
1085  Dune::GeometryType::BasicType gt, SolverCategory::Category st = SolverCategory::sequential,
1088  {
1089  public:
1090 
1091  // export types
1092  typedef T Grid;
1093  typedef typename T::LeafGridView GV;
1094  typedef typename T::ctype ctype;
1095  static const int dim = T::dimension;
1096  static const int dimworld = T::dimensionworld;
1097  typedef N NT;
1098  typedef QkDGLocalFiniteElementMap<ctype,NT,degree,dim> FEM;
1100  typedef typename CONB::CON CON;
1101  typedef VBET VBE;
1105  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
1107 
1108  // constructor making the grid function space an all that is needed
1109  DGQkSpace (const GV& gridview) : gv(gridview), conb()
1110  {
1111  femp = std::shared_ptr<FEM>(new FEM());
1112  gfsp = std::shared_ptr<GFS>(new GFS(gv,*femp));
1113  // initialize ordering
1114  gfsp->update();
1115  ccp = std::shared_ptr<CC>(new CC());
1116  }
1117 
1118  FEM& getFEM() { return *femp; }
1119  const FEM& getFEM() const { return *femp; }
1120 
1121  // return gfs reference
1122  GFS& getGFS () { return *gfsp; }
1123 
1124  // return gfs reference const version
1125  const GFS& getGFS () const {return *gfsp;}
1126 
1127  // return gfs reference
1128  CC& getCC () { return *ccp;}
1129 
1130  // return gfs reference const version
1131  const CC& getCC () const { return *ccp;}
1132 
1133  template<class BCTYPE>
1134  void assembleConstraints (const BCTYPE& bctype)
1135  {
1136  ccp->clear();
1137  constraints(bctype,*gfsp,*ccp);
1138  }
1139 
1141  {
1142  ccp->clear();
1143  }
1144 
1145  void setConstrainedDOFS (DOF& x, NT nt) const
1146  {
1147  set_constrained_dofs(*ccp,nt,x);
1148  conb.make_consistent(*gfsp,x);
1149  }
1150 
1151  void setNonConstrainedDOFS (DOF& x, NT nt) const
1152  {
1153  set_nonconstrained_dofs(*ccp,nt,x);
1154  conb.make_consistent(*gfsp,x);
1155  }
1156 
1157  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
1158  {
1159  copy_constrained_dofs(*ccp,xin,xout);
1160  conb.make_consistent(*gfsp,xout);
1161  }
1162 
1163  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
1164  {
1165  copy_nonconstrained_dofs(*ccp,xin,xout);
1166  conb.make_consistent(*gfsp,xout);
1167  }
1168 
1169  private:
1170  GV gv; // need this object here because FEM and GFS store a const reference !!
1171  CONB conb;
1172  std::shared_ptr<FEM> femp;
1173  std::shared_ptr<GFS> gfsp;
1174  std::shared_ptr<CC> ccp;
1175  };
1176 
1177 
1178  // Discontinuous space using QK with Gauss Lobatto points (use only for cube elements)
1179  template<typename T, typename N, unsigned int degree,
1180  Dune::GeometryType::BasicType gt, SolverCategory::Category st = SolverCategory::sequential,
1181  //typename VBET=ISTL::VectorBackend<ISTL::Blocking::fixed,Dune::QkStuff::QkSize<degree,T::dimension>::value> >
1182  typename VBET=ISTL::VectorBackend<> >
1184  {
1185  public:
1186 
1187  // export types
1188  typedef T Grid;
1189  typedef typename T::LeafGridView GV;
1190  typedef typename T::ctype ctype;
1191  static const int dim = T::dimension;
1192  static const int dimworld = T::dimensionworld;
1193  typedef N NT;
1194  typedef QkDGLocalFiniteElementMap<ctype,NT,degree,dim,QkDGBasisPolynomial::lobatto> FEM;
1196  typedef typename CONB::CON CON;
1197  typedef VBET VBE;
1201  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
1203 
1204  // constructor making the grid function space an all that is needed
1205  DGQkGLSpace (const GV& gridview) : gv(gridview), conb()
1206  {
1207  femp = std::shared_ptr<FEM>(new FEM());
1208  gfsp = std::shared_ptr<GFS>(new GFS(gv,*femp));
1209  // initialize ordering
1210  gfsp->update();
1211  ccp = std::shared_ptr<CC>(new CC());
1212  }
1213 
1214  FEM& getFEM() { return *femp; }
1215  const FEM& getFEM() const { return *femp; }
1216 
1217  // return gfs reference
1218  GFS& getGFS () { return *gfsp; }
1219 
1220  // return gfs reference const version
1221  const GFS& getGFS () const {return *gfsp;}
1222 
1223  // return gfs reference
1224  CC& getCC () { return *ccp;}
1225 
1226  // return gfs reference const version
1227  const CC& getCC () const { return *ccp;}
1228 
1229  template<class BCTYPE>
1230  void assembleConstraints (const BCTYPE& bctype)
1231  {
1232  ccp->clear();
1233  constraints(bctype,*gfsp,*ccp);
1234  }
1235 
1237  {
1238  ccp->clear();
1239  }
1240 
1241  void setConstrainedDOFS (DOF& x, NT nt) const
1242  {
1243  set_constrained_dofs(*ccp,nt,x);
1244  conb.make_consistent(*gfsp,x);
1245  }
1246 
1247  void setNonConstrainedDOFS (DOF& x, NT nt) const
1248  {
1249  set_nonconstrained_dofs(*ccp,nt,x);
1250  conb.make_consistent(*gfsp,x);
1251  }
1252 
1253  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
1254  {
1255  copy_constrained_dofs(*ccp,xin,xout);
1256  conb.make_consistent(*gfsp,xout);
1257  }
1258 
1259  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
1260  {
1261  copy_nonconstrained_dofs(*ccp,xin,xout);
1262  conb.make_consistent(*gfsp,xout);
1263  }
1264 
1265  private:
1266  GV gv; // need this object here because FEM and GFS store a const reference !!
1267  CONB conb;
1268  std::shared_ptr<FEM> femp;
1269  std::shared_ptr<GFS> gfsp;
1270  std::shared_ptr<CC> ccp;
1271  };
1272 
1273 
1274  // Discontinuous space using Legendre polynomials (use only for cube elements)
1275  template<typename T, typename N, unsigned int degree,
1276  Dune::GeometryType::BasicType gt, SolverCategory::Category st = SolverCategory::sequential,
1277  //typename VBET=ISTL::VectorBackend<ISTL::Blocking::fixed,Dune::QkStuff::QkSize<degree,T::dimension>::value> >
1278  typename VBET=ISTL::VectorBackend<> >
1280  {
1281  public:
1282 
1283  // export types
1284  typedef T Grid;
1285  typedef typename T::LeafGridView GV;
1286  typedef typename T::ctype ctype;
1287  static const int dim = T::dimension;
1288  static const int dimworld = T::dimensionworld;
1289  typedef N NT;
1290  typedef QkDGLocalFiniteElementMap<ctype,NT,degree,dim,QkDGBasisPolynomial::legendre> FEM;
1292  typedef typename CONB::CON CON;
1293  typedef VBET VBE;
1297  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
1299 
1300  // constructor making the grid function space an all that is needed
1301  DGLegendreSpace (const GV& gridview) : gv(gridview), conb()
1302  {
1303  femp = std::shared_ptr<FEM>(new FEM());
1304  gfsp = std::shared_ptr<GFS>(new GFS(gv,*femp));
1305  // initialize ordering
1306  gfsp->update();
1307  ccp = std::shared_ptr<CC>(new CC());
1308  }
1309 
1310  FEM& getFEM() { return *femp; }
1311  const FEM& getFEM() const { return *femp; }
1312 
1313  // return gfs reference
1314  GFS& getGFS () { return *gfsp; }
1315 
1316  // return gfs reference const version
1317  const GFS& getGFS () const {return *gfsp;}
1318 
1319  // return gfs reference
1320  CC& getCC () { return *ccp;}
1321 
1322  // return gfs reference const version
1323  const CC& getCC () const { return *ccp;}
1324 
1325  template<class BCTYPE>
1326  void assembleConstraints (const BCTYPE& bctype)
1327  {
1328  ccp->clear();
1329  constraints(bctype,*gfsp,*ccp);
1330  }
1331 
1333  {
1334  ccp->clear();
1335  }
1336 
1337  void setConstrainedDOFS (DOF& x, NT nt) const
1338  {
1339  set_constrained_dofs(*ccp,nt,x);
1340  conb.make_consistent(*gfsp,x);
1341  }
1342 
1343  void setNonConstrainedDOFS (DOF& x, NT nt) const
1344  {
1345  set_nonconstrained_dofs(*ccp,nt,x);
1346  conb.make_consistent(*gfsp,x);
1347  }
1348 
1349  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
1350  {
1351  copy_constrained_dofs(*ccp,xin,xout);
1352  conb.make_consistent(*gfsp,xout);
1353  }
1354 
1355  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
1356  {
1357  copy_nonconstrained_dofs(*ccp,xin,xout);
1358  conb.make_consistent(*gfsp,xout);
1359  }
1360 
1361  private:
1362  GV gv; // need this object here because FEM and GFS store a const reference !!
1363  CONB conb;
1364  std::shared_ptr<FEM> femp;
1365  std::shared_ptr<GFS> gfsp;
1366  std::shared_ptr<CC> ccp;
1367  };
1368 
1369 
1370  // Discontinuous P0 space
1371  template<typename T, typename N,
1372  Dune::GeometryType::BasicType gt, SolverCategory::Category st = SolverCategory::sequential,
1373  typename VBET=ISTL::VectorBackend<> >
1374  class P0Space
1375  {
1376  public:
1377 
1378  // export types
1379  typedef T Grid;
1380  typedef typename T::LeafGridView GV;
1381  typedef typename T::ctype ctype;
1382  static const int dim = T::dimension;
1383  static const int dimworld = T::dimensionworld;
1384  typedef N NT;
1387  typedef typename CONB::CON CON;
1388  typedef VBET VBE;
1392  typedef typename GFS::template ConstraintsContainer<N>::Type CC;
1394 
1395  // constructor making the grid function space an all that is needed
1396  P0Space (const GV& gridview) : gv(gridview), conb()
1397  {
1398  femp = std::shared_ptr<FEM>(new FEM(Dune::GeometryType(gt,dim)));
1399  gfsp = std::shared_ptr<GFS>(new GFS(gv,*femp));
1400  // initialize ordering
1401  gfsp->update();
1402  ccp = std::shared_ptr<CC>(new CC());
1403  }
1404 
1405  FEM& getFEM() { return *femp; }
1406  const FEM& getFEM() const { return *femp; }
1407 
1408  // return gfs reference
1409  GFS& getGFS () { return *gfsp; }
1410 
1411  // return gfs reference const version
1412  const GFS& getGFS () const {return *gfsp;}
1413 
1414  // return gfs reference
1415  CC& getCC () { return *ccp;}
1416 
1417  // return gfs reference const version
1418  const CC& getCC () const { return *ccp;}
1419 
1420  template<class BCTYPE>
1421  void assembleConstraints (const BCTYPE& bctype)
1422  {
1423  ccp->clear();
1424  constraints(bctype,*gfsp,*ccp);
1425  }
1426 
1428  {
1429  ccp->clear();
1430  }
1431 
1432  void setConstrainedDOFS (DOF& x, NT nt) const
1433  {
1434  set_constrained_dofs(*ccp,nt,x);
1435  conb.make_consistent(*gfsp,x);
1436  }
1437 
1438  void setNonConstrainedDOFS (DOF& x, NT nt) const
1439  {
1440  set_nonconstrained_dofs(*ccp,nt,x);
1441  conb.make_consistent(*gfsp,x);
1442  }
1443 
1444  void copyConstrainedDOFS (const DOF& xin, DOF& xout) const
1445  {
1446  copy_constrained_dofs(*ccp,xin,xout);
1447  conb.make_consistent(*gfsp,xout);
1448  }
1449 
1450  void copyNonConstrainedDOFS (const DOF& xin, DOF& xout) const
1451  {
1452  copy_nonconstrained_dofs(*ccp,xin,xout);
1453  conb.make_consistent(*gfsp,xout);
1454  }
1455 
1456  private:
1457  GV gv; // need this object here because FEM and GFS store a const reference !!
1458  CONB conb;
1459  std::shared_ptr<FEM> femp;
1460  std::shared_ptr<GFS> gfsp;
1461  std::shared_ptr<CC> ccp;
1462  };
1463 
1464 
1465  // how can we most easily specify a grid function
1466  // pass a function space as parameter
1467  template<typename FS, typename Functor>
1469  : public GridFunctionBase<GridFunctionTraits<typename FS::GV, typename FS::NT,
1470  1,FieldVector<typename FS::NT,1> >
1471  ,UserFunction<FS,Functor> >
1472  {
1473  public:
1474  typedef GridFunctionTraits<typename FS::GV, typename FS::NT,
1475  1,FieldVector<typename FS::NT,1> > Traits;
1476 
1478  UserFunction (const FS& fs_, const Functor& f_)
1479  : fs(fs_), f(f_)
1480  {}
1481 
1483  inline void evaluate (const typename Traits::ElementType& e,
1484  const typename Traits::DomainType& x,
1485  typename Traits::RangeType& y) const
1486  {
1487  typename Traits::DomainType x_ = e.geometry().global(x);
1488  std::vector<double> x__(x.size());
1489  for (size_t i=0; i<x.size(); ++i) x__[i]=x_[i];
1490  y = f(x__);
1491  }
1492 
1493  inline const typename FS::GV& getGridView () const
1494  {
1495  return fs.getGFS().gridView();
1496  }
1497 
1498  private:
1499  const FS fs; // store a copy of the function space
1500  const Functor f;
1501  };
1502 
1503 
1504  template<typename FS, typename LOP, SolverCategory::Category st = SolverCategory::sequential>
1506  {
1507  public:
1508  // export types
1510  typedef Dune::PDELab::GridOperator<typename FS::GFS,typename FS::GFS,LOP,MBE,
1511  typename FS::NT,typename FS::NT,typename FS::NT,
1512  typename FS::CC,typename FS::CC> GO;
1513  typedef typename GO::Jacobian MAT;
1514 
1515  GalerkinGlobalAssembler (const FS& fs, LOP& lop, const std::size_t nonzeros)
1516  {
1517  gop = std::shared_ptr<GO>(new GO(fs.getGFS(),fs.getCC(),fs.getGFS(),fs.getCC(),lop,MBE(nonzeros)));
1518  }
1519 
1520  // return grid reference
1521  GO& getGO ()
1522  {
1523  return *gop;
1524  }
1525 
1526  // return grid reference const version
1527  const GO& getGO () const
1528  {
1529  return *gop;
1530  }
1531 
1533  {
1534  return *gop;
1535  }
1536 
1538  {
1539  return gop.operator->();
1540  }
1541 
1542  const GO& operator*() const
1543  {
1544  return *gop;
1545  }
1546 
1547  const GO* operator->() const
1548  {
1549  return gop.operator->();
1550  }
1551 
1552  private:
1553  std::shared_ptr<GO> gop;
1554  };
1555 
1556 
1557  template<typename FS, typename LOP, SolverCategory::Category st = SolverCategory::sequential>
1559  {
1560  public:
1561  // export types
1563  typedef Dune::PDELab::GridOperator<typename FS::GFS,typename FS::GFS,LOP,MBE,
1564  typename FS::NT,typename FS::NT,typename FS::NT,
1565  typename FS::CC,typename FS::CC> GO;
1566  typedef typename GO::Jacobian MAT;
1567 
1568  GalerkinGlobalAssemblerNewBackend (const FS& fs, LOP& lop, const MBE& mbe)
1569  {
1570  gop = std::shared_ptr<GO>(new GO(fs.getGFS(),fs.getCC(),fs.getGFS(),fs.getCC(),lop,mbe));
1571  }
1572 
1573  // return grid reference
1574  GO& getGO ()
1575  {
1576  return *gop;
1577  }
1578 
1579  // return grid reference const version
1580  const GO& getGO () const
1581  {
1582  return *gop;
1583  }
1584 
1586  {
1587  return *gop;
1588  }
1589 
1591  {
1592  return gop.operator->();
1593  }
1594 
1595  const GO& operator*() const
1596  {
1597  return *gop;
1598  }
1599 
1600  const GO* operator->() const
1601  {
1602  return gop.operator->();
1603  }
1604 
1605  private:
1606  std::shared_ptr<GO> gop;
1607  };
1608 
1609 
1610  // variant with two different function spaces
1611  template<typename FSU, typename FSV, typename LOP, SolverCategory::Category st>
1613  {
1614  public:
1615  // export types
1617  typedef Dune::PDELab::GridOperator<typename FSU::GFS,typename FSV::GFS,LOP,MBE,
1618  typename FSU::NT,typename FSU::NT,typename FSU::NT,
1619  typename FSU::CC,typename FSV::CC> GO;
1620  typedef typename GO::Jacobian MAT;
1621 
1622  GlobalAssembler (const FSU& fsu, const FSV& fsv, LOP& lop, const std::size_t nonzeros)
1623  {
1624  gop = std::shared_ptr<GO>(new GO(fsu.getGFS(),fsu.getCC(),fsv.getGFS(),fsv.getCC(),lop,MBE(nonzeros)));
1625  }
1626 
1627  // return grid reference
1628  GO& getGO ()
1629  {
1630  return *gop;
1631  }
1632 
1633  // return grid reference const version
1634  const GO& getGO () const
1635  {
1636  return *gop;
1637  }
1638 
1640  {
1641  return *gop;
1642  }
1643 
1645  {
1646  return gop.operator->();
1647  }
1648 
1649  const GO& operator*() const
1650  {
1651  return *gop;
1652  }
1653 
1654  const GO* operator->() const
1655  {
1656  return gop.operator->();
1657  }
1658 
1659  private:
1660  std::shared_ptr<GO> gop;
1661  };
1662 
1663 
1664  template<typename GO1, typename GO2, bool implicit = true>
1666  {
1667  public:
1668  // export types
1671  typedef typename GO::Jacobian MAT;
1672 
1673  OneStepGlobalAssembler (GO1& go1, GO2& go2)
1674  {
1675  gop = std::shared_ptr<GO>(new GO(*go1,*go2));
1676  }
1677 
1678  // return grid reference
1679  GO& getGO ()
1680  {
1681  return *gop;
1682  }
1683 
1684  // return grid reference const version
1685  const GO& getGO () const
1686  {
1687  return *gop;
1688  }
1689 
1691  {
1692  return *gop;
1693  }
1694 
1696  {
1697  return gop.operator->();
1698  }
1699 
1700  const GO& operator*() const
1701  {
1702  return *gop;
1703  }
1704 
1705  const GO* operator->() const
1706  {
1707  return gop.operator->();
1708  }
1709 
1710  private:
1711  std::shared_ptr<GO> gop;
1712  };
1713 
1714 
1715  // packaging of the CG_AMG_SSOR solver: default version is sequential
1716  template<typename FS, typename ASS, SolverCategory::Category st = SolverCategory::sequential>
1718  {
1719  public:
1720  // types exported
1722 
1723  ISTLSolverBackend_CG_AMG_SSOR (const FS& fs, const ASS& ass, unsigned maxiter_=5000,
1724  int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
1725  {
1726  lsp = std::shared_ptr<LS>(new LS(maxiter_,verbose_,reuse_,usesuperlu_));
1727  }
1728 
1729  LS& getLS () {return *lsp;}
1730  const LS& getLS () const { return *lsp;}
1731  LS& operator*(){return *lsp;}
1732  LS* operator->() { return lsp.operator->(); }
1733  const LS& operator*() const{return *lsp;}
1734  const LS* operator->() const{ return lsp.operator->();}
1735 
1736  private:
1737  std::shared_ptr<LS> lsp;
1738  };
1739 
1740  // packaging of the CG_AMG_SSOR solver: nonoverlapping version
1741  template<typename FS, typename ASS>
1742  class ISTLSolverBackend_CG_AMG_SSOR<FS,ASS, SolverCategory::nonoverlapping>
1743  {
1744  public:
1745  // types exported
1747 
1748  ISTLSolverBackend_CG_AMG_SSOR (const FS& fs, const ASS& ass, unsigned maxiter_=5000,
1749  int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
1750  {
1751  lsp = std::shared_ptr<LS>(new LS(fs.getGFS(),maxiter_,verbose_,reuse_,usesuperlu_));
1752  }
1753 
1754  LS& getLS () {return *lsp;}
1755  const LS& getLS () const { return *lsp;}
1756  LS& operator*(){return *lsp;}
1757  LS* operator->() { return lsp.operator->(); }
1758  const LS& operator*() const{return *lsp;}
1759  const LS* operator->() const{ return lsp.operator->();}
1760 
1761  private:
1762  std::shared_ptr<LS> lsp;
1763  };
1764 
1765  // packaging of the CG_AMG_SSOR solver: overlapping version
1766  template<typename FS, typename ASS>
1767  class ISTLSolverBackend_CG_AMG_SSOR<FS,ASS, SolverCategory::overlapping>
1768  {
1769  public:
1770  // types exported
1772 
1773  ISTLSolverBackend_CG_AMG_SSOR (const FS& fs, const ASS& ass, unsigned maxiter_=5000,
1774  int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
1775  {
1776  lsp = std::shared_ptr<LS>(new LS(fs.getGFS(),maxiter_,verbose_,reuse_,usesuperlu_));
1777  }
1778 
1779  LS& getLS () {return *lsp;}
1780  const LS& getLS () const { return *lsp;}
1781  LS& operator*(){return *lsp;}
1782  LS* operator->() { return lsp.operator->(); }
1783  const LS& operator*() const{return *lsp;}
1784  const LS* operator->() const{ return lsp.operator->();}
1785 
1786  private:
1787  std::shared_ptr<LS> lsp;
1788  };
1789 
1790  // packaging of the CG_SSOR solver: default version is sequential
1791  template<typename FS, typename ASS, SolverCategory::Category st = SolverCategory::sequential>
1793  {
1794  public:
1795  // types exported
1797 
1798  ISTLSolverBackend_CG_SSOR (const FS& fs, const ASS& ass, unsigned maxiter_=5000,
1799  int steps_=5, int verbose_=1)
1800  {
1801  lsp = std::shared_ptr<LS>(new LS(maxiter_,verbose_));
1802  }
1803 
1804  LS& getLS () {return *lsp;}
1805  const LS& getLS () const { return *lsp;}
1806  LS& operator*(){return *lsp;}
1807  LS* operator->() { return lsp.operator->(); }
1808  const LS& operator*() const{return *lsp;}
1809  const LS* operator->() const{ return lsp.operator->();}
1810 
1811  private:
1812  std::shared_ptr<LS> lsp;
1813  };
1814 
1815  // packaging of the CG_SSOR solver: nonoverlapping version
1816  template<typename FS, typename ASS>
1817  class ISTLSolverBackend_CG_SSOR<FS,ASS,SolverCategory::nonoverlapping>
1818  {
1819  public:
1820  // types exported
1822 
1823  ISTLSolverBackend_CG_SSOR (const FS& fs, const ASS& ass, unsigned maxiter_=5000,
1824  int steps_=5, int verbose_=1)
1825  {
1826  lsp = std::shared_ptr<LS>(new LS(fs.getGFS(),maxiter_,steps_,verbose_));
1827  }
1828 
1829  LS& getLS () {return *lsp;}
1830  const LS& getLS () const { return *lsp;}
1831  LS& operator*(){return *lsp;}
1832  LS* operator->() { return lsp.operator->(); }
1833  const LS& operator*() const{return *lsp;}
1834  const LS* operator->() const{ return lsp.operator->();}
1835 
1836  private:
1837  std::shared_ptr<LS> lsp;
1838  };
1839 
1840  // packaging of the CG_SSOR solver: overlapping version
1841  template<typename FS, typename ASS>
1842  class ISTLSolverBackend_CG_SSOR<FS,ASS,SolverCategory::overlapping>
1843  {
1844  public:
1845  // types exported
1847 
1848  ISTLSolverBackend_CG_SSOR (const FS& fs, const ASS& ass, unsigned maxiter_=5000,
1849  int steps_=5, int verbose_=1)
1850  {
1851  lsp = std::shared_ptr<LS>(new LS(fs.getGFS(),fs.getCC(),maxiter_,steps_,verbose_));
1852  }
1853 
1854  LS& getLS () {return *lsp;}
1855  const LS& getLS () const { return *lsp;}
1856  LS& operator*(){return *lsp;}
1857  LS* operator->() { return lsp.operator->(); }
1858  const LS& operator*() const{return *lsp;}
1859  const LS* operator->() const{ return lsp.operator->();}
1860 
1861  private:
1862  std::shared_ptr<LS> lsp;
1863  };
1864 
1865 
1866  // packaging of a default solver that should always work
1867  // in the sequential case : BCGS SSOR
1868  template<typename FS, typename ASS, SolverCategory::Category st = SolverCategory::sequential>
1870  {
1871  public:
1872  // types exported
1874 
1875  ISTLSolverBackend_IterativeDefault (const FS& fs, const ASS& ass, unsigned maxiter_=5000, int verbose_=1)
1876  {
1877  lsp = std::shared_ptr<LS>(new LS(maxiter_,verbose_));
1878  }
1879 
1880  LS& getLS () {return *lsp;}
1881  const LS& getLS () const { return *lsp;}
1882  LS& operator*(){return *lsp;}
1883  LS* operator->() { return lsp.operator->(); }
1884  const LS& operator*() const{return *lsp;}
1885  const LS* operator->() const{ return lsp.operator->();}
1886 
1887  private:
1888  std::shared_ptr<LS> lsp;
1889  };
1890 
1891  // in the nonoverlapping case : BCGS SSORk
1892  template<typename FS, typename ASS>
1893  class ISTLSolverBackend_IterativeDefault<FS,ASS,SolverCategory::nonoverlapping>
1894  {
1895  public:
1896  // types exported
1898 
1899  ISTLSolverBackend_IterativeDefault (const FS& fs, const ASS& ass, unsigned maxiter_=5000, int verbose_=1)
1900  {
1901  lsp = std::shared_ptr<LS>(new LS(ass.getGO(),maxiter_,3,verbose_));
1902  }
1903 
1904  LS& getLS () {return *lsp;}
1905  const LS& getLS () const { return *lsp;}
1906  LS& operator*(){return *lsp;}
1907  LS* operator->() { return lsp.operator->(); }
1908  const LS& operator*() const{return *lsp;}
1909  const LS* operator->() const{ return lsp.operator->();}
1910 
1911  private:
1912  std::shared_ptr<LS> lsp;
1913  };
1914 
1915  // in the overlapping case : BCGS SSORk
1916  template<typename FS, typename ASS>
1917  class ISTLSolverBackend_IterativeDefault<FS,ASS,SolverCategory::overlapping>
1918  {
1919  public:
1920  // types exported
1922 
1923  ISTLSolverBackend_IterativeDefault (const FS& fs, const ASS& ass, unsigned maxiter_=5000, int verbose_=1)
1924  {
1925  lsp = std::shared_ptr<LS>(new LS(fs.getGFS(),fs.getCC(),maxiter_,3,verbose_));
1926  }
1927 
1928  LS& getLS () {return *lsp;}
1929  const LS& getLS () const { return *lsp;}
1930  LS& operator*(){return *lsp;}
1931  LS* operator->() { return lsp.operator->(); }
1932  const LS& operator*() const{return *lsp;}
1933  const LS* operator->() const{ return lsp.operator->();}
1934 
1935  private:
1936  std::shared_ptr<LS> lsp;
1937  };
1938 
1939  // packaging of a default solver that should always work
1940  // in the sequential case : BCGS SSOR
1941  template<typename FS, typename ASS, SolverCategory::Category st = SolverCategory::sequential>
1943  {
1944  public:
1945  // types exported
1947 
1948  ISTLSolverBackend_ExplicitDiagonal (const FS& fs, const ASS& ass, unsigned maxiter_=5000, int verbose_=1)
1949  {
1950  lsp = std::shared_ptr<LS>(new LS());
1951  }
1952 
1953  LS& getLS () {return *lsp;}
1954  const LS& getLS () const { return *lsp;}
1955  LS& operator*(){return *lsp;}
1956  LS* operator->() { return lsp.operator->(); }
1957  const LS& operator*() const{return *lsp;}
1958  const LS* operator->() const{ return lsp.operator->();}
1959 
1960  private:
1961  std::shared_ptr<LS> lsp;
1962  };
1963 
1964  // packaging of a default solver that should always work
1965  // in the sequential case : BCGS SSOR
1966  template<typename FS, typename ASS>
1967  class ISTLSolverBackend_ExplicitDiagonal<FS,ASS,SolverCategory::overlapping>
1968  {
1969  public:
1970  // types exported
1972 
1973  ISTLSolverBackend_ExplicitDiagonal (const FS& fs, const ASS& ass, unsigned maxiter_=5000, int verbose_=1)
1974  {
1975  lsp = std::shared_ptr<LS>(new LS(fs.getGFS()));
1976  }
1977 
1978  LS& getLS () {return *lsp;}
1979  const LS& getLS () const { return *lsp;}
1980  LS& operator*(){return *lsp;}
1981  LS* operator->() { return lsp.operator->(); }
1982  const LS& operator*() const{return *lsp;}
1983  const LS* operator->() const{ return lsp.operator->();}
1984 
1985  private:
1986  std::shared_ptr<LS> lsp;
1987  };
1988 
1989  // packaging of a default solver that should always work
1990  // in the sequential case : BCGS SSOR
1991  template<typename FS, typename ASS>
1992  class ISTLSolverBackend_ExplicitDiagonal<FS,ASS,SolverCategory::nonoverlapping>
1993  {
1994  public:
1995  // types exported
1997 
1998  ISTLSolverBackend_ExplicitDiagonal (const FS& fs, const ASS& ass, unsigned maxiter_=5000, int verbose_=1)
1999  {
2000  lsp = std::shared_ptr<LS>(new LS(fs.getGFS()));
2001  }
2002 
2003  LS& getLS () {return *lsp;}
2004  const LS& getLS () const { return *lsp;}
2005  LS& operator*(){return *lsp;}
2006  LS* operator->() { return lsp.operator->(); }
2007  const LS& operator*() const{return *lsp;}
2008  const LS* operator->() const{ return lsp.operator->();}
2009 
2010  private:
2011  std::shared_ptr<LS> lsp;
2012  };
2013 
2014 
2015 } // end namespace PDELab
2016  } // end namespace Dune
2017 
2018 #endif // DUNE_PDELAB_BOILERPLATE_PDELAB_HH
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1056
LS & getLS()
Definition: pdelab.hh:1729
void maskForeignDOFs(X &x) const
Mask out all DOFs not owned by the current process with 0.
Definition: parallelhelper.hh:115
HangingNodesDirichletConstraints< Grid, HangingNodesConstraintsAssemblers::CubeGridQ1Assembler, BCType > CON
Definition: pdelab.hh:473
extend conforming constraints class by processor boundary
Definition: conforming.hh:101
const LS & operator*() const
Definition: pdelab.hh:1808
void clearConstraints()
Definition: pdelab.hh:1332
T::LeafGridView GV
Definition: pdelab.hh:1380
periodic boundary intersection (neighbor() == true && boundary() == true)
T & getGrid()
Definition: pdelab.hh:157
const Grid & operator*() const
Definition: pdelab.hh:309
N NT
Definition: pdelab.hh:898
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1163
const FS::GV & getGridView() const
Definition: pdelab.hh:1493
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1259
N NT
Definition: pdelab.hh:1097
GO * operator->()
Definition: pdelab.hh:1644
const FEM & getFEM() const
Definition: pdelab.hh:1024
ISTLSolverBackend_ExplicitDiagonal(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1)
Definition: pdelab.hh:1998
void assembleConstraints(const BCTYPE &bctype)
Definition: pdelab.hh:1230
FEM & getFEM()
Definition: pdelab.hh:616
Dune::PDELab::ISTLBackend_NOVLP_CG_AMG_SSOR< typename ASS::GO > LS
Definition: pdelab.hh:1746
T & getGrid()
Definition: pdelab.hh:349
T::LeafGridView GV
Definition: pdelab.hh:584
FEM & getFEM()
Definition: pdelab.hh:1310
const FEM & getFEM() const
Definition: pdelab.hh:621
ISTLSolverBackend_CG_SSOR(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
Definition: pdelab.hh:1848
DGCONBase< st > CONB
Definition: pdelab.hh:1386
void clearConstraints()
Definition: pdelab.hh:656
T Grid
Definition: pdelab.hh:329
Nonoverlapping parallel CG solver preconditioned with AMG smoothed by SSOR.
Definition: novlpistlsolverbackend.hh:1045
QkDGLocalFiniteElementMap< ctype, NT, degree, dim > FEM
Definition: pdelab.hh:1098
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:1201
Dune::PDELab::ISTLBackend_NOVLP_BCGS_SSORk< typename ASS::GO > LS
Definition: pdelab.hh:1897
std::shared_ptr< T > getSharedPtr()
Definition: pdelab.hh:343
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1062
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
Definition: pdelab.hh:821
Dune::PDELab::OneStepGridOperator< typename GO1::GO, typename GO2::GO, implicit > GO
Definition: pdelab.hh:1670
const FEM & getFEM() const
Definition: pdelab.hh:1311
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1438
GFS & getGFS()
Definition: pdelab.hh:1218
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:908
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1241
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:962
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:1295
void clearConstraints()
Definition: pdelab.hh:1140
LS & getLS()
Definition: pdelab.hh:1953
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:720
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1151
Dune::PDELab::GridOperator< typename FS::GFS, typename FS::GFS, LOP, MBE, typename FS::NT, typename FS::NT, typename FS::NT, typename FS::CC, typename FS::CC > GO
Definition: pdelab.hh:1512
GFS & getGFS()
Definition: pdelab.hh:627
CGCONBase(Grid &grid, const BCType &bctype, const GV &gv)
Definition: pdelab.hh:501
LS & operator*()
Definition: pdelab.hh:1731
void make_consistent(const GFS &gfs, DOF &x) const
Definition: pdelab.hh:853
CONB::CON CON
Definition: pdelab.hh:1196
ISTLSolverBackend_IterativeDefault(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1)
Definition: pdelab.hh:1923
const GFS & getGFS() const
Definition: pdelab.hh:1030
Dune::FieldVector< GV::Grid::ctype, GV::dimension > DomainType
domain type in dim-size coordinates
Definition: function.hh:48
GO::Jacobian MAT
Definition: pdelab.hh:1671
CONB::CON CON
Definition: pdelab.hh:1292
FEM & getFEM()
Definition: pdelab.hh:1118
Definition: gridoperator/onestep.hh:16
T Grid
Definition: pdelab.hh:1188
T::LeafGridView GV
Definition: pdelab.hh:1093
const T & operator*() const
Definition: pdelab.hh:370
Definition: pdelab.hh:988
FEM & getFEM()
Definition: pdelab.hh:1405
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:1007
DGQkGLSpace(const GV &gridview)
Definition: pdelab.hh:1205
Dune::PDELab::ISTLBackend_OVLP_ExplicitDiagonal< typename FS::GFS > LS
Definition: pdelab.hh:1971
const LS & getLS() const
Definition: pdelab.hh:1805
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: ovlpistlsolverbackend.hh:1010
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1444
GO & getGO()
Definition: pdelab.hh:1521
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:596
const GFS & getGFS() const
Definition: pdelab.hh:930
Definition: l2orthonormal.hh:257
const GO & operator*() const
Definition: pdelab.hh:1542
CGCONBase< Grid, degree, gt, mt, SolverCategory::nonoverlapping, BCType > CONB
Definition: pdelab.hh:709
ISTLSolverBackend_IterativeDefault(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1)
Definition: pdelab.hh:1899
ISTLSolverBackend_CG_AMG_SSOR(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: pdelab.hh:1773
GO::Jacobian MAT
Definition: pdelab.hh:1566
const GO & getGO() const
Definition: pdelab.hh:1634
LS & operator*()
Definition: pdelab.hh:1882
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:1294
void clearConstraints()
Definition: pdelab.hh:1236
GO & operator*()
Definition: pdelab.hh:1532
NoConstraints CON
Definition: pdelab.hh:828
ISTLSolverBackend_CG_AMG_SSOR(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: pdelab.hh:1748
LS * operator->()
Definition: pdelab.hh:1807
std::shared_ptr< T > getSharedPtr()
Definition: pdelab.hh:151
void assembleConstraints(const BCTYPE &bctype)
Definition: pdelab.hh:1039
const GO * operator->() const
Definition: pdelab.hh:1600
const T * operator->() const
Definition: pdelab.hh:183
LS & getLS()
Definition: pdelab.hh:1880
DGCONBase< st > CONB
Definition: pdelab.hh:1291
Definition: istl/descriptors.hh:50
VBET VBE
Definition: pdelab.hh:1388
GO * operator->()
Definition: pdelab.hh:1590
const FEM & getFEM() const
Definition: pdelab.hh:924
CONB::CON CON
Definition: pdelab.hh:1005
Definition: pdelab.hh:1468
CC & getCC()
Definition: pdelab.hh:1224
void assembleConstraints(const BCTYPE &bctype)
Definition: pdelab.hh:939
const GFS & getGFS() const
Definition: pdelab.hh:1125
const GFS & getGFS() const
Definition: pdelab.hh:1412
Sequential conjugate gradient solver preconditioned with AMG smoothed by SSOR.
Definition: seqistlsolverbackend.hh:686
Dune::PDELab::NonOverlappingEntitySet< GV > ES
Definition: pdelab.hh:704
traits class holding the function signature, same as in local function
Definition: function.hh:175
Definition: parallelhelper.hh:53
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:968
Parallel P0 constraints for nonoverlapping grids with ghosts.
Definition: p0ghost.hh:16
LS & operator*()
Definition: pdelab.hh:1955
LS * operator->()
Definition: pdelab.hh:1732
Definition: pdelab.hh:1183
Dune::PDELab::ISTLBackend_NOVLP_ExplicitDiagonal< typename FS::GFS > LS
Definition: pdelab.hh:1996
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:1200
T & operator*()
Definition: pdelab.hh:360
T::ctype ctype
Definition: pdelab.hh:895
CC & getCC()
Definition: pdelab.hh:1033
CONB::CON CON
Definition: pdelab.hh:593
VBET VBE
Definition: pdelab.hh:1101
const Grid * operator->() const
Definition: pdelab.hh:314
const FEM & getFEM() const
Definition: pdelab.hh:1406
ISTLSolverBackend_IterativeDefault(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1)
Definition: pdelab.hh:1875
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:909
std::shared_ptr< Grid > getSharedPtr()
Definition: pdelab.hh:282
Grid & operator*()
Definition: pdelab.hh:299
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:599
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1253
OneStepGlobalAssembler(GO1 &go1, GO2 &go2)
Definition: pdelab.hh:1673
const FEM & getFEM() const
Definition: pdelab.hh:1215
Dune::PDELab::GridOperator< typename FS::GFS, typename FS::GFS, LOP, MBE, typename FS::NT, typename FS::NT, typename FS::NT, typename FS::CC, typename FS::CC > GO
Definition: pdelab.hh:1565
T * operator->()
Definition: pdelab.hh:173
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1450
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:1298
P0ParallelGhostConstraints CON
Definition: pdelab.hh:845
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
T Grid
Definition: pdelab.hh:92
HangingNodesDirichletConstraints< Grid, HangingNodesConstraintsAssemblers::SimplexGridP1Assembler, BCType > CON
Definition: pdelab.hh:447
ISTLBackend_SEQ_CG_SSOR LS
Definition: pdelab.hh:1796
StructuredGrid(Dune::GeometryType::BasicType meshtype, array< double, dimworld > lower_left, array< double, dimworld > upper_right, array< unsigned int, dim > cells)
Definition: pdelab.hh:118
GO & operator*()
Definition: pdelab.hh:1585
MeshType
Definition: pdelab.hh:433
DGCONBase< st > CONB
Definition: pdelab.hh:1195
const GO * operator->() const
Definition: pdelab.hh:1654
ISTL::BCRSMatrixBackend MBE
Definition: pdelab.hh:1509
T Grid
Definition: pdelab.hh:583
const CC & getCC() const
Definition: pdelab.hh:645
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:1392
DGPkSpace(const GV &gridview)
Definition: pdelab.hh:914
Definition: genericdatahandle.hh:622
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:956
GO::Jacobian MAT
Definition: pdelab.hh:1513
T::LeafGridView GV
Definition: pdelab.hh:894
VBET VBE
Definition: pdelab.hh:1293
void set_nonconstrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
Definition: constraints.hh:962
CONB::CON CON
Definition: pdelab.hh:1100
const LS & operator*() const
Definition: pdelab.hh:1884
GO & getGO()
Definition: pdelab.hh:1628
T & operator*()
Definition: pdelab.hh:168
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:1102
convert a grid function space and a coefficient vector into a grid function
Definition: gridfunctionspaceutilities.hh:53
DGCONBase< st > CONB
Definition: pdelab.hh:1099
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:1391
Standard grid operator implementation.
Definition: gridoperator.hh:57
const GFS & getGFS() const
Definition: pdelab.hh:1221
LS & operator*()
Definition: pdelab.hh:1806
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:1106
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: novlpistlsolverbackend.hh:633
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1343
const GO & getGO() const
Definition: pdelab.hh:1685
Dune::PDELab::ISTLBackend_CG_AMG_SSOR< typename ASS::GO > LS
Definition: pdelab.hh:1771
P0Space(const GV &gridview)
Definition: pdelab.hh:1396
Overlapping parallel CGS solver with SSOR preconditioner.
Definition: ovlpistlsolverbackend.hh:727
DGQkOPBSpace(const GV &gridview)
Definition: pdelab.hh:1014
Solver to be used for explicit time-steppers with (block-)diagonal mass matrix.
Definition: seqistlsolverbackend.hh:503
VBET VBE
Definition: pdelab.hh:1006
Overlapping parallel conjugate gradient solver preconditioned with AMG smoothed by SSOR...
Definition: ovlpistlsolverbackend.hh:1258
CC & getCC()
Definition: pdelab.hh:1415
GFS & getGFS()
Definition: pdelab.hh:1122
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:1198
CGCONBase< Grid, degree, gt, mt, st, BCType > CONB
Definition: pdelab.hh:590
Definition: pdelab.hh:888
CC & getCC()
Definition: pdelab.hh:1128
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:661
VBET VBE
Definition: pdelab.hh:595
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:602
Overlapping parallel BiCGStab solver with SSOR preconditioner.
Definition: ovlpistlsolverbackend.hh:661
Definition: partitionviewentityset.hh:34
ISTLSolverBackend_CG_SSOR(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
Definition: pdelab.hh:1798
Definition: pdelab.hh:440
wrap a GridFunction so it can be used with the VTKWriter from dune-grid.
Definition: vtkexport.hh:22
const GFS & getGFS() const
Definition: pdelab.hh:1317
ISTLSolverBackend_CG_AMG_SSOR(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1, bool reuse_=false, bool usesuperlu_=true)
Definition: pdelab.hh:1723
StructuredGrid(Dune::GeometryType::BasicType meshtype, unsigned int cells, int overlap=1)
Definition: pdelab.hh:205
OPBLocalFiniteElementMap< ctype, NT, degree, dim, gt, N, Dune::PB::BasisType::Qk > FEM
Definition: pdelab.hh:1002
QkDGLocalFiniteElementMap< ctype, NT, degree, dim, QkDGBasisPolynomial::legendre > FEM
Definition: pdelab.hh:1290
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:1202
T * operator->()
Definition: pdelab.hh:365
GO & operator*()
Definition: pdelab.hh:1639
DGLegendreSpace(const GV &gridview)
Definition: pdelab.hh:1301
T::ctype ctype
Definition: pdelab.hh:1381
LS * operator->()
Definition: pdelab.hh:1883
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:950
const Entity & e
Definition: localfunctionspace.hh:111
Grid & getGrid()
Definition: pdelab.hh:288
ISTL::BCRSMatrixBackend MBE
Definition: pdelab.hh:1669
const GO & operator*() const
Definition: pdelab.hh:1649
N NT
Definition: pdelab.hh:598
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperator.hh:75
StructuredGrid(Dune::GeometryType::BasicType meshtype, array< double, dimworld > lower_left, array< double, dimworld > upper_right, array< unsigned int, dim > cells, int overlap=1)
Definition: pdelab.hh:221
Dune::PDELab::ISTLBackend_SEQ_ExplicitDiagonal LS
Definition: pdelab.hh:1946
GO * operator->()
Definition: pdelab.hh:1695
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:1008
void evaluate(const typename Traits::ElementType &e, const typename Traits::DomainType &x, typename Traits::RangeType &y) const
Evaluate the GridFunction at given position.
Definition: pdelab.hh:1483
GridFunctionSpace< ES, FEM, CON, VBE > GFS
Definition: pdelab.hh:715
const CC & getCC() const
Definition: pdelab.hh:1036
FEM & getFEM()
Definition: pdelab.hh:923
const Grid & getGrid() const
Definition: pdelab.hh:294
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:907
Definition: l2orthonormal.hh:257
T::LeafGridView GV
Definition: pdelab.hh:994
StructuredGrid(Dune::GeometryType::BasicType meshtype, array< double, dimworld > lower_left, array< double, dimworld > upper_right, array< unsigned int, dim > cells, array< bool, dim > periodic, int overlap=1)
Definition: pdelab.hh:251
Nonoverlapping parallel BiCGSTAB solver preconditioned by block SSOR.
Definition: novlpistlsolverbackend.hh:825
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1349
QkLocalFiniteElementMap< GV, C, R, degree > FEM
Definition: pdelab.hh:416
Definition: pdelab.hh:1374
FEMB::FEM FEM
Definition: pdelab.hh:592
Traits::Jacobian Jacobian
Definition: gridoperator/onestep.hh:55
void assembleConstraints(const BCType &bctype)
Definition: pdelab.hh:650
const GO * operator->() const
Definition: pdelab.hh:1705
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1145
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:1393
void clearConstraints()
Definition: pdelab.hh:945
Parallel P0 constraints for overlapping grids.
Definition: p0.hh:15
Definition: pdelab.hh:1087
T::ctype ctype
Definition: pdelab.hh:585
N NT
Definition: pdelab.hh:1289
T::ctype ctype
Definition: pdelab.hh:1190
A grid function space.
Definition: gridfunctionspace.hh:169
void clearConstraints()
Definition: pdelab.hh:1045
GO & operator*()
Definition: pdelab.hh:1690
static const int dim
Definition: pdelab.hh:94
Dirichlet Constraints construction.
Definition: conforming.hh:36
GO::Jacobian MAT
Definition: pdelab.hh:1620
DGQkSpace(const GV &gridview)
Definition: pdelab.hh:1109
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:1199
QkDGLocalFiniteElementMap< ctype, NT, degree, dim, QkDGBasisPolynomial::lobatto > FEM
Definition: pdelab.hh:1194
const LS & getLS() const
Definition: pdelab.hh:1881
const LS & getLS() const
Definition: pdelab.hh:1730
CC & getCC()
Definition: pdelab.hh:639
const CC & getCC() const
Definition: pdelab.hh:1323
Dune::PDELab::GridOperator< typename FSU::GFS, typename FSV::GFS, LOP, MBE, typename FSU::NT, typename FSU::NT, typename FSU::NT, typename FSU::CC, typename FSV::CC > GO
Definition: pdelab.hh:1619
Hanging Node constraints construction.
Definition: hangingnode.hh:318
ISTLBackend_OVLP_CG_SSORk< typename FS::GFS, typename FS::CC > LS
Definition: pdelab.hh:1846
void assembleConstraints(const BCTYPE &bctype)
Definition: pdelab.hh:1421
const LS * operator->() const
Definition: pdelab.hh:1734
T::LeafGridView GV
Definition: pdelab.hh:1285
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1432
const CON & getCON() const
Definition: pdelab.hh:834
CC & getCC()
Definition: pdelab.hh:933
CGFEMBase< GV, ctype, N, degree, dim, gt > FEMB
Definition: pdelab.hh:589
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:679
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:910
const GO & getGO() const
Definition: pdelab.hh:1527
T::LeafGridView GV
Definition: pdelab.hh:1189
ISTLSolverBackend_ExplicitDiagonal(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1)
Definition: pdelab.hh:1973
Definition: pdelab.hh:391
Definition: pdelab.hh:1665
GFS & getGFS()
Definition: pdelab.hh:1409
void make_consistent(const GFS &gfs, DOF &x) const
Definition: pdelab.hh:836
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1157
const GO * operator->() const
Definition: pdelab.hh:1547
T::ctype ctype
Definition: pdelab.hh:330
UserFunction(const FS &fs_, const Functor &f_)
constructor
Definition: pdelab.hh:1478
GridFunctionSpace< GV, FEM, CON, VBE > GFS
Definition: pdelab.hh:1389
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:1390
Definition: pdelab.hh:88
Definition: pdelab.hh:1505
GFS & getGFS()
Definition: pdelab.hh:1314
GFS & getGFS()
Definition: pdelab.hh:1027
StructuredGrid(Dune::GeometryType::BasicType meshtype, unsigned int cells)
Definition: pdelab.hh:98
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:1010
const CC & getCC() const
Definition: pdelab.hh:1418
FEM & getFEM()
Definition: pdelab.hh:1214
Definition: pdelab.hh:1279
T Grid
Definition: pdelab.hh:1379
const GO & operator*() const
Definition: pdelab.hh:1700
const T * operator->() const
Definition: pdelab.hh:375
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:667
LS * operator->()
Definition: pdelab.hh:1956
Definition: noconstraints.hh:16
PkLocalFiniteElementMap< GV, C, R, degree > FEM
Definition: pdelab.hh:398
CC & getCC()
Definition: pdelab.hh:1320
const LS & getLS() const
Definition: pdelab.hh:1954
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:600
GridFunctionTraits< typename FS::GV, typename FS::NT, 1, FieldVector< typename FS::NT, 1 > > Traits
Definition: pdelab.hh:1475
const T & getGrid() const
Definition: pdelab.hh:355
GO * operator->()
Definition: pdelab.hh:1537
T::ctype ctype
Definition: pdelab.hh:995
OPBLocalFiniteElementMap< ctype, NT, degree, dim, gt > FEM
Definition: pdelab.hh:902
const LS & operator*() const
Definition: pdelab.hh:1957
Backend for sequential BiCGSTAB solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:261
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:1009
ISTLBackend_SEQ_BCGS_SSOR LS
Definition: pdelab.hh:1873
void constraints(const GFS &gfs, CG &cg, const bool verbose=false)
construct constraints
Definition: constraints.hh:751
T Grid
Definition: pdelab.hh:893
const T & operator*() const
Definition: pdelab.hh:178
ISTLSolverBackend_CG_SSOR(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int steps_=5, int verbose_=1)
Definition: pdelab.hh:1823
GalerkinGlobalAssembler(const FS &fs, LOP &lop, const std::size_t nonzeros)
Definition: pdelab.hh:1515
Definition: pdelab.hh:434
const CC & getCC() const
Definition: pdelab.hh:936
GO & getGO()
Definition: pdelab.hh:1679
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:798
YaspGrid< dim > Grid
Definition: pdelab.hh:200
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:1297
T::ctype ctype
Definition: pdelab.hh:1094
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:601
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
T Grid
Definition: pdelab.hh:1284
CONB::CON CON
Definition: pdelab.hh:905
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:673
Definition: pdelab.hh:325
Definition: pdelab.hh:1612
GlobalAssembler(const FSU &fsu, const FSV &fsv, LOP &lop, const std::size_t nonzeros)
Definition: pdelab.hh:1622
LS & getLS()
Definition: pdelab.hh:1804
static const int dimworld
Definition: pdelab.hh:95
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:1296
VBET VBE
Definition: pdelab.hh:906
FEM & getFEM()
Definition: pdelab.hh:1023
void assembleConstraints(const BCTYPE &bctype)
Definition: pdelab.hh:1326
void setNonConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1247
N NT
Definition: pdelab.hh:998
void copyConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:792
N NT
Definition: pdelab.hh:1384
GalerkinGlobalAssemblerNewBackend(const FS &fs, LOP &lop, const MBE &mbe)
Definition: pdelab.hh:1568
const GFS & getGFS() const
Definition: pdelab.hh:633
void set_constrained_dofs(const CG &cg, typename XG::ElementType x, XG &xg)
construct constraints from given boundary condition function
Definition: constraints.hh:798
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1355
Dune::PDELab::ISTLBackend_SEQ_CG_AMG_SSOR< typename ASS::GO > LS
Definition: pdelab.hh:1721
N NT
Definition: pdelab.hh:1193
Backend for sequential conjugate gradient solver with SSOR preconditioner.
Definition: seqistlsolverbackend.hh:348
const CC & getCC() const
Definition: pdelab.hh:1227
const GO & operator*() const
Definition: pdelab.hh:1595
Dune::PDELab::ISTL::BCRSMatrixBackend MBE
Definition: pdelab.hh:1562
GV::Traits::template Codim< 0 >::Entity ElementType
codim 0 entity
Definition: function.hh:117
const LS & operator*() const
Definition: pdelab.hh:1733
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:1104
Backend using (possibly nested) ISTL BCRSMatrices.
Definition: bcrsmatrixbackend.hh:190
const LS * operator->() const
Definition: pdelab.hh:1809
const LS * operator->() const
Definition: pdelab.hh:1885
Backend::Vector< GFS, N > DOF
Definition: pdelab.hh:1103
void clearConstraints()
Definition: pdelab.hh:1427
const LS * operator->() const
Definition: pdelab.hh:1958
GFS::template ConstraintsContainer< N >::Type CC
Definition: pdelab.hh:1105
void make_consistent(const GFS &gfs, DOF &x) const
Definition: pdelab.hh:870
Dune::PDELab::DiscreteGridFunction< GFS, DOF > DGF
Definition: pdelab.hh:719
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1050
T Grid
Definition: pdelab.hh:993
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > &>::type native(T &t)
Definition: backend/interface.hh:192
ISTLBackend_NOVLP_CG_SSORk< typename ASS::GO > LS
Definition: pdelab.hh:1821
CONB::CON CON
Definition: pdelab.hh:1387
VBET VBE
Definition: pdelab.hh:1197
CGCONBase(Grid &grid, const BCType &bctype, const GV &gv)
Definition: pdelab.hh:527
leaf of a function tree
Definition: function.hh:575
const CON & getCON() const
Definition: pdelab.hh:851
Definition: pdelab.hh:1792
DGCONBase< st > CONB
Definition: pdelab.hh:1004
T Grid
Definition: pdelab.hh:1092
void copyNonConstrainedDOFS(const DOF &xin, DOF &xout) const
Definition: pdelab.hh:1068
CGFEMBase< ES, ctype, N, degree, dim, gt > FEMB
Definition: pdelab.hh:708
UnstructuredGrid(std::string filename, bool verbose=true, bool insert_boundary_segments=true)
Definition: pdelab.hh:335
GO & getGO()
Definition: pdelab.hh:1574
void copy_constrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:938
void assembleConstraints(const BCTYPE &bctype)
Definition: pdelab.hh:1134
Grid::ctype ctype
Definition: pdelab.hh:201
Grid * operator->()
Definition: pdelab.hh:304
Definition: pdelab.hh:435
DGCONBase< st > CONB
Definition: pdelab.hh:904
const T & getGrid() const
Definition: pdelab.hh:163
const CON & getCON() const
Definition: pdelab.hh:868
P0ParallelConstraints CON
Definition: pdelab.hh:862
void setConstrainedDOFS(DOF &x, NT nt) const
Definition: pdelab.hh:1337
ISTLBackend_OVLP_BCGS_SSORk< typename FS::GFS, typename FS::CC > LS
Definition: pdelab.hh:1921
GFS & getGFS()
Definition: pdelab.hh:927
CGSpace(Grid &grid, const BCType &bctype)
Definition: pdelab.hh:605
ISTLSolverBackend_ExplicitDiagonal(const FS &fs, const ASS &ass, unsigned maxiter_=5000, int verbose_=1)
Definition: pdelab.hh:1948
Nonoverlapping parallel CG solver preconditioned by block SSOR.
Definition: novlpistlsolverbackend.hh:847
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:911
const CC & getCC() const
Definition: pdelab.hh:1131
const FEM & getFEM() const
Definition: pdelab.hh:1119
const GO & getGO() const
Definition: pdelab.hh:1580
T::ctype ctype
Definition: pdelab.hh:93
ISTL::BCRSMatrixBackend MBE
Definition: pdelab.hh:1616
T::ctype ctype
Definition: pdelab.hh:1286
Dune::PDELab::P0LocalFiniteElementMap< ctype, NT, dim > FEM
Definition: pdelab.hh:1385
Definition: pdelab.hh:579
VTKGridFunctionAdapter< DGF > VTKF
Definition: pdelab.hh:1011