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