OctreeNodeInterface.h

Go to the documentation of this file.

 
#pragma once
 
#include <array>
#include <set>
#include <vector>
 
#include "autopas/cells/FullParticleCell.h"
#include "autopas/containers/octree/OctreeDirection.h"
#include "autopas/utils/inBox.h"
 
namespace autopas {
template <typename Particle_T>
class OctreeLeafNode;
 
template <class Particle_T>
class OctreeNodeInterface {
 public:
  OctreeNodeInterface(const std::array<double, 3> &boxMin, const std::array<double, 3> &boxMax,
                      OctreeNodeInterface<Particle_T> *parent, const int unsigned treeSplitThreshold,
                      const double interactionLength, const double cellSizeFactor)
      : _boxMin(boxMin),
        _boxMax(boxMax),
        _parent(parent),
        _treeSplitThreshold(treeSplitThreshold),
        _interactionLength(interactionLength),
        _cellSizeFactor(cellSizeFactor) {}
 
  virtual ~OctreeNodeInterface() = default;
 
  OctreeNodeInterface(const OctreeNodeInterface<Particle_T> &) = default;
 
  virtual std::unique_ptr<OctreeNodeInterface<Particle_T>> insert(const Particle_T &p) = 0;
 
  virtual bool deleteParticle(Particle_T &particle) = 0;
 
  virtual void collectAllParticles(std::vector<Particle_T *> &ps) const = 0;
 
  virtual void appendAllLeafBoxes(
      std::vector<std::pair<std::array<double, 3>, std::array<double, 3>>> &boxes) const = 0;
 
  virtual void appendAllLeaves(std::vector<OctreeLeafNode<Particle_T> *> &leaves) const = 0;
 
  virtual void clearChildren(std::unique_ptr<OctreeNodeInterface<Particle_T>> &ref) = 0;
 
  virtual size_t size() const = 0;
 
  virtual size_t getNumberOfParticles(IteratorBehavior behavior = IteratorBehavior::owned) const = 0;
 
  virtual OctreeNodeInterface<Particle_T> *SON(octree::Octant O) = 0;
 
  virtual bool hasChildren() = 0;
 
  virtual OctreeNodeInterface<Particle_T> *getChild(int index) = 0;
 
  virtual std::set<OctreeLeafNode<Particle_T> *> getLeavesInRange(const std::array<double, 3> &min,
                                                                  const std::array<double, 3> &max) = 0;
 
  bool isInside(const std::array<double, 3> &point) {
    using namespace autopas::utils;
    return inBox(point, _boxMin, _boxMax);
  }
 
  static bool volumeExistsOnAxis(const int axis, const std::array<double, 3> &aMin, const std::array<double, 3> &aMax,
                                 const std::array<double, 3> &bMin, const std::array<double, 3> &bMax) {
    const bool o1 = aMin[axis] < bMax[axis];
    const bool o2 = bMin[axis] < aMax[axis];
    return o1 and o2;
  }
 
  bool enclosesVolumeWithOtherOnAxis(const int axis, const OctreeNodeInterface<Particle_T> *other) {
    return volumeExistsOnAxis(axis, this->getBoxMin(), this->getBoxMax(), other->getBoxMin(), other->getBoxMax());
  }
 
  bool overlapsBox(const std::array<double, 3> &otherMin, const std::array<double, 3> &otherMax) {
    bool result = true;
    for (auto d = 0; d < 3; ++d) {
      result &= (this->_boxMin[d] <= otherMax[d]) and (this->_boxMax[d] >= otherMin[d]);
    }
    return result;
  }
 
  static double getEnclosedVolumeWith(const std::array<double, 3> &aMin, const std::array<double, 3> &aMax,
                                      const std::array<double, 3> &bMin, const std::array<double, 3> &bMax) {
    auto product = 1.0;
    int count = 0;
    for (auto d = 0; d < 3; ++d) {
      auto minOnAxis = std::max(aMin[d], bMin[d]);
      auto maxOnAxis = std::min(aMax[d], bMax[d]);
      auto dim = (maxOnAxis - minOnAxis);
      if (dim > 0) {
        ++count;
      } else {
        break;
      }
      product *= dim;
    }
    auto result = count == 3 ? product : 0;
    return result;
  }
 
  double getEnclosedVolumeWith(const std::array<double, 3> &otherMin, const std::array<double, 3> &otherMax) {
    return getEnclosedVolumeWith(this->getBoxMin(), this->getBoxMax(), otherMin, otherMax);
  }
 
