multi_grid.h

Go to the documentation of this file.

/*
 * Copyright (c) 2009-2015:  G-CSC, Goethe University Frankfurt
 * Author: Sebastian Reiter
 * 
 * This file is part of UG4.
 * 
 * UG4 is free software: you can redistribute it and/or modify it under the
 * terms of the GNU Lesser General Public License version 3 (as published by the
 * Free Software Foundation) with the following additional attribution
 * requirements (according to LGPL/GPL v3 §7):
 * 
 * (1) The following notice must be displayed in the Appropriate Legal Notices
 * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
 * 
 * (2) The following notice must be displayed at a prominent place in the
 * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
 * 
 * (3) The following bibliography is recommended for citation and must be
 * preserved in all covered files:
 * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
 *   parallel geometric multigrid solver on hierarchically distributed grids.
 *   Computing and visualization in science 16, 4 (2013), 151-164"
 * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
 *   flexible software system for simulating pde based models on high performance
 *   computers. Computing and visualization in science 16, 4 (2013), 165-179"
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
 * GNU Lesser General Public License for more details.
 */
 
#ifndef __H__LIB_GRID__MULTI_GRID__
#define __H__LIB_GRID__MULTI_GRID__
 
#include <vector>
#include <cassert>
#include "grid/grid.h"
#include "tools/subset_handler_grid.h"
#include "common_attachments.h"
#include "common/util/array_util.h"
#include "multi_grid_child_info.h"
#include "algorithms/attachment_util.h"
 
namespace ug
{
 
 
class MultiGrid : public Grid, public GridObserver
{
  friend struct MGVertexInfo;
  friend struct MGEdgeInfo;
  friend struct MGFaceInfo;
  friend struct MGVolumeInfo;
 
  protected:
    typedef MGVertexInfo VertexInfo;
    typedef MGEdgeInfo EdgeInfo;
    typedef MGFaceInfo FaceInfo;
    typedef MGVolumeInfo VolumeInfo;
 
  public:
    //  methods from Grid, that would be hidden if not explicitly
    //  declared as required.
    using Grid::begin;
    using Grid::end;
    using Grid::num;
    using Grid::get_grid_objects;
    using Grid::create;
    using Grid::create_by_cloning;
 
 
    MultiGrid();
 
    MultiGrid(uint options);
    
    virtual ~MultiGrid();
 
    void enable_hierarchical_insertion(bool bEnable);
    inline bool hierarchical_insertion_enabled()  {return m_bHierarchicalInsertion;}
 
  //  element creation
 
    template<class TGeomObj>
    typename geometry_traits<TGeomObj>::iterator
    create(size_t level);
 
    template <class TGeomObj>
    typename geometry_traits<TGeomObj>::iterator
    create(const typename geometry_traits<TGeomObj>::Descriptor& descriptor,
        size_t level);
    VertexIterator create_by_cloning(Vertex* pCloneMe,
                      int level);
 
    EdgeIterator create_by_cloning(Edge* pCloneMe,
                       const EdgeVertices& ev,
                       int level);
 
    FaceIterator create_by_cloning(Face* pCloneMe,
                     const FaceVertices& fv,
                     int level);
 
    VolumeIterator create_by_cloning(Volume* pCloneMe,
                     const VolumeVertices& vv,
                     int level);
 
    inline size_t num_levels() const  {return (size_t)m_hierarchy.num_subsets();}
 
    inline size_t top_level() const;
 
    inline void level_required(int lvl);
 
    template <class TElem> inline
    size_t num(int level) const   {return m_hierarchy.num<TElem>(level);}
 
    template <class TElem> inline
    typename geometry_traits<TElem>::iterator
    begin(int level)
    {
      //assert(level < (int)num_levels() && "ERROR in MultiGrid::begin(...): requested level too high!");
      if(level >= (int)num_levels())
        return end<TElem>();
      return m_hierarchy.begin<TElem>(level);
    }
 
