Traverse.h

// ============================================================================
// ==                                                                        ==
// == AMDiS - Adaptive multidimensional simulations                          ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
// ==  crystal growth group                                                  ==
// ==                                                                        ==
// ==  Stiftung caesar                                                       ==
// ==  Ludwig-Erhard-Allee 2                                                 ==
// ==  53175 Bonn                                                            ==
// ==  germany                                                               ==
// ==                                                                        ==
// ============================================================================
// ==                                                                        ==
// ==  http://www.caesar.de/cg/AMDiS                                         ==
// ==                                                                        ==
// ============================================================================

/** \file Traverse.h */

/** \defgroup Traverse Traverse module
 * @{ <img src="traverse.png"> @}
 *
 * \brief
 * Contains classes used for mesh traversal.
 */

#ifndef AMDIS_TRAVERSE_H
#define AMDIS_TRAVERSE_H

#include "Flag.h"
#include "Global.h"
#include "ElInfo.h"
#include "ElInfoStack.h"
#include "MemoryManager.h"
#include "OpenMP.h"
#include <vector>
#include <deque>
#include <stack>

namespace AMDiS {

  class MacroElement;
  class Mesh;

  // ===========================================================================
  // ===== class TraverseStack =================================================
  // ===========================================================================

  /** \ingroup Traverse
   * \brief
   * Mesh refinement and coarsening routines are examples of functions which 
   * need a non-recursive access to the mesh elements.
   * The implementation of the non-recursive mesh traversal routines uses a 
   * stack to save the tree path from a macro element to the current element. 
   * TraverseStack holds such information. Before calling the non-recursive mesh
   * traversal routines, such a stack must be allocated. The stack is 
   * initialized by each call to \ref traverseFirst().
   */
  class TraverseStack 
  {
  public:
    MEMORY_MANAGED(TraverseStack);

    /** \brief
     * Creates an empty TraverseStack
     */
    TraverseStack() 
      : traverse_mel(NULL),
        stack_size(0),
        stack_used(0),
        save_stack_used(0),
        myThreadId_(0),
        maxThreads_(1)
    {
#ifdef _OPENMP
#pragma omp critical
#endif
      id = ElInfo::traverseIdCounter++;
    }

    /** \brief
     * Destructor
     */
    ~TraverseStack() 
    {
      for (int i = 0; i < static_cast<int>(elinfo_stack.size()); i++) {
	DELETE elinfo_stack[i];
      }
      for (int i = 0; i < static_cast<int>(save_elinfo_stack.size()); i++) {
	DELETE save_elinfo_stack[i];
      }
    }

  public:
    /** \brief
     * Returns the first ElInfo of a non recursive traversal which fullfills the
     * selection criterion specified by level and fill_flag and initializes the
     * stack. After a call of traverseFirst the next Elements are traversed via 
     * \ref traverseNext(). 
     */
    ElInfo* traverseFirst(Mesh *mesh, int level, Flag fill_flag);

    /** \brief
     * Returns the next ElInfo in a traversal initiated by \ref traverseFirst()
     * If NULL is returned, the traversal is finished.
     */
    ElInfo* traverseNext(ElInfo* elinfo_old);

    /// Returns the neighbour-th neighbour of elInfoOld
    ElInfo* traverseNeighbour(ElInfo* elInfoOld, int neighbour);

    /// Returns the neighbour-th neighbour of elInfoOld
    ElInfo* traverseNeighbour3d(ElInfo* elInfoOld, int neighbour);

    /// Returns the neighbour-th neighbour of elInfoOld
    ElInfo* traverseNeighbour2d(ElInfo* elInfoOld, int neighbour);

    /// Not yet implemented
    ElInfo* traverseMultiGridLevel();

    /// Preorder traversal of all elements
    ElInfo* traverseEveryElementPreorder();

    /// Inorder traversal of all elements
    ElInfo* traverseEveryElementInorder();

    /// Postorder traversal of all elements
    ElInfo* traverseEveryElementPostorder();

    /// Only for 3d: Calls update of all ElInfo3d onjects in \ref elinfo_stack
    void update();

    ///
    void getCoordsInElem(const ElInfo *upperElInfo, DimMat<double> *coords);

    /// Is used for parallel mesh traverse.
    inline void setMyThreadId(int myThreadId) {
      myThreadId_ = myThreadId;
    }

