vtkStructuredGridConnectivity.h

Go to the documentation of this file.

/*=========================================================================
 
 Program:   Visualization Toolkit
 Module:    vtkStructuredGridConnectivity.h
 
 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
 All rights reserved.
 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
 
 This software is distributed WITHOUT ANY WARRANTY; without even
 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
 PURPOSE.  See the above copyright notice for more information.
 
 =========================================================================*/
#ifndef vtkStructuredGridConnectivity_h
#define vtkStructuredGridConnectivity_h
 
#define VTK_NO_OVERLAP 0
#define VTK_NODE_OVERLAP 1
#define VTK_EDGE_OVERLAP 2
#define VTK_PARTIAL_OVERLAP 3
 
// VTK include directives
#include "vtkAbstractGridConnectivity.h"
#include "vtkFiltersGeometryModule.h" // For export macro
#include "vtkStructuredData.h"        // For data description definitions
#include "vtkStructuredNeighbor.h"    // For Structured Neighbor object definition
 
// C++ include directives
#include <cassert>  // For assert()
#include <iostream> // For cout
#include <map>      // For STL map
#include <utility>  // For STL pair and overloaded relational operators
#include <vector>   // For STL vector
 
// Forward Declarations
class vtkIdList;
class vtkUnsignedCharArray;
class vtkPointData;
class vtkCellData;
class vtkPoints;
 
class VTKFILTERSGEOMETRY_EXPORT vtkStructuredGridConnectivity : public vtkAbstractGridConnectivity
{
public:
  static vtkStructuredGridConnectivity* New();
  vtkTypeMacro(vtkStructuredGridConnectivity, vtkAbstractGridConnectivity);
  void PrintSelf(ostream& os, vtkIndent indent) override;
 
 
  vtkSetVector6Macro(WholeExtent, int);
  vtkGetVector6Macro(WholeExtent, int);
 
 
  vtkGetMacro(DataDimension, int);
 
  void SetNumberOfGrids(const unsigned int N) override;
 
  virtual void RegisterGrid(const int gridID, int extents[6], vtkUnsignedCharArray* nodesGhostArray,
    vtkUnsignedCharArray* cellGhostArray, vtkPointData* pointData, vtkCellData* cellData,
    vtkPoints* gridNodes);
 
  void GetGridExtent(const int gridID, int extent[6]);
 
  void SetGhostedGridExtent(const int gridID, int ext[6]);
 
  void GetGhostedGridExtent(const int gridID, int ext[6]);
 
  void ComputeNeighbors() override;
 
  int GetNumberOfNeighbors(const int gridID)
  {
    return (static_cast<int>(this->Neighbors[gridID].size()));
  }
 
  vtkStructuredNeighbor GetGridNeighbor(const int gridID, const int nei);
 
  vtkIdList* GetNeighbors(const int gridID, int* extents);
 
  void FillGhostArrays(
    const int gridID, vtkUnsignedCharArray* nodesArray, vtkUnsignedCharArray* cellsArray) override;
 
  void CreateGhostLayers(const int N = 1) override;
 
protected:
  vtkStructuredGridConnectivity();
  ~vtkStructuredGridConnectivity() override;
 
  bool InBounds(const int idx, const int Lo, const int Hi) { return ((idx >= Lo) && (idx <= Hi)); }
 
  bool StrictlyInsideBounds(const int idx, const int Lo, const int Hi)
  {
    return ((idx > Lo) && (idx < Hi));
  }
 
  bool IsSubset(int A[2], int B[2])
  {
    return (this->InBounds(A[0], B[0], B[1]) && this->InBounds(A[1], B[0], B[1]));
  }
 
  int Cardinality(int S[2]) { return (S[1] - S[0] + 1); }
 
 
  int GetNumberOfNodesPerCell(const int dim)
  {
    int numNodes = 0;
    switch (dim)
    {
      case 1:
        numNodes = 2; // line cell
        break;
      case 2:
        numNodes = 4; // quad cell
        break;
      case 3:
        numNodes = 8; // hex cell
        break;
      default:
        assert("ERROR: code should not reach here!" && false);
    } // END switch
    return (numNodes);
  }
 