    template <class TElem> inline
    typename geometry_traits<TElem>::iterator
    end(int level)
    {
      //assert(level < (int)num_levels() && "ERROR in MultiGrid::end(...): requested level too high!");
      if(level >= (int)num_levels())
        return end<TElem>();
      return m_hierarchy.end<TElem>(level);
    }
 
    template <class TElem> inline
    typename geometry_traits<TElem>::const_iterator
    begin(int level) const
    {
      //assert(level < (int)num_levels() && "ERROR in MultiGrid::begin(...): requested level too high!");
      if(level >= (int)num_levels())
        return end<TElem>();
      return m_hierarchy.begin<TElem>(level);
    }
 
    template <class TElem> inline
    typename geometry_traits<TElem>::const_iterator
    end(int level) const
    {
      //assert(level < (int)num_levels() && "ERROR in MultiGrid::end(...): requested level too high!");
      if(level >= (int)num_levels())
        return end<TElem>();
      return m_hierarchy.end<TElem>(level);
    }
 
  //  geometric-object-collection
    inline GridObjectCollection
    get_grid_objects(int level)
    {return m_hierarchy.get_grid_objects_in_subset(level);}
    
  //  multi-level-geometric-object-collection
    virtual GridObjectCollection get_grid_objects()
    {return m_hierarchy.get_grid_objects();}
    
    template <class TElem> inline
    int get_level(TElem* elem) const
    {return m_hierarchy.get_subset_index(elem);}
 
    GridObject* get_parent(GridObject* parent) const;
    inline GridObject* get_parent(Vertex* o) const  {return get_info(o).m_pParent;}
    inline GridObject* get_parent(Edge* o) const  {return get_info(o).m_pParent;}
    inline GridObject* get_parent(Face* o) const    {return m_aaParentFACE[o];}
    inline GridObject* get_parent(Volume* o) const    {return m_aaParentVOL[o];}
 
  //  number of children
    template <class TElem> inline
    bool has_children(TElem* elem) const
    {return get_info(elem).has_children();}
 
  //  CHILD QUANTITIES
 
    template <class TChild, class TElem>
    inline size_t num_children(TElem* elem) const {return num_children(elem, TChild());}
 
    template <class TChild>
    size_t num_children(GridObject* elem) const;
 
    template <class TElem>
    inline size_t num_children_total(TElem* elem) const;
 
    template <class TElem>
    inline size_t num_child_vertices(TElem* elem) const {return get_info(elem).num_child_vertices();}
 
 
    template <class TElem>
    inline size_t num_child_edges(TElem* elem) const  {return get_info(elem).num_child_edges();}
    inline size_t num_child_edges(Vertex*) const  {return 0;}
 
    template <class TElem>
    inline size_t num_child_faces(TElem* elem) const  {return get_info(elem).num_child_faces();}
    inline size_t num_child_faces(Vertex*) const  {return 0;}
    inline size_t num_child_faces(Edge*) const    {return 0;}
 
    inline size_t num_child_volumes(Volume* elem) const {return get_info(elem).num_child_volumes();}
    template <class TElem>
    inline size_t num_child_volumes(TElem*) const   {return 0;}
  //  CHILD ACCESS
 
    template <class TChild, class TElem>
    inline TChild* get_child(TElem* elem, size_t ind) const {return get_child(elem, ind, TChild());}
 
    template <class TChild>
    TChild* get_child(GridObject* elem, size_t ind) const;
    template <class TElem>
    inline Vertex* get_child_vertex(TElem* elem) const  {return get_info(elem).child_vertex();}
 
 
    template <class TElem>
    inline Edge* get_child_edge(TElem* elem, size_t ind) const  {return get_info(elem).child_edge(ind);}
    inline Edge* get_child_edge(Vertex*, size_t) const    {return NULL;}
 
    template <class TElem>
    inline Face* get_child_face(TElem* elem, size_t ind) const  {return get_info(elem).child_face(ind);}
    inline Face* get_child_face(Vertex*, size_t) const    {return NULL;}
    inline Face* get_child_face(Edge*, size_t) const      {return NULL;}
 