    /// Is used for parallel mesh traverse.
    inline void setMaxThreads(int maxThreads) {
      maxThreads_ = maxThreads;
    }

    ///
    int getStackData(std::vector<ElInfo*> &elInfos, std::vector<int> &infos) {
      elInfos = elinfo_stack;
      infos = info_stack;

      return stack_used;
    }

  private:
    /// Enlargement of the stack
    void enlargeTraverseStack();

    /// Used by \ref traverseFirst() \ref traverseNext()
    ElInfo* traverseLeafElement();

    /// Used by \ref traverseFirst() \ref traverseNext()
    ElInfo* traverseLeafElementLevel();

    /// Used by \ref traverseFirst() \ref traverseNext()
    ElInfo* traverseElementLevel();

    /** \brief
     * Avoids copy of a traverse stack. If copy should be possible
     * the operator must be implemented (deep copy not flat copy!)
     */
    void operator=(const TraverseStack& /*rhs*/) {
      FUNCNAME("TraverseStack::operator=()");
      ERROR_EXIT("not implemented");
    }

    /** \brief
     * Avoids copy of a traverse stack. If copy should be possible
     * the operator must be implemented (deep copy not flat copy!)
     */
    TraverseStack(const TraverseStack&) {
      FUNCNAME("TraverseStack::TraverseStack()");
      ERROR_EXIT("not implemented");
    }

  private:
    /// Iterator to the current MacroElement
    std::deque<MacroElement*>::const_iterator currentMacro;

    /// Mesh which is currently traversed
    Mesh* traverse_mesh;

    /** \brief
     * Traverse level. Used if CALL_EL_LEVEL or CALL_LEAF_EL_LEVEL are set in
     * \ref traverse_fill_flag
     */
    int traverse_level;

    /** \brief
     * Flags for traversal. Specifies which elements are visited and which
     * information are filled to the ElInfo objects.
     */
    Flag traverse_fill_flag;
  
    /// current macro element
    const MacroElement *traverse_mel;

    /// Current size of the stack
    int stack_size;

    /// Used size of the stack
    int stack_used;

    ///
    std::vector<ElInfo*> elinfo_stack;

    /** \brief
     * Stores for all level, which children of this element was already 
     * visited. If 0, no children were visited (this status occurs only
     * in intermediate computations). If 1, only the first was vistied.
     * If 2, both children of the element on this level were already 
     * visited.
     */
    std::vector<int> info_stack;

    ///
    const MacroElement *save_traverse_mel;

    ///
    std::vector<ElInfo*> save_elinfo_stack;

    ///
    std::vector<unsigned char> save_info_stack;

    ///
    int save_stack_used;

    ///
    int id;

    /** \brief
     * Is used for parallel mesh traverse. The thread with the id
     * myThreadId is only allowed to access coarse elements, which id
     * satisfies: elId % maxThreads = myThreadId.
     */
    int myThreadId_;

    /** \brief
     * Is used for parallel mesh traverse. The thread with the id
     * myThreadId is only allowed to access coarse elements, which id
     * satisfies: elId % maxThreads = myThreadId.
     */
    int maxThreads_;
  };

  // ============================================================================
  // ===== class Traverse =======================================================
  // ============================================================================

  /** \ingroup Traverse
   * \brief
   * Class for recursive mesh traversal. Used in \ref Mesh::traverse()
   */
  class Traverse 
  {
  public:
    MEMORY_MANAGED(Traverse);

    /** \brief
     * Creates a Traverse object for Mesh m with fill flags f and level l.
     * el_fct is a pointer to a the function that will be called for every 
     * visited element.
     */
    Traverse(Mesh* m,
	     const Flag f,
	     const unsigned char l,
	     int (*ef)(ElInfo*))
      : mesh(m),
        flag(f),
        level(l),
        el_fct(ef) 
    {
      TEST_EXIT(m)("No traverse without mesh!\n");
      id = ElInfo::traverseIdCounter++;
    }

    /// Performs the recursive traversal
    int recursive(ElInfoStack*);
  
  private:
    Mesh *mesh;

    Flag flag; 

    int level;

    int (*el_fct)(ElInfo*);

    int id;

    Traverse();

    Traverse(const Traverse&);
  };

}

#endif  // AMDIS_TRAVERSE_H