  void FillNodesGhostArray(const int gridID, const int dataDescription, int GridExtent[6],
    int RealExtent[6], vtkUnsignedCharArray* nodesArray);
 
  void FillCellsGhostArray(const int dataDescription, const int numNodesPerCell, int dims[3],
    int CellExtent[6], vtkUnsignedCharArray* nodesArray, vtkUnsignedCharArray* cellsArray);
 
  void SearchNeighbors(const int gridID, const int i, const int j, const int k, vtkIdList* neiList);
 
  void MarkNodeProperty(const int gridID, const int i, const int j, const int k, int ext[6],
    int RealExtent[6], unsigned char& pfield);
 
  void MarkCellProperty(unsigned char& pfield, unsigned char* nodeGhostFields, const int numNodes);
 
  void GetRealExtent(const int gridID, int GridExtent[6], int RealExtent[6]);
 
  bool IsGhostNode(int GridExtent[6], int RealExtent[6], const int i, const int j, const int k);
 
  bool IsNodeOnBoundaryOfExtent(const int i, const int j, const int k, int ext[6]);
 
  bool IsNodeOnSharedBoundary(
    const int gridID, int RealExtent[6], const int i, const int j, const int k);
 
  bool IsNodeOnBoundary(const int i, const int j, const int k);
 
  bool IsNodeInterior(const int i, const int j, const int k, int GridExtent[6]);
 
  bool IsNodeWithinExtent(const int i, const int j, const int k, int GridExtent[6])
  {
    bool status = false;
 
    switch (this->DataDescription)
    {
      case VTK_X_LINE:
        if ((GridExtent[0] <= i) && (i <= GridExtent[1]))
        {
          status = true;
        }
        break;
      case VTK_Y_LINE:
        if ((GridExtent[2] <= j) && (j <= GridExtent[3]))
        {
          status = true;
        }
        break;
      case VTK_Z_LINE:
        if ((GridExtent[4] <= k) && (k <= GridExtent[5]))
        {
          status = true;
        }
        break;
      case VTK_XY_PLANE:
        if ((GridExtent[0] <= i) && (i <= GridExtent[1]) && (GridExtent[2] <= j) &&
          (j <= GridExtent[3]))
        {
          status = true;
        }
        break;
      case VTK_YZ_PLANE:
        if ((GridExtent[2] <= j) && (j <= GridExtent[3]) && (GridExtent[4] <= k) &&
          (k <= GridExtent[5]))
        {
          status = true;
        }
        break;
      case VTK_XZ_PLANE:
        if ((GridExtent[0] <= i) && (i <= GridExtent[1]) && (GridExtent[4] <= k) &&
          (k <= GridExtent[5]))
        {
          status = true;
        }
        break;
      case VTK_XYZ_GRID:
        if ((GridExtent[0] <= i) && (i <= GridExtent[1]) && (GridExtent[2] <= j) &&
          (j <= GridExtent[3]) && (GridExtent[4] <= k) && (k <= GridExtent[5]))
        {
          status = true;
        }
        break;
      default:
        std::cout << "Data description is: " << this->DataDescription << "\n";
        std::cout.flush();
        assert("pre: Undefined data-description!" && false);
    } // END switch
 
    return (status);
  }
 
  void SetNeighbors(
    const int i, const int j, int i2jOrientation[3], int j2iOrientation[3], int overlapExtent[6]);
 
  void DetermineNeighborOrientation(
    const int idx, int A[2], int B[2], int overlap[2], int orient[3]);
 
  void DetectNeighbors(
    const int i, const int j, int ex1[6], int ex2[6], int orientation[3], int ndim);
 
  int IntervalOverlap(int A[2], int B[2], int overlap[2]);
 
  int DoPartialOverlap(int s[2], int S[2], int overlap[2]);
 
  int PartialOverlap(
    int A[2], const int CardinalityOfA, int B[2], const int CardinalityOfB, int overlap[2]);
 