    inline Volume* get_child_volume(Volume* elem, size_t ind) const {return get_info(elem).child_volume(ind);}
    template <class TElem>
    inline Volume* get_child_volume(TElem*, size_t) const {return NULL;}
 
    template <class TElem>
    void clear_child_connections(TElem* parent);
 
 
    template <class TElem>
    void associate_parent(TElem* elem, GridObject* parent);
 
    template <class TElem>
    char parent_type(TElem* elem) const;
 
 
    template <class TElem>
    void set_parent_type(TElem* elem, char type);
 
    void check_edge_elem_infos(int level) const;
    void check_face_elem_infos(int level) const;
    void check_volume_elem_infos(int level) const;
    
 
    SubsetHandler& get_hierarchy_handler()    {return m_hierarchy;}
    
  //  Don't invoke the following methods directly!
  //  They are intended for internal feedback only.
 
  //  grid callbacks
    virtual void elements_to_be_cleared(Grid* grid);
 
    virtual void vertex_created(Grid* grid, Vertex* vrt,
                  GridObject* pParent = NULL,
                  bool replacesParent = false);
 
    virtual void edge_created(Grid* grid, Edge* e,
                  GridObject* pParent = NULL,
                  bool replacesParent = false);
 
    virtual void face_created(Grid* grid, Face* f,
                  GridObject* pParent = NULL,
                  bool replacesParent = false);
 
    virtual void volume_created(Grid* grid, Volume* vol,
                  GridObject* pParent = NULL,
                  bool replacesParent = false);
 
    virtual void vertex_to_be_erased(Grid* grid, Vertex* vrt,
                     Vertex* replacedBy = NULL);
 
    virtual void edge_to_be_erased(Grid* grid, Edge* e,
                     Edge* replacedBy = NULL);
 
    virtual void face_to_be_erased(Grid* grid, Face* f,
                     Face* replacedBy = NULL);
 
    virtual void volume_to_be_erased(Grid* grid, Volume* vol,
                     Volume* replacedBy = NULL);
 
  protected:
 
  //  Note: VertexInfo and EdgeInfo are stored directly, FaceInfo and
  //  VolumeInfo are stored dynamically.
    typedef Attachment<GridObject*> AParent;
    typedef Attachment<VertexInfo>      AVertexInfo;
    typedef Attachment<EdgeInfo>      AEdgeInfo;
    typedef Attachment<FaceInfo*>     AFaceInfo;
    typedef Attachment<VolumeInfo*>     AVolumeInfo;
    typedef Attachment<char>        AParentType;
 
 
  //  initialization
    void init();
    
  //  create levels
    void create_levels(int numLevels);
 
  //  info-access
    inline VertexInfo& get_info(Vertex* v);
    inline EdgeInfo& get_info(Edge* e);
    inline FaceInfo& get_info(Face* f);
    inline VolumeInfo& get_info(Volume* v);
 
  //  const info-access
    inline const VertexInfo& get_info(Vertex* v) const;
    inline const EdgeInfo& get_info(Edge* e) const;
    inline const FaceInfo& get_info(Face* f) const;
    inline const VolumeInfo& get_info(Volume* v) const;
 
  //  elem creation
    template <class TElem>
    inline void element_created(TElem* elem)  {element_created<TElem, TElem>(elem, NULL);}
 
    template <class TElem, class TParent>
    void element_created(TElem* elem, TParent* pParent);
 
    template <class TElem, class TParent>
    void element_created(TElem* elem, TParent* pParent, TElem* pReplaceMe);
 
    template <class TElem>
    void element_to_be_erased(TElem* elem);
 
    template <class TElem, class TParent>
    void element_to_be_erased(TElem* elem, TParent* pParent);
 
    //template <class TElem>
    //void element_to_be_replaced(TElem* elemOld, TElem* elemNew);
 
 
 
    template <class TElem>
    inline size_t num_children(TElem* elem, const Vertex&) const
      {return num_child_vertices(elem);}
 
    template <class TElem>
    inline size_t num_children(TElem* elem, const Edge&) const
      {return num_child_edges(elem);}
 
