SmrCell.h

Go to the documentation of this file.

// -*- C++ -*-

#if !defined(__geom_mesh_simplicial_SmrCell_h__)
#define __geom_mesh_simplicial_SmrCell_h__

#if defined(DEBUG_geom) && !defined(DEBUG_SmrCell)
#define DEBUG_SmrCell
#endif

#include "../../defs.h"

#include "../simplex/topology.h"

BEGIN_NAMESPACE_GEOM


template<class SMR>
class SmrCell {
  //
  // Enumerations.
  //

public:

  enum {M = SMR::M};

  //
  // Public types.
  //

public:

  typedef SMR Mesh;

  typedef typename Mesh::CellIterator CellIterator;
  typedef typename Mesh::CellConstIterator CellConstIterator;

  typedef typename Mesh::NodeIterator NodeIterator;
  typedef typename Mesh::NodeConstIterator NodeConstIterator;
  typedef typename Mesh::Vertex Vertex;
  typedef typename Mesh::Number Number;

  typedef Simplex<M,NodeIterator> NodeIteratorSimplex;
  typedef Simplex<M,CellIterator> CellIteratorSimplex;

  //
  // Private types.
  //

private:

  //
  // Public types.
  //

public:

  typedef typename NodeIteratorSimplex::Face Face;
      
  //
  // Data
  //

private:
  
  NodeIteratorSimplex _nodes;
  CellIteratorSimplex _neighbors;
  mutable int _identifier;
  CellIterator _self;

public:

  //--------------------------------------------------------------------------

  SmrCell() :
    _nodes(),
    _neighbors(),
    _identifier(-1),
    _self(0) {
    _nodes.setEach(NodeIterator(0));
    _neighbors.setEach(CellIterator(0));
  }

  SmrCell(const NodeIterator v0, const NodeIterator v1, 
          const CellIterator c0 = CellIterator(0), 
          const CellIterator c1 = CellIterator(0),
          const int identifier = -1) :
    _nodes(v0, v1),
    _neighbors(c0, c1),
    _identifier(identifier),
    _self(0) {
    LOKI_STATIC_CHECK(M == 1, TheSimplexDimensionMustBe1);
  }
  
  SmrCell(const NodeIterator v0, const NodeIterator v1, 
          const NodeIterator v2, 
          const CellIterator c0 = CellIterator(0), 
          const CellIterator c1 = CellIterator(0), 
          const CellIterator c2 = CellIterator(0),
          const int identifier = -1) :
    _nodes(v0, v1, v2),
    _neighbors(c0, c1, c2),
    _identifier(identifier),
    _self(0) {
    LOKI_STATIC_CHECK(M == 2, TheSimplexDimensionMustBe2);
  }
  
  SmrCell(const NodeIterator v0, const NodeIterator v1, 
          const NodeIterator v2, const NodeIterator v3, 
          const CellIterator c0 = CellIterator(0), 
          const CellIterator c1 = CellIterator(0), 
          const CellIterator c2 = CellIterator(0), 
          const CellIterator c3 = CellIterator(0),
          const int identifier = -1) :
    _nodes(v0, v1, v2, v3),
    _neighbors(c0, c1, c2, c3),
    _identifier(identifier),
    _self(0) {
    LOKI_STATIC_CHECK(M == 3, TheSimplexDimensionMustBe3);
  }
  
  SmrCell(const SmrCell& other) :
    _nodes(other._nodes),
    _neighbors(other._neighbors),
    _identifier(other._identifier),
    _self(other._self)
  {}

  ~SmrCell()
  {}

  //--------------------------------------------------------------------------
      
  SmrCell& 
  operator=(const SmrCell& other) {
    if (&other != this) {
      _nodes = other._nodes;
      _neighbors = other._neighbors;
      _identifier = other._identifier;
      _self = other._self;
    }
    return *this;
  }

  //--------------------------------------------------------------------------


  int
  getIdentifier() const {
    return _identifier;
  }

  CellConstIterator
  getSelf() const {
    return _self;
  }

  NodeConstIterator
  getNode(const int m) const {
    return _nodes[m];
  }

  const NodeIteratorSimplex&
  getNodes() const {
    return _nodes;
  }

  int
  getIndex(NodeConstIterator node) const {
    return _nodes.getVertexIndex(node);
  }

  int
  getIndex(NodeIterator node) const {
    return _nodes.getVertexIndex(node);
  }


  int
  getIndex(const Face& f) const {
    for (int m = 0; m != M + 1; ++m) {
      if (! f.hasVertex(_nodes[m])) {
        return m;
      }
    }
    assert(false);
    return -1;
  }

  bool
  hasNode(NodeConstIterator node) const {
    return _nodes.hasVertex(node);
  }


  bool
  hasNode(NodeConstIterator node, int* m) const {
    return _nodes.hasVertex(node, m);
  }

  bool
  hasIncidentBoundaryFace(NodeConstIterator node) const {
#ifdef DEBUG_SmrCell
    assert(hasNode(node));
#endif
    return hasIncidentBoundaryFace(getIndex(node));
  }

  bool
  hasIncidentBoundaryFace(const int n) const {
#ifdef DEBUG_SmrCell
    assert(0 <= n && n < M + 1);
#endif    
    // For each face incident to the node.
    for (int m = 0; m != M+1; ++m) {
      // Exclude the face opposite this node.
      if (m != n) {
        // If the face is on the boundary.
        if (getNeighbor(m) == 0) {
          return true;
        }
      }
    }
    // We did not encounter any incident boundary faces.
    return false;
  }

  void
  getCentroid(Vertex* centroid) const {
    // CONTINUE: This is the arithmetic mean.  Look up the formula for 
    // the centroid.
    *centroid = _nodes[0]->getVertex();
    for (int m = 1; m != M+1; ++m) {
      *centroid += _nodes[m]->getVertex();
    }
    *centroid /= (M + 1);
  }

  CellConstIterator
  getNeighbor(const int m) const {
    return _neighbors[m];
  }

  bool
  isFaceOnBoundary(const int m) const {
    return getNeighbor(m) == CellConstIterator(0);
  }

  int
  getNumberOfNeighbors() const {
    return M + 1 - int(std::count(_neighbors.begin(), _neighbors.end(), 
                                  CellConstIterator(0)));
  }

  int
  getIndex(CellConstIterator c) const {
    return _neighbors.getVertexIndex(c);
  }

  bool
  hasNeighbor(CellConstIterator c) const {
    return _neighbors.hasVertex(c);
  }


  bool
  hasNeighbor(CellConstIterator c, int* m) const {
    return _neighbors.hasVertex(c, m);
  }

  int
  getMirrorIndex(const int m) const {
    if (_neighbors[m] == 0) {
      return -1;
    }
    return _neighbors[m]->getIndex(_self);
  }


  Number
  computeMinimumEdgeLength(int* a, int* b) const {
    Number x = std::numeric_limits<Number>::max();
    Number d;
    int j;
    // For each edge.
    for (int i = 0; i != M; ++i) {
      for (j = i+1; j != M + 1; ++j) {
        d = geom::computeDistance(_nodes[i]->getVertex(), 
                                  _nodes[j]->getVertex());
        if (d < x) {
          x = d;
          *a = i;
          *b = j;
        }
      }
    }
    return x;
  }
  

  Number
  computeMaximumEdgeLength(int* a, int* b) const {
    Number x = 0;
    Number d;
    int j;
    // For each edge.
    for (int i = 0; i != M; ++i) {
      for (j = i+1; j != M + 1; ++j) {
        d = geom::computeDistance(_nodes[i]->getVertex(), 
                                  _nodes[j]->getVertex());
        if (d > x) {
          x = d;
          *a = i;
          *b = j;
        }
      }
    }
    return x;
  }
  
  Number
  computeMinimumEdgeLength() const {
    int a, b;
    return computeMinimumEdgeLength(&a, &b);
  }
  
  Number
  computeMaximumEdgeLength() const {
    int a, b;
    return computeMaximumEdgeLength(&a, &b);
  }

  void
  getSimplex(Simplex<M,Vertex>* simplex) {
    for (int m = 0; m != M + 1; ++m) {
      (*simplex)[m] = getNode(m)->getVertex();
    }
  }
 
  //--------------------------------------------------------------------------


  void
  setIdentifier(const int identifier) const { 
    _identifier = identifier;
  }

  CellIterator
  getSelf() {
    return _self;
  }

  void
  setSelf(const CellIterator self) {
    _self = self;
  }

  NodeIterator
  getNode(const int m) {
    return _nodes[m];
  }

  void
  setNode(const int m, const NodeIterator node) {
    _nodes[m] = node;
  }

  CellIterator
  getNeighbor(const int m) {
    return _neighbors[m];
  }

  void
  setNeighbor(const int m, const CellIterator c) {
    _neighbors[m] = c;
  }

  NodeIterator
  getMirrorVertex(const int m) {
    if (_neighbors[m] == 0) {
      return NodeIterator(0);
    }
    return _neighbors[m]->getNode(getMirrorIndex(m));
  }

  void
  getFace(const int m, Face* f) {
    _nodes.getFace(m, f);
  }

  void
  unlink() {
    // For each vertex.
    for (typename NodeIteratorSimplex::Iterator i = _nodes.getBeginning();
          i != _nodes.getEnd(); ++i) {
      // Remove the vertices' link to this cell.
      (*i)->removeCell(getSelf());
      // Remove the link to the vertex.
      *i = 0;
    }
    // For each neigbor.
    for (typename CellIteratorSimplex::Iterator i = _neighbors.getBeginning();
          i != _neighbors.getEnd(); ++i) {
      // If there is a neighbor.
      if (*i != 0) {
        // Remove the neighbor's link to this cell.
        (*i)->removeNeighbor(getSelf());
        // Remove the link to the neighbor.
        *i = 0;
      }
    }
  }


  void
  removeNeighbor(const CellIterator c) {
    // For each neigbor.
    typename CellIteratorSimplex::Iterator i = _neighbors.getBeginning();
    for (; i != _neighbors.getEnd(); ++i) {
      // If this is the specified neighbor.
      if (*i == c) {
        // Remove the link to the neighbor.
        *i = 0;
        break;
      }
    }
    // The specified cell must be a neighbor.
    assert(i != _neighbors.getEnd());
  }

  void
  negate() {
    _nodes.negate();
    _neighbors.negate();
  }

  //--------------------------------------------------------------------------

  void
  put(std::ostream& out) const {
    // Write the simplex identifier.
    out << "Id = "<< getIdentifier() << " Vertices = ";
    // Write the incident nodes.
    for (int m = 0; m != M+1; ++m) {
      assert(_nodes[m] != 0);
      out << _nodes[m]->getIdentifier() << " ";
    }
    out << " Neighbors = ";
    // Write the neighbors.
    for (int m = 0; m != M+1; ++m) {
      if (_neighbors[m] != 0) {
        out << _neighbors[m]->getIdentifier() << " ";
      }
      else {
        out << "-1 ";
      }
    }
    out << "\n";
  }

};



template<class SMR>
bool
doesCellHaveIncidentFaceOnBoundary(const typename SMR::CellConstIterator& c, 
                                   int i, int j);



template<class SMR>
bool
isOnBoundary(const typename SMR::CellConstIterator& c, int i, int j);



template<class SMR>
inline
bool
isOnBoundary(const typename SMR::Edge& edge) {
  return isOnBoundary<SMR>(edge.first, edge.second, edge.third);
}



template<class SMR>
inline
bool
hasFace(typename SMR::CellConstIterator cell, 
        typename SMR::Cell::Face& face,
        int* faceIndex) {
  return hasFace(cell->getNodes(), face, faceIndex);
}



template<class SMR>
inline
bool
hasFace(typename SMR::CellConstIterator cell, 
        typename SMR::Cell::Face& face) {
  int faceIndex;
  return hasFace<SMR>(cell, face, &faceIndex);
}



template<class SMR>
inline
bool
hasFace(typename SMR::CellConstIterator cell, 
        typename SMR::NodeConstIterator a,
        typename SMR::NodeConstIterator b,
        typename SMR::NodeConstIterator c,
        int* faceIndex) {
  return hasFace(cell->getNodes(), a, b, c, faceIndex);
}



template<class SMR>
inline
bool
hasFace(typename SMR::CellConstIterator cell, 
        typename SMR::NodeConstIterator a,
        typename SMR::NodeConstIterator b,
        typename SMR::NodeConstIterator c) {
  int dummy;
  return hasFace<SMR>(cell, a, b, c, &dummy);
}



template<class SMR>
inline
bool
hasFace(typename SMR::CellIterator cell, 
        typename SMR::NodeIterator a,
        typename SMR::NodeIterator b,
        typename SMR::NodeIterator c,
        int* faceIndex) {
  return hasFace(cell->getNodes(), a, b, c, faceIndex);
}



template<class SMR>
inline
bool
hasFace(typename SMR::CellIterator cell, 
        typename SMR::NodeIterator a,
        typename SMR::NodeIterator b,
        typename SMR::NodeIterator c) {
  int dummy;
  return hasFace<SMR>(cell, a, b, c, &dummy);
}


// Below are the even and odd permutations of the nodes.  We use even 
// permutations to get the next node and odd permutations to get the 
// previous node.
//
// Even Odd
// 0123 0132
// 0231 0213
// 0312 0321
// 1032 1032
// 1203 1230
// 1320 1302
// 2013 2031
// 2130 2103
// 2301 2310
// 3021 3012
// 3102 3120
// 3210 3201


// CONTINUE: Make this a singleton.  Static data may be a bad idea.
// Also, I should move this to Simplex.
inline
int
getNextNodeIndex(const int i, const int j) {
  const int Size = 12;
  static int even[Size][3] = 
    {{0,1,2},
     {0,2,3},
     {0,3,1},
     {1,0,3},
     {1,2,0},
     {1,3,2},
     {2,0,1},
     {2,1,3},
     {2,3,0},
     {3,0,2},
     {3,1,0},
     {3,2,1}};
  for (int n = 0; n != Size; ++n) {
    if (even[n][0] == i && even[n][1] == j) {
      return even[n][2];
    }
  }
  // Here we check that i and j have sensible values.
  assert(false);
  return -1;
}


// CONTINUE: Make this a singleton.  Static data may be a bad idea.
inline
int
getPreviousNodeIndex(const int i, const int j) {
  const int Size = 12;
  static int odd[Size][3] = 
    {{0,1,3},
     {0,2,1},
     {0,3,2},
     {1,0,2},
     {1,2,3},
     {1,3,0},
     {2,0,3},
     {2,1,0},
     {2,3,1},
     {3,0,1},
     {3,1,2},
     {3,2,0}};
  for (int n = 0; n != Size; ++n) {
    if (odd[n][0] == i && odd[n][1] == j) {
      return odd[n][2];
    }
  }
  // Here we check that i and j have sensible values.
  assert(false);
  return -1;
}


template <typename SMR>
inline
int
getNextNodeIndex(const typename SMR::CellConstIterator cell, 
                 const typename SMR::NodeConstIterator a, 
                 const typename SMR::NodeConstIterator b) {
  assert(cell->hasNode(a) && cell->hasNode(b));
  return getNextNodeIndex(cell->getIndex(a), cell->getIndex(b));
}


template <typename SMR>
inline
int
getPreviousNodeIndex(const typename SMR::CellConstIterator cell, 
                     const typename SMR::NodeConstIterator a, 
                     const typename SMR::NodeConstIterator b) {
  assert(cell->hasNode(a) && cell->hasNode(b));
  return getPreviousNodeIndex(cell->getIndex(a), cell->getIndex(b));
}


END_NAMESPACE_GEOM

#define __geom_mesh_simplicial_SmrCell_ipp__
#include "SmrCell.ipp"
#undef __geom_mesh_simplicial_SmrCell_ipp__

#endif

---