  void EstablishNeighbors(const int i, const int j);
 
  void AcquireDataDescription();
 
  bool HasBlockConnection(const int gridID, const int blockDirection);
 
  void RemoveBlockConnection(const int gridID, const int blockDirection);
 
  void AddBlockConnection(const int gridID, const int blockDirection);
 
  void ClearBlockConnections(const int gridID);
 
  int GetNumberOfConnectingBlockFaces(const int gridID);
 
  void SetBlockTopology(const int gridID);
 
  void GetIJKBlockOrientation(
    const int i, const int j, const int k, int ext[6], int orientation[3]);
 
  int Get1DOrientation(const int idx, const int ExtentLo, const int ExtentHi, const int OnLo,
    const int OnHi, const int NotOnBoundary);
 
  void CreateGhostedExtent(const int gridID, const int N);
 
  void GetGhostedExtent(
    int* ghostedExtent, int GridExtent[6], const int minIdx, const int maxIdx, const int N);
 
  void CreateGhostedMaskArrays(const int gridID);
 
  void InitializeGhostData(const int gridID);
 
  void AllocatePointData(vtkPointData* RPD, const int N, vtkPointData* PD);
 
  void AllocateCellData(vtkCellData* RCD, const int N, vtkCellData* CD);
 
  void TransferRegisteredDataToGhostedData(const int gridID);
 
  void ComputeNeighborSendAndRcvExtent(const int gridID, const int N);
 
  virtual void TransferGhostDataFromNeighbors(const int gridID);
 
  void TransferLocalNeighborData(const int gridID, const vtkStructuredNeighbor& Neighbor);
 
  void CopyCoordinates(
    vtkPoints* source, vtkIdType sourceIdx, vtkPoints* target, vtkIdType targetIdx);
 
  void CopyFieldData(
    vtkFieldData* source, vtkIdType sourceIdx, vtkFieldData* target, vtkIdType targetIdx);
 
  int GetNeighborIndex(const int gridIdx, const int NeighborGridIdx);
 
  void PrintExtent(int extent[6]);
 
  int DataDimension;
  int DataDescription;
  int WholeExtent[6];
 
  std::vector<int> GridExtents;
  std::vector<int> GhostedExtents;
  std::vector<unsigned char> BlockTopology;
  std::vector<std::vector<vtkStructuredNeighbor>> Neighbors;
  std::map<std::pair<int, int>, int> NeighborPair2NeighborListIndex;
 
private:
  vtkStructuredGridConnectivity(const vtkStructuredGridConnectivity&) = delete;
  void operator=(const vtkStructuredGridConnectivity&) = delete;
};
 
//=============================================================================
//  INLINE METHODS
//=============================================================================
 