    template <class TElem>
    inline size_t num_children(TElem* elem, const Face&) const
      {return num_child_faces(elem);}
 
    template <class TElem>
    inline size_t num_children(TElem* elem, const Volume&) const
      {return num_child_volumes(elem);}
 
    template <class TElem>
    inline Vertex* get_child(TElem* elem, size_t ind, const Vertex&) const
      {return get_child_vertex(elem);}
 
    template <class TElem>
    inline Edge* get_child(TElem* elem, size_t ind, const Edge&) const
      {return get_child_edge(elem, ind);}
 
    template <class TElem>
    inline Face* get_child(TElem* elem, size_t ind, const Face&) const
      {return get_child_face(elem, ind);}
 
    template <class TElem>
    inline Volume* get_child(TElem* elem, size_t ind, const Volume&) const
      {return get_child_volume(elem, ind);}
 
    inline void set_parent(Vertex* o, GridObject* p)  {get_info(o).m_pParent = p;}
    inline void set_parent(Edge* o, GridObject* p)    {get_info(o).m_pParent = p;}
    inline void set_parent(Face* o, GridObject* p)      {m_aaParentFACE[o] = p;}
    inline void set_parent(Volume* o, GridObject* p)    {m_aaParentVOL[o] = p;}
 
    template <class TParent, class TChild>
    void add_child(TParent* p, TChild* c);
 
    template <class TChild>
    void add_child(GridObject* p, TChild* c);
 
    template <class TParent, class TChild>
    void remove_child(TParent* p, TChild* c);
 
    template <class TChild>
    void remove_child(GridObject* p, TChild* c);
 
    inline void create_child_info(Vertex* o){}
    inline void create_child_info(Edge* o)  {}
    inline void create_child_info(Face* o)    {if(!m_aaFaceInf[o]) m_aaFaceInf[o] = new FaceInfo();}
    inline void create_child_info(Volume* o)  {if(!m_aaVolInf[o]) m_aaVolInf[o] = new VolumeInfo();}
 
    inline void release_child_info(Vertex* o) {}
    inline void release_child_info(Edge* o)   {}
    inline void release_child_info(Face* o)     {if(m_aaFaceInf[o]) delete m_aaFaceInf[o]; m_aaFaceInf[o] = NULL;}
    inline void release_child_info(Volume* o)   {if(m_aaVolInf[o]) delete m_aaVolInf[o]; m_aaVolInf[o] = NULL;}
  //  hierarchy
    SubsetHandler m_hierarchy;
    bool m_bHierarchicalInsertion;
 
  //  parent attachment
    AParent   m_aParent;
 
  //  info attachments
    AVertexInfo m_aVertexInfo;
    AEdgeInfo m_aEdgeInfo;
    AFaceInfo m_aFaceInfo;
    AVolumeInfo m_aVolumeInfo;
    AParentType m_aParentType;
 
  //  parent access - only required for faces and volumes.
    Grid::FaceAttachmentAccessor<AParent>   m_aaParentFACE;
    Grid::VolumeAttachmentAccessor<AParent>   m_aaParentVOL;
 
  //  element info access
    Grid::VertexAttachmentAccessor<AVertexInfo> m_aaVrtInf;
    Grid::EdgeAttachmentAccessor<AEdgeInfo>   m_aaEdgeInf;
    Grid::FaceAttachmentAccessor<AFaceInfo>   m_aaFaceInf;
    Grid::VolumeAttachmentAccessor<AVolumeInfo> m_aaVolInf;
 
    MultiElementAttachmentAccessor<AParentType> m_aaParentType;
};
 
 
 
 
template <class TGrid>
class MGWrapper
{
  public:
    MGWrapper(TGrid& grid);
    
    inline uint num_levels() const;
 
    template <class TElem> inline
    uint num(int level) const;
 
    template <class TElem> inline
    typename geometry_traits<TElem>::iterator
    begin(int level);
 
    template <class TElem> inline
    typename geometry_traits<TElem>::iterator
    end(int level);
};
 
}// end of namespace
 
//  include implementation
#include "multi_grid_impl.hpp"
 
#endif