  OctreeNodeInterface<Particle_T> *EQ_FACE_NEIGHBOR(const octree::Face I) {
    OctreeNodeInterface<Particle_T> *param, *P = this;
    if (ADJ(I, SONTYPE(P))) {
      param = FATHER(P)->EQ_FACE_NEIGHBOR(I);
    } else {
      param = FATHER(P);
    }
    return param->SON(REFLECT(I, SONTYPE(P)));
  }
 
  OctreeNodeInterface<Particle_T> *EQ_EDGE_NEIGHBOR(const octree::Edge I) {
    OctreeNodeInterface<Particle_T> *param, *P = this;
    if (ADJ(I, SONTYPE(P))) {
      param = FATHER(P)->EQ_EDGE_NEIGHBOR(I);
    } else if (COMMON_FACE(I, SONTYPE(P)) != octree::O) {
      param = FATHER(P)->EQ_FACE_NEIGHBOR(COMMON_FACE(I, SONTYPE(P)));
    } else {
      param = FATHER(P);
    }
    return param->SON(REFLECT(I, SONTYPE(P)));
  }
 
  OctreeNodeInterface<Particle_T> *EQ_VERTEX_NEIGHBOR(const octree::Vertex I) {
    OctreeNodeInterface<Particle_T> *param, *P = this;
    if (ADJ(I, SONTYPE(P))) {
      param = FATHER(P)->EQ_VERTEX_NEIGHBOR(I);
    } else if (COMMON_EDGE(I, SONTYPE(P)) != octree::OO) {
      param = FATHER(P)->EQ_EDGE_NEIGHBOR(COMMON_EDGE(I, SONTYPE(P)));
    } else if (COMMON_FACE(I, SONTYPE(P)) != octree::O) {
      param = FATHER(P)->EQ_FACE_NEIGHBOR(COMMON_FACE(I, SONTYPE(P)));
    } else {
      param = FATHER(P);
    }
    return param->SON(REFLECT(I, SONTYPE(P)));
  }
 
  OctreeNodeInterface<Particle_T> *GTEQ_FACE_NEIGHBOR(octree::Face I);
 
  OctreeNodeInterface<Particle_T> *GTEQ_EDGE_NEIGHBOR(octree::Edge I);
 
  OctreeNodeInterface<Particle_T> *GTEQ_VERTEX_NEIGHBOR(octree::Vertex I);
 
  virtual std::vector<OctreeLeafNode<Particle_T> *> getLeavesFromDirections(
      const std::vector<octree::Vertex> &directions) = 0;
 
  std::vector<OctreeLeafNode<Particle_T> *> getNeighborLeaves(const octree::Any direction) {
    auto opposite = octree::getOppositeDirection(direction);
    auto directions = octree::getAllowedDirections(opposite);
    auto neighborLeaves = getLeavesFromDirections(directions);
    return neighborLeaves;
  }
 
  std::set<OctreeLeafNode<Particle_T> *> getNeighborLeaves() {
    std::set<OctreeLeafNode<Particle_T> *> result;
 
    // Get all face neighbors
    for (auto face : octree::Tables::faces) {
      OctreeNodeInterface<Particle_T> *neighbor = GTEQ_FACE_NEIGHBOR(face);
      if (neighbor) {
        auto leaves = neighbor->getNeighborLeaves(face);
        result.insert(leaves.begin(), leaves.end());
      }
    }
 
    // Get all edge neighbors
    for (auto edge : octree::Tables::edges) {
      OctreeNodeInterface<Particle_T> *neighbor = GTEQ_EDGE_NEIGHBOR(edge);
      if (neighbor) {
        auto leaves = neighbor->getNeighborLeaves(edge);
        result.insert(leaves.begin(), leaves.end());
      }
    }
 
    // Get all face neighbors
    for (auto vertex : octree::Tables::vertices) {
      OctreeNodeInterface<Particle_T> *neighbor = GTEQ_VERTEX_NEIGHBOR(vertex);
      if (neighbor) {
        auto leaves = neighbor->getNeighborLeaves(vertex);
        result.insert(leaves.begin(), leaves.end());
      }
    }
 
    return result;
  }
 
  [[nodiscard]] const std::array<double, 3> &getBoxMin() const { return _boxMin; }
 
  [[nodiscard]] const std::array<double, 3> &getBoxMax() const { return _boxMax; }
 
  OctreeNodeInterface<Particle_T> *getParent() const { return _parent; }
 
 protected:
  bool hasParent() { return _parent != nullptr; }
 
  OctreeNodeInterface<Particle_T> *_parent;
 
  std::array<double, 3> _boxMin;
 
  std::array<double, 3> _boxMax;
 
  int unsigned _treeSplitThreshold;
 
  double _interactionLength;
 
  double _cellSizeFactor;
};
 
template <class Particle_T>
inline bool GRAY(OctreeNodeInterface<Particle_T> *node) {
  // According to Samet: "All non-leaf nodes are said to be GRAY"
  return node->hasChildren();
}
 
template <class Particle_T>
inline OctreeNodeInterface<Particle_T> *FATHER(const OctreeNodeInterface<Particle_T> *node) {
  return node->getParent();
}
 
template <class Particle_T>
static octree::Octant SONTYPE(const OctreeNodeInterface<Particle_T> *node) {
  octree::Octant result = octree::OOO;
  if (FATHER(node)) {
    for (octree::Vertex test : octree::Tables::vertices) {
      if (FATHER(node)->SON(test) == node) {
        result = test;
        break;
      }
    }
    if (result == octree::OOO) {
      throw std::runtime_error("[OctreeNodeInterface::SONTYPE()] Unable to determine SONTYPE");
    }
  }
  return result;
}
 
template <class Particle_T>
OctreeNodeInterface<Particle_T> *OctreeNodeInterface<Particle_T>::GTEQ_FACE_NEIGHBOR(const octree::Face I) {
  // Check precondition
  if (not isFace(I)) {
    throw std::runtime_error("[OctreeNodeInterface::GTEQ_FACE_NEIGHBOR()] Received invalid face.");
  }
 
  auto null = [](OctreeNodeInterface<Particle_T> *T) { return T == nullptr; };
 
  // Find a common ancestor
  OctreeNodeInterface<Particle_T> *Q, *P = this;
  if ((not null(FATHER(P))) and ADJ(I, SONTYPE(P))) {
    Q = FATHER(P)->GTEQ_FACE_NEIGHBOR(I);
  } else {
    Q = FATHER(P);
  }
 
  if ((not null(Q)) and GRAY(Q)) {
    // Follow the reflected path to locate the neighbor
    return Q->SON(REFLECT(I, SONTYPE(P)));
  } else {
    return Q;
  }
}
 
template <class Particle_T>
OctreeNodeInterface<Particle_T> *OctreeNodeInterface<Particle_T>::GTEQ_EDGE_NEIGHBOR(const octree::Edge I) {
  // Check precondition
  if (not isEdge(I)) {
    throw std::runtime_error("[OctreeNodeInterface::GTEQ_EDGE_NEIGHBOR()] Received invalid edge.");
  }
 
  auto null = [](OctreeNodeInterface<Particle_T> *T) { return T == nullptr; };
 
  // Find a common ancestor
  OctreeNodeInterface<Particle_T> *Q, *P = this;
  if (null(FATHER(P))) {
    Q = nullptr;
  } else if (ADJ(I, SONTYPE(P))) {
    Q = FATHER(P)->GTEQ_EDGE_NEIGHBOR(I);
  } else if (octree::Face common = COMMON_FACE(I, SONTYPE(P)); common != octree::O) {
    Q = FATHER(P)->GTEQ_FACE_NEIGHBOR(common);
  } else {
    Q = FATHER(P);
  }
 
  if ((not null(Q)) and GRAY(Q)) {
    // Follow opposite path to locate the neighbor
    return Q->SON(REFLECT(I, SONTYPE(P)));
  } else {
    return Q;
  }
}
 
template <class Particle_T>
OctreeNodeInterface<Particle_T> *OctreeNodeInterface<Particle_T>::GTEQ_VERTEX_NEIGHBOR(const octree::Vertex I) {
  // Check precondition
  if (not isVertex(I)) {
    throw std::runtime_error("[OctreeNodeInterface::GTEQ_VERTEX_NEIGHBOR()] Received invalid vertex.");
  }
 
  // Find a common ancestor
  OctreeNodeInterface<Particle_T> *Q, *P = this;
  if (not FATHER(P)) {
    Q = nullptr;
  } else if (ADJ(I, SONTYPE(P))) {
    Q = FATHER(P)->GTEQ_VERTEX_NEIGHBOR(I);
  } else if (octree::Edge commonEdge = COMMON_EDGE(I, SONTYPE(P)); commonEdge != octree::OO) {
    Q = FATHER(P)->GTEQ_EDGE_NEIGHBOR(commonEdge);
  } else if (octree::Face commonFace = COMMON_FACE(I, SONTYPE(P)); commonFace != octree::O) {
    Q = FATHER(P)->GTEQ_FACE_NEIGHBOR(commonFace);
  } else {
    Q = FATHER(P);
  }
 
  if (Q and GRAY(Q)) {
    // Follow opposite path to locate the neighbor
    return Q->SON(REFLECT(I, SONTYPE(P)));
  } else {
    return Q;
  }
}
}  // namespace autopas