//------------------------------------------------------------------------------
inline int vtkStructuredGridConnectivity::GetNeighborIndex(
  const int gridIdx, const int NeighborGridIdx)
{
  assert("pre: Grid index is out-of-bounds!" && (gridIdx >= 0) &&
    (gridIdx < static_cast<int>(this->NumberOfGrids)));
  assert("pre: Neighbor grid index is out-of-bounds!" && (NeighborGridIdx >= 0) &&
    (NeighborGridIdx < static_cast<int>(this->NumberOfGrids)));
 
  std::pair<int, int> gridPair = std::make_pair(gridIdx, NeighborGridIdx);
  assert("pre: Neighboring grid pair does not exist in hash!" &&
    (this->NeighborPair2NeighborListIndex.find(gridPair) !=
      this->NeighborPair2NeighborListIndex.end()));
 
  return (this->NeighborPair2NeighborListIndex[gridPair]);
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::GetGhostedExtent(
  int* ghostedExtent, int GridExtent[6], const int minIdx, const int maxIdx, const int N)
{
  assert("pre: Number of ghost layers must be N >= 1" && (N >= 1));
  assert("pre: ghosted extent pointer is nullptr" && ghostedExtent != nullptr);
 
  ghostedExtent[minIdx] = GridExtent[minIdx] - N;
  ghostedExtent[maxIdx] = GridExtent[maxIdx] + N;
 
  // Clamp the ghosted extent to be within the WholeExtent
  ghostedExtent[minIdx] = (ghostedExtent[minIdx] < this->WholeExtent[minIdx])
    ? this->WholeExtent[minIdx]
    : ghostedExtent[minIdx];
  ghostedExtent[maxIdx] = (ghostedExtent[maxIdx] > this->WholeExtent[maxIdx])
    ? this->WholeExtent[maxIdx]
    : ghostedExtent[maxIdx];
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::SetGhostedGridExtent(const int gridID, int ext[6])
{
  assert("pre: gridID is out-of-bounds" && (gridID >= 0) &&
    (gridID < static_cast<int>(this->NumberOfGrids)));
  assert("pre: ghosted-extents vector has not been allocated" &&
    (this->NumberOfGrids == this->GhostedExtents.size() / 6));
 
  for (int i = 0; i < 6; ++i)
  {
    this->GhostedExtents[gridID * 6 + i] = ext[i];
  }
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::GetGridExtent(const int gridID, int ext[6])
{
  assert("pre: gridID out-of-bounds!" &&
    (gridID >= 0 && gridID < static_cast<int>(this->NumberOfGrids)));
  for (int i = 0; i < 6; ++i)
  {
    ext[i] = this->GridExtents[gridID * 6 + i];
  }
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::GetGhostedGridExtent(const int gridID, int ext[6])
{
  assert("pre: gridID out-of-bounds!" &&
    (gridID >= 0 && gridID < static_cast<int>(this->NumberOfGrids)));
 
  if (this->GhostedExtents.empty())
  {
    ext[0] = ext[2] = ext[4] = -1;
    ext[1] = ext[3] = ext[5] = 0;
    vtkErrorMacro("No ghosted extents found for registered grid extends!!!");
    return;
  }
 
  assert("GhostedExtents are not aligned with registered grid extents" &&
    (this->GhostedExtents.size() == this->GridExtents.size()));
  for (int i = 0; i < 6; ++i)
  {
    ext[i] = this->GhostedExtents[gridID * 6 + i];
  }
}
 
//------------------------------------------------------------------------------
inline bool vtkStructuredGridConnectivity::IsNodeOnBoundaryOfExtent(
  const int i, const int j, const int k, int ext[6])
{
  if (!this->IsNodeWithinExtent(i, j, k, ext))
  {
    return false;
  }
 
  bool status = false;
  switch (this->DataDescription)
  {
    case VTK_X_LINE:
      if (i == ext[0] || i == ext[1])
      {
        status = true;
      }
      break;
    case VTK_Y_LINE:
      if (j == ext[2] || j == ext[3])
      {
        status = true;
      }
      break;
    case VTK_Z_LINE:
      if (k == ext[4] || k == ext[5])
      {
        status = true;
      }
      break;
    case VTK_XY_PLANE:
      if ((i == ext[0] || i == ext[1]) || (j == ext[2] || j == ext[3]))
      {
        status = true;
      }
      break;
    case VTK_YZ_PLANE:
      if ((j == ext[2] || j == ext[3]) || (k == ext[4] || k == ext[5]))
      {
        status = true;
      }
      break;
    case VTK_XZ_PLANE:
      if ((i == ext[0] || i == ext[1]) || (k == ext[4] || k == ext[5]))
      {
        status = true;
      }
      break;
    case VTK_XYZ_GRID:
      if ((i == ext[0] || i == ext[1]) || (j == ext[2] || j == ext[3]) ||
        (k == ext[4] || k == ext[5]))
      {
        status = true;
      }
      break;
    default:
      std::cout << "Data description is: " << this->DataDescription << "\n";
      std::cout.flush();
      assert("pre: Undefined data-description!" && false);
  } // END switch
 
  return (status);
}
 
//------------------------------------------------------------------------------
inline bool vtkStructuredGridConnectivity::IsNodeInterior(
  const int i, const int j, const int k, int GridExtent[6])
{
  bool status = false;
 
  switch (this->DataDescription)
  {
    case VTK_X_LINE:
      if ((GridExtent[0] < i) && (i < GridExtent[1]))
      {
        status = true;
      }
      break;
    case VTK_Y_LINE:
      if ((GridExtent[2] < j) && (j < GridExtent[3]))
      {
        status = true;
      }
      break;
    case VTK_Z_LINE:
      if ((GridExtent[4] < k) && (k < GridExtent[5]))
      {
        status = true;
      }
      break;
    case VTK_XY_PLANE:
      if ((GridExtent[0] < i) && (i < GridExtent[1]) && (GridExtent[2] < j) && (j < GridExtent[3]))
      {
        status = true;
      }
      break;
    case VTK_YZ_PLANE:
      if ((GridExtent[2] < j) && (j < GridExtent[3]) && (GridExtent[4] < k) && (k < GridExtent[5]))
      {
        status = true;
      }
      break;
    case VTK_XZ_PLANE:
      if ((GridExtent[0] < i) && (i < GridExtent[1]) && (GridExtent[4] < k) && (k < GridExtent[5]))
      {
        status = true;
      }
      break;
    case VTK_XYZ_GRID:
      if ((GridExtent[0] < i) && (i < GridExtent[1]) && (GridExtent[2] < j) &&
        (j < GridExtent[3]) && (GridExtent[4] < k) && (k < GridExtent[5]))
      {
        status = true;
      }
      break;
    default:
      std::cout << "Data description is: " << this->DataDescription << "\n";
      std::cout.flush();
      assert("pre: Undefined data-description!" && false);
  } // END switch
 
  return (status);
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::DetermineNeighborOrientation(
  const int idx, int A[2], int B[2], int overlap[2], int orient[3])
{
  assert("pre: idx is out-of-bounds" && (idx >= 0) && (idx < 3));
 
  // A. Non-overlapping cases
  if (overlap[0] == overlap[1])
  {
    if (A[1] == B[0])
    {
      orient[idx] = vtkStructuredNeighbor::HI;
    }
    else if (A[0] == B[1])
    {
      orient[idx] = vtkStructuredNeighbor::LO;
    }
    else
    {
      orient[idx] = vtkStructuredNeighbor::UNDEFINED;
      assert("ERROR: Code should not reach here!" && false);
    }
  } // END non-overlapping cases
  // B. Sub-set cases
  else if (this->IsSubset(A, B))
  {
    if ((A[0] == B[0]) && (A[1] == B[1]))
    {
      orient[idx] = vtkStructuredNeighbor::ONE_TO_ONE;
    }
    else if (this->StrictlyInsideBounds(A[0], B[0], B[1]) &&
      this->StrictlyInsideBounds(A[1], B[0], B[1]))
    {
      orient[idx] = vtkStructuredNeighbor::SUBSET_BOTH;
    }
    else if (A[0] == B[0])
    {
      orient[idx] = vtkStructuredNeighbor::SUBSET_HI;
    }
    else if (A[1] == B[1])
    {
      orient[idx] = vtkStructuredNeighbor::SUBSET_LO;
    }
    else
    {
      orient[idx] = vtkStructuredNeighbor::UNDEFINED;
      assert("ERROR: Code should not reach here!" && false);
    }
  }
  // C. Super-set cases
  else if (this->IsSubset(B, A))
  {
    orient[idx] = vtkStructuredNeighbor::SUPERSET;
  }
  // D. Partially-overlapping (non-subset) cases
  else if (!(this->IsSubset(A, B) || this->IsSubset(A, B)))
  {
    if (this->InBounds(A[0], B[0], B[1]))
    {
      orient[idx] = vtkStructuredNeighbor::LO;
    }
    else if (this->InBounds(A[1], B[0], B[1]))
    {
      orient[idx] = vtkStructuredNeighbor::HI;
    }
    else
    {
      orient[idx] = vtkStructuredNeighbor::UNDEFINED;
      assert("ERROR: Code should not reach here!" && false);
    }
  }
  else
  {
    orient[idx] = vtkStructuredNeighbor::UNDEFINED;
    assert("ERROR: Code should not reach here!" && false);
  }
}
 
//------------------------------------------------------------------------------
inline int vtkStructuredGridConnectivity::Get1DOrientation(const int idx, const int ExtentLo,
  const int ExtentHi, const int OnLo, const int OnHi, const int NotOnBoundary)
{
  if (idx == ExtentLo)
  {
    return OnLo;
  }
  else if (idx == ExtentHi)
  {
    return OnHi;
  }
  return NotOnBoundary;
}
 
//------------------------------------------------------------------------------
inline bool vtkStructuredGridConnectivity::HasBlockConnection(
  const int gridID, const int blockDirection)
{
  // Sanity check
  assert("pre: gridID is out-of-bounds" && (gridID >= 0) &&
    (gridID < static_cast<int>(this->NumberOfGrids)));
  assert("pre: BlockTopology has not been properly allocated" &&
    (this->NumberOfGrids == this->BlockTopology.size()));
  assert("pre: blockDirection is out-of-bounds" && (blockDirection >= 0) && (blockDirection < 6));
  bool status = false;
  if (this->BlockTopology[gridID] & (1 << blockDirection))
  {
    status = true;
  }
  return (status);
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::RemoveBlockConnection(
  const int gridID, const int blockDirection)
{
  // Sanity check
  assert("pre: gridID is out-of-bounds" && (gridID >= 0) &&
    (gridID < static_cast<int>(this->NumberOfGrids)));
  assert("pre: BlockTopology has not been properly allocated" &&
    (this->NumberOfGrids == this->BlockTopology.size()));
  assert("pre: blockDirection is out-of-bounds" && (blockDirection >= 0) && (blockDirection < 6));
 
  this->BlockTopology[gridID] &= ~(1 << blockDirection);
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::AddBlockConnection(
  const int gridID, const int blockDirection)
{
  // Sanity check
  assert("pre: gridID is out-of-bounds" && (gridID >= 0) &&
    (gridID < static_cast<int>(this->NumberOfGrids)));
  assert("pre: BlockTopology has not been properly allocated" &&
    (this->NumberOfGrids == this->BlockTopology.size()));
  assert("pre: blockDirection is out-of-bounds" && (blockDirection >= 0) && (blockDirection < 6));
  this->BlockTopology[gridID] |= (1 << blockDirection);
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::ClearBlockConnections(const int gridID)
{
  // Sanity check
  assert("pre: gridID is out-of-bounds" && (gridID >= 0) &&
    (gridID < static_cast<int>(this->NumberOfGrids)));
  assert("pre: BlockTopology has not been properly allocated" &&
    (this->NumberOfGrids == this->BlockTopology.size()));
  for (int i = 0; i < 6; ++i)
  {
    this->RemoveBlockConnection(gridID, i);
  } // END for all block directions
}
 
//------------------------------------------------------------------------------
inline int vtkStructuredGridConnectivity::GetNumberOfConnectingBlockFaces(const int gridID)
{
  // Sanity check
  assert("pre: gridID is out-of-bounds" && (gridID >= 0) &&
    (gridID < static_cast<int>(this->NumberOfGrids)));
  assert("pre: BlockTopology has not been properly allocated" &&
    (this->NumberOfGrids == this->BlockTopology.size()));
 
  int count = 0;
  for (int i = 0; i < 6; ++i)
  {
    if (this->HasBlockConnection(gridID, i))
    {
      ++count;
    }
  }
  assert("post: count must be in [0,5]" && (count >= 0 && count <= 6));
  return (count);
}
 
//------------------------------------------------------------------------------
inline void vtkStructuredGridConnectivity::SetNumberOfGrids(const unsigned int N)
{
  if (N == 0)
  {
    vtkErrorMacro("Number of grids cannot be 0.");
    return;
  }
 
  this->NumberOfGrids = N;
  this->AllocateUserRegisterDataStructures();
 
  this->GridExtents.resize(6 * N, -1);
  this->Neighbors.resize(N);
  this->BlockTopology.resize(N);
}
#endif /* vtkStructuredGridConnectivity_h */