[dune-pdelab] Finite Volume Method with slope reconstruction

[Christian Engwer]

Hi Gregor,

> Would I have to call lfs.child(k).bind(*it) for each component?

no, you just bind the root space.

> What does lfs_cache do? Do I need a separate cache for each child of lfs?

lfs_cache caches the data, you always have containers for the whole
local functions space, thus you only need one cache for the root
space. When accessing the entries you have to take case to use the
accumulate methods, so that the indices are correctly computed.

Christian

> 
> Cheers,
> 
> Gregor
> 
> 
> On 02.05.2017 10:20, Marian Piatkowski wrote:
> > Hi Gregor
> > 
> > 
> > On 04/28/2017 03:13 PM, Gregor Corbin wrote:
> > > Hello dear Dune developers,
> > > 
> > > I have a dune-pdelab code that implements a first-order finite
> > > volume method for a hyperbolic system. Currently I am trying to
> > > implement a slope reconstruction to achieve second order accuracy in
> > > space. Therefore, before each stage of the time-stepping scheme, I
> > > need to iterate over the grid and compute the slopes on each cell.
> > > My current approach is to write a limiter class and give it to the
> > > explicit one step method like this:
> > > 
> > > oneStepMethod.apply(t,dt,x_old,x_new,LinearReconstruction);
> > 
> > That's how it is exactly done in one of our Dune modules for the
> > lecture, dune-funcep.
> > Maybe you could have a look at the code in our Git repository
> > - https://parcomp-git.iwr.uni-heidelberg.de/Teaching/dune-funcep
> > especially the file
> > - https://parcomp-git.iwr.uni-heidelberg.de/Teaching/dune-funcep/blob/master/src/common/limiter.hh
> > 
> > I think the second class in this file is just what you are looking for.
> > 
> > Regards,
> > 
> > Marian
> > 
> > > 
> > > In the prestage method of this class, which gets only the solution
> > > vector, i iterate over all the cells on the grid. For each cell I
> > > want to extract the solution on that cell and on it's neighbors.
> > > But I have no idea how to do this.
> > > 
> > > I appreciate any help on this problem. Is there maybe even a better
> > > way to implement slope limiting in pdelab?
> > > 
> > > Cheers and thanks in advance,
> > > Gregor Corbin
> > > 
> > > Below is a sketch for the code of the limiter class:
> > > 
> > > template <typename GridView, typename GFS, int n_components>
> > > class LinearReconstruction
> > > {
> > >   enum {dim = GridView::dimension};
> > >   enum {n_comp = n_components};
> > > 
> > > public:
> > >   LinearReconstruction(const GridView &grid_view, const GFS &gfs):
> > >   grid_view_(grid_view),
> > >   gfs_(gfs),
> > >   mapper_(grid_view_),
> > >   slopes_(mapper_.size(),0.0)
> > >   {}
> > > 
> > >   template<typename X>
> > >   void prestage(X& x)
> > >   {
> > >     // reconstruct slopes
> > >     for (const auto& cell : elements(grid_view_) )
> > >     {
> > >        // I would like to do something like in the local operator
> > > methods,
> > >        // but how do I get the lfsu_s ???
> > >        // is there a method to get it from the cell and the grid
> > > function space?
> > > 
> > >        // lfsu_s = gfs_.localFunctionSpace(cell); <- does something
> > > similar exist?
> > >       Dune::FieldVector<double,n_comp> u_s(0.0);
> > >       for (int k = 0; k < n_comp; ++k)
> > >         u_s[k] = x(lfsu_s.child(k),0);
> > > 
> > >       // get the solutions from all the neighbors
> > >       std::vector<Dune::FieldVector<double,n_comp> >
> > > u_n(pow(2,dim),0.0); // only valid for structured cubic grids
> > >       int i = 0;
> > >       for (const auto& is : intersections(grid_view_,cell))
> > >       {
> > >         const auto &neighbor = is.outside();
> > >          // lfsu_n = gfs_.localFunctionSpace(neighbor); <- does
> > > something similar exist?
> > > 
> > >          for (int k = 0; k < n_comp; ++k)
> > >               u_n.at(i)[k] = x(lfsu_n.child(k),0);
> > >          ++i;
> > >       }
> > > 
> > >       // use the info to reconstruct slopes
> > >       slopes_.at(mapper_.index(cell)) =  // some function(u_s, u_n)
> > >     }
> > >   }
> > > 
> > >   template<typename V>
> > >   void poststage(V& v)
> > >   {}
> > > private:
> > >   const GridView &grid_view_;
> > >    const GFS &gfs_;
> > >   using Mapper =
> > > Dune::MultipleCodimMultipleGeomTypeMapper<GridView,Dune::MCMGElementLayout>;
> > > 
> > >   const Mapper mapper_;
> > >   std::vector<Dune::FieldMatrix<double,n_comp,dim> > slopes_;
> > > };
> > > 
> > > 
> > > 
> > > 
> > > 
> > > _______________________________________________
> > > dune-pdelab mailing list
> > > dune-pdelab at dune-project.org
> > > http://lists.dune-project.org/mailman/listinfo/dune-pdelab
> > 
> > 
> > _______________________________________________
> > dune-pdelab mailing list
> > dune-pdelab at dune-project.org
> > http://lists.dune-project.org/mailman/listinfo/dune-pdelab
> 
> 
> _______________________________________________
> dune-pdelab mailing list
> dune-pdelab at dune-project.org
> http://lists.dune-project.org/mailman/listinfo/dune-pdelab
> 

-- 
Prof. Dr. Christian Engwer 
Institut für Numerische und Angewandte Mathematik
Fachbereich Mathematik und Informatik der Universität Münster
Einsteinstrasse 62
48149 Münster

E-Mail  christian.engwer at uni-muenster.de
Telefon +49 251 83-35067
FAX     +49 251 83-32729