uibhmTopology.h

Go to the documentation of this file.

/************************************************************************
**
** Copyright (C) 2014 by Carlos Augusto Teixera Mendes
** All rights reserved.
**
** This file is part of the "GeMA" software. It's use should respect
** the terms in the license agreement that can be found together
** with this source code.
** It is provided AS IS, with NO WARRANTY OF ANY KIND,
** INCLUDING THE WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR
** A PARTICULAR PURPOSE.
**
************************************************************************/

#ifndef _UIBHM_TOPOLOGY_H_
#define _UIBHM_TOPOLOGY_H_

#include <vector>
#include <assert.h>

#include "gmpGemaMeshConfig.h"

#include <gmLog.h>

class GmCellMesh;
template <class T> class UibhmQuery;
template <class T> class UibhmSurfaceQuery;
template <class T> class UibhmSolidQuery;

#define GMP_UIBHM_SIDE_BITS 4

#define GM_UIBHM_SIDE_MASK (~(~0u >> GMP_UIBHM_SIDE_BITS))

#define GMP_UIBHM_MAX_SIDE_VALUE ((1 << GMP_UIBHM_SIDE_BITS) - 1)

#define GMP_UIBHM_BORDER_SIDE  GMP_UIBHM_MAX_SIDE_VALUE

#define GMP_UIBHM_MID_SIDE (GMP_UIBHM_MAX_SIDE_VALUE - 1)

#define GMP_UIBHM_MAX_NUM_SIDES GMP_UIBHM_MID_SIDE

#define GMP_UIBHM_CELL_BITS (32 - GMP_UIBHM_SIDE_BITS)

#define GM_UIBHM_CELL_MASK (~GM_UIBHM_SIDE_MASK)

#define GMP_UIBHM_MAX_NUM_CELLS (1 << GMP_UIBHM_CELL_BITS)

#define GMP_UIBHM_INVALID_VALUE (~0u)


// We don't want the topology structure to be able to store less cells than the mesh 
// representation
static_assert(GMP_UIBHM_CELL_BITS >= GMP_GEMAMESH_CELL_ID_BITS, "Unexpected number of cell id bits");


class GMP_GEMAMESH_API_EXPORT UibhmTopologyBase
{
public:
  virtual ~UibhmTopologyBase() {}

  virtual bool init(bool buildEdgeMap) = 0;
  virtual bool addNewCells() = 0;
  virtual bool cellIsIndexed(int cellId) = 0;
};


template <class Geometry> 
class GMP_GEMAMESH_API_EXPORT UibhmTopology : public UibhmTopologyBase, public Geometry
{
public:

  using HalfSide = quint32;

  using HSVT = typename Geometry::template VectorType<HalfSide>;

  using IVT = typename Geometry::template VectorType<int>;

  UibhmTopology(GmCellMesh* mesh)
    : Geometry(mesh), _nhalf(0), _lastAddedCells(0), _hasOffset(false), _needsEH(false) {}

  virtual bool init(bool buildEdgeMap);
  virtual bool addNewCells();

  virtual bool cellIsIndexed(int cellId) { return cellId < _lastAddedCells; }

  void printParameters(const GmLogCategory& logger, GmLogLevel level) const;
  void print          (const GmLogCategory& logger, GmLogLevel level) const;

  size_t printMemoryStatistics(const GmLogCategory& logger, GmLogLevel level) const;
  size_t usedMemory() const;

  static HalfSide hsMake(int cellId, int side)
  {
    assert(cellId >= 0 && cellId < GMP_UIBHM_MAX_NUM_CELLS);
    assert(side   >= 0 && side   < GMP_UIBHM_MAX_NUM_SIDES);
    return (side << GMP_UIBHM_CELL_BITS) | cellId;
  }

  static HalfSide hsMakeSpecial(int id, int side)
  {
    assert(id >= 0 && id < GMP_UIBHM_MAX_NUM_CELLS);
    assert(side == GMP_UIBHM_BORDER_SIDE || side == GMP_UIBHM_MID_SIDE);
    return (side << GMP_UIBHM_CELL_BITS) | id;
  }

  static int hsCell(HalfSide hs) { return hs & GM_UIBHM_CELL_MASK; }

  static int hsSide(HalfSide hs) { return hs >> GMP_UIBHM_CELL_BITS; }

  void hsBreak(HalfSide hs, int* cellId, int* side) const
  {
    // TODO: Check if returning by a std:: pair is more efficient
    *cellId = hsCell(hs);
    *side   = hsSide(hs);
    assert(*cellId >= 0 && *cellId < this->numCells());
    assert(*side   >= 0 && *side   < this->numSides(*cellId));
  }

  static bool hsIsValid(HalfSide hs) { return hsSide(hs) < GMP_UIBHM_MAX_NUM_SIDES; }

  static bool hsIsMidElement(HalfSide hs) { return hsSide(hs) == GMP_UIBHM_MID_SIDE; }

  int linearHs(int cellId, int side) const
  {
    int lhs;
    if constexpr(Geometry::IsHybrid)
      lhs = (_hasOffset ? _Of[cellId] : this->co(cellId)) + side;
    else
      lhs = cellId * this->_typeTemplate->numSides() + side;
    assert(lhs >= 0 && lhs < _nhalf);
    return lhs;
  }

  int linearHs(HalfSide hs) const
  {
    int cellId, side;
    hsBreak(hs, &cellId, &side);
    return linearHs(cellId, side);
  }

  HalfSide oppositeHs(int cellId, int side) const { return _O[linearHs(cellId, side)]; }

  HalfSide oppositeHs(HalfSide hs) const { return _O[linearHs(hs)]; }

  bool isBorderHs(int cellId, int side) const { return hsSide(oppositeHs(cellId, side)) == GMP_UIBHM_BORDER_SIDE; }

  bool isBorderHs(HalfSide hs) const { assert(hsIsValid(hs));  return hsSide(oppositeHs(hs)) == GMP_UIBHM_BORDER_SIDE; }

  bool isNonInterfaceBorderHs(HalfSide hs) const
  {
    assert(hsIsValid(hs));
    if(this->solid())
      return isBorderHs(hs); // TODO: Change this?
    else
    {
      HalfSide o = oppositeHs(hs);
      return (hsSide(o) == GMP_UIBHM_BORDER_SIDE) && (o != GMP_UIBHM_INVALID_VALUE);
    }
  }

  bool isInterfaceBorderHs(HalfSide hs) const
  {
    assert(hsIsValid(hs) && !this->solid());
    return oppositeHs(hs) == GMP_UIBHM_INVALID_VALUE;
  }

  HalfSide firstHs(int cellId) const
  {
    assert(cellId >= 0 && cellId < this->numCells());
    return hsMake(cellId, 0);
  }

  HalfSide nextHs(HalfSide hs) const
  {
    int cellId, side;
    hsBreak(hs, &cellId, &side);
    return hsMake(cellId, (side + 1) % this->numSides(cellId)); 
  }

  HalfSide prevHs(HalfSide hs) const
  {
    int cellId, side;
    hsBreak(hs, &cellId, &side);
    int ncs = this->numSides(cellId);
    return hsMake(cellId, (side + ncs - 1) % ncs);
  }

  HalfSide adjacentHf(HalfSide hf, int faceEdge) const
  {
    assert(this->solid());
    int cellId, face;
    hsBreak(hf, &cellId, &face);
    assert(faceEdge >= 0 && faceEdge < this->tt(cellId)->numFaceVertices(face));
    return hsMake(cellId, this->tt(cellId)->adjacentFaces(face)[faceEdge]);
  }

  HalfSide radialHf(HalfSide hf, int faceEdge) const { return oppositeHs(adjacentHf(hf, faceEdge)); }

 protected:

  int  buildOffsetTable();
  int  updateOffsetTable(int numAddedCells);

  bool buildOppositeTable();
  bool updateOppositeTable(int numAddedCells);

  bool buildVertexHalfSideTable();
  bool updateVertexHalfSideTable(int numAddedCells);

  bool buildEdgeHalfFaceMap();
  bool updateEdgeHalfFaceMap(int numAddedCells);

  bool buildBordersTable();
  bool updateBordersTable(int numAddedCells);

  template <class Key, class MakeKey> void oppositeTableTraverseCells(MakeKey makeKey, int first, int last, QMap <Key, HalfSide>& adjMap);

  void vertexHalfSideTableTraverseCells(int first, int last);
  bool edgeHalfFaceMapTraverseCells(int first, int last);

  bool numNodesDifferentFromNumEdges() const;

  // The auxiliar map used to build the opposite table must be indexed by a key
  // that uniquely identifies a side. Since we want to support large meshes we 
  // can NOT use a simple coding method like just doing v1*(nv^2) + v2*nv + v3
  // for a face identified by v1, v2 & v3 (with a maximum of nv vertices), since 
  // that will result on an overflow, even if storing the result in a 64 bits 
  // integer, with around 2.6 million nodes.
  // For surface meshes we encode the adge with a size_t as proposed by 
  // Geometry::encodeEdge, while for solid meshes we use a vector with 3 numbers
  // filled with the smallest vertices from the face. Using only 3 vertices is possible
  // since we do NOT support non-manifold meshes and so if two elements share a face,
  // both faces should have the exact same vertices (although in different orders) and
  // there can not be a "partial" intersection between faces (other than at an edge or 
  // a vertex).
  struct V3
  {
    int v[3];
    bool operator<(const struct V3& other) const
    {
      return (v[0] < other.v[0]) ||
             (v[0] == other.v[0] && v[1] < other.v[1]) ||
             (v[0] == other.v[0] && v[1] == other.v[1] && v[2] < other.v[2]);
    }
    bool operator==(const struct V3& other) const
    {
      return v[0] == other.v[0] && v[1] == other.v[1] && v[2] == other.v[2];
    }
  };

  static size_t makeSurfaceKey(QVarLengthArray<int, 10>& nodeList)
  {
    assert(nodeList.size() == 2);
    return UibhmTopology<Geometry>::encodeEdge(nodeList[0], nodeList[1]);
  };

  static V3 makeSolidKey(QVarLengthArray<int, 10>& nodeList)
  {
    assert(nodeList.size() > 2);
    qSort(nodeList);
    return V3{ nodeList[0], nodeList[1], nodeList[2] };
  };

  // Let the query classes access our internal structure
  friend class UibhmQuery<Geometry>;
  friend class UibhmSurfaceQuery<Geometry>;
  friend class UibhmSolidQuery<Geometry>;

  int  _nhalf;     
  bool _hasOffset; 
  bool _needsEH;   
  int  _lastAddedCells; 

  HSVT _O;

  IVT _Of; 

  HSVT _VH;

  QVector<HalfSide> _B;

  QMap<size_t, HalfSide> _EH;
};


#endif