Commit 8a217f00 authored by Quentin Khan's avatar Quentin Khan

Remove adaptive development from master

parent 90d89c0d
#ifndef SCALFMM_BOX_HPP_
#define SCALFMM_BOX_HPP_
#include <ostream>
#include "FZCurve.hpp"
/* Implements a N dimensions box
*
* \author Quentin Khan <quentin.khan@inria.fr>
*
* The box is described by two opposite corners : the maximum and
* the minimum one. All the class transformations maintain this
* predicate.
*
* \tparam Real Floating number representation.
* \tparam Dim Space dimension count.
* \tparam SpaceFillingCurve A templatize implementation of a space filling curve
*
*/
template<class _Position, template<std::size_t> class SpaceFillingCurve = FZCurve>
class FBox {
public:
/// Position type
using position_t = _Position;
/// Floating number representation
using Real = typename position_t::Real;
/// Space dimension
constexpr static const std::size_t Dim = position_t::Dim;
/// Space filling curve type
using space_filling_curve_t = SpaceFillingCurve<Dim>;
private:
position_t _c1; ///< Minimum corner
position_t _c2; ///< Maximum corner
position_t _center; ///< Center
/** Rearranges the corners to ensure the maximum-minimum predicate */
void rearrange_corners () {
for(std::size_t i = 0; i < Dim; ++i){
if(_c1[i] > _c2[i]) {
std::swap(_c1[i], _c2[i]);
}
}
_center = (_c1 + _c2) / 2;
}
public:
/// Accessor for the minimum corner
const position_t& c1 = _c1;
/// Accessor for the maximum corner
const position_t& c2 = _c2;
/** Builds an empty box at the origin */
FBox() = default;
/** Copies an existing box */
FBox(const FBox& other) : _c1(other._c1), _c2(other._c2), _center(other._center) {}
/** Copies an other box */
FBox& operator=(const FBox& other) {
this->_c1 = other._c1;
this->_c2 = other._c2;
this->_center = other._center;
return *this;
}
/** Builds a cube from the lower corner and its side length
*
* \param min_corner The lowest corner
* \param side_length The cube's side length
**/
FBox(const position_t& min_corner, Real side_length) :
_c1(min_corner),
_c2(min_corner)
{
if(side_length < 0) {
side_length = - side_length;
}
for(auto&& d : _c2) {
d += side_length;
}
_center = (_c1 + _c2) / 2;
}
/** Builds a box from two corners
*
* The maximum and minimum corners are deduced from the given corners.
*
* \param corner_1 The first corner.
* \param corner_2 The second corner.
*/
FBox(const position_t& corner_1, const position_t& corner_2) : _c1(corner_1), _c2(corner_2) {
rearrange_corners();
}
/** Changes the box corners
*
* The maximum and minimum corners are deduced from the given corners.
*
* \param corner_1 The first corner.
* \param corner_2 The second corner.
*/
void set(const position_t& corner_1, const position_t& corner_2) {
_c1 = corner_1;
_c2 = corner_2;
rearrange_corners();
}
/** Checks whether a position is within the box bounds
*
* \param p The position to check.
*/
bool contains(const position_t& p) const {
for(std::size_t i = 0; i < Dim; i++) {
if(p[i] < _c1[i] || p[i] > _c2[i]) {
return false;
}
}
return true;
}
/** Checks whether an object's position is within the box bounds
*
* The object must implement a `position_t position() const;` method.
*
* \tparam T The object type.
* \param obj The object which position to check.
*/
template<class T>
bool contains(const T& obj) const {
return contains(obj.position());
}
/** Accessor for the box center */
const position_t& center() const {
return _center;
}
/** Accessor for the box corners
*
* The corners are numbered using a space filling curve.
*
* \param idx The corner index.
* \return The idx'th corner.
*/
position_t corner(std::size_t idx) const {
position_t c;
std::size_t i = 0;
for(bool choice : space_filling_curve_t::position(idx)) {
c[i] = choice ? _c2[i] : _c1[i];
++i;
}
return c;
}
/** Setter for the corners
*
* Moves a corner to a new position and modifies the relevant other
* ones. The corners are numbered using a space filling curve.
*
* \param idx The moved corner index.
* \param pos The new corner position.
*/
void corner(std::size_t idx, const position_t& pos) {
std::size_t i = 0;
for(bool choice : space_filling_curve_t::position(idx)) {
if(choice) {
_c2[i] = pos[i];
} else {
_c1[i] = pos[i];
}
++i;
}
rearrange_corners();
}
/** Sums the corners of two boxes */
FBox operator+(const FBox& other) const {
return FBox(_c1 + other._c1, _c2 + other._c2);
}
/** Tests two boxes equality */
bool operator==(const FBox& other) const {
return c1 == other.c1 && c2 == other.c2;
}
/** Tests two boxes inequality */
bool operator!=(const FBox& other) const {
return ! this->operator==(other);
}
friend std::ostream& operator<<(std::ostream& os, const FBox& box) {
return os << "[" << box.c1 << "," << box.c2 << "]";
}
};
#endif
#ifndef _SCALFMM_IN_ORDER_NODE_ITERATOR_HPP_
#define _SCALFMM_IN_ORDER_NODE_ITERATOR_HPP_
#include "FNodeIterator.hpp"
namespace scalfmm {
namespace tests {
struct test_InOrderNodeIterator;
}
}
template<class Node>
class FInOrderNodeIterator : public FNodeIterator<FInOrderNodeIterator<Node>, Node> {
using base_t = FNodeIterator<FInOrderNodeIterator<Node>, Node>;
using Direction = typename base_t::Direction;
friend scalfmm::tests::test_InOrderNodeIterator;
friend base_t;
public:
using base_t::base_t;
using base_t::operator++;
using base_t::operator--;
FInOrderNodeIterator& operator++() {
return move_iterator<Direction::forwards>();
}
FInOrderNodeIterator& operator--() {
return move_iterator<Direction::backwards>();
}
private:
/** Indicates the index of the middle child */
template<Direction dir>
typename base_t::index_t mid_idx() const {
if(dir == Direction::forwards) { // Compile time optimisation will remove unused branch
return (base_t::child_count/2)-1;
} else {
return (base_t::child_count/2);
}
}
/** Moves the cursor to the first node
*
* Resets the iterator to a valid state.
*/
template<Direction dir>
FInOrderNodeIterator& move_to_boundary() {
base_t::template goto_root<dir>();
while(! this->_cur->is_leaf()) {
base_t::template move_down<dir>();
}
return *this;
}
template<Direction dir>
FInOrderNodeIterator& leaf_next() {
// Leave current leaf
if(!this->_cur_state.empty()) {
base_t::move_up();
}
// If we visited half the parent's children, stay on the parent
if (!this->_cur_state.empty() && mid_idx<dir>() == this->_cur_state.back()) {
return *this;
}
// If we visited all the children of the node, we move up
while(!this->_cur_state.empty() && base_t::template at_last_child<dir>()) {
base_t::move_up();
}
// Moving too far up means we got to the end of the tree
if(this->_cur_state.empty()) {
base_t::template set_end<dir>();
return *this;
}
// Once we moved up enough, we move down to the next leaf
while(! this->_cur->is_leaf()) {
base_t::template move_down<dir>();
}
return *this;
}
template<Direction dir>
FInOrderNodeIterator& node_next() {
// If we are on an internal node, we move down to the children
if(dir == Direction::forwards) { // Compile time optimisation will remove unused branch
this->_cur_state.back() = static_cast<char>(mid_idx<dir>()+1);
} else {
this->_cur_state.back() = static_cast<char>(mid_idx<dir>()-1);
}
this->_cur = this->_cur->getChild((std::size_t)this->_cur_state.back());
this->_cur_state.push_back(base_t::template first_idx<dir>());
while(! this->_cur->is_leaf()) {
base_t::template move_down<dir>();
}
return *this;
}
template<Direction dir>
FInOrderNodeIterator& move_iterator() {
// If before the beginning, move to the first node
if(base_t::template is_other_end<dir>()) {
return move_to_boundary<dir>();
}
// Don't move if already past the end
if(base_t::template is_end<dir>()) {
return *this;
}
if(this->_cur->is_leaf()) {
return leaf_next<dir>();
} else {
return node_next<dir>();
}
}
};
#endif
#ifndef SCALFMM_NODE_HPP_
#define SCALFMM_NODE_HPP_
#include <cmath>
#include <functional>
#include <memory>
#include <stdexcept>
#include <string>
#include <vector>
#include <unordered_set>
#include "FBox.hpp"
#include "Utils/FPoint.hpp"
#include "Utils/FConstFuncs.hpp"
namespace scalfmm {
namespace tests {
struct test_Node;
struct test_NodeIterator;
struct test_InOrderNodeIterator;
struct test_PreOrderNodeIterator;
struct test_PostOrderNodeIterator;
}
}
struct NodeEmptyData {};
/** \brief Tree node implementation
*
*/
template<class _Tree, class _ParticleContainer, class NodeData>
class FNode {
public:
/// Tree type this class belongs to
using tree_t = _Tree;
/// Type used to represent real numbers
using Real = typename tree_t::Real;
/// Space dimension count
constexpr static const std::size_t Dim = _Tree::Dim;
/// Child count if the node is not a leaf
constexpr static std::size_t child_count = Fpow(2,Dim);
/// Node position type
using position_t = FPoint<Real, Dim>;
/// Node bounding box type
using box_t = FBox<position_t>;
/// Interaction lists type
using interaction_list_t = std::unordered_set<FNode*>;
/// Children array type
using child_node_array_t = std::array<FNode*, child_count>;
/// Particle container type
using particle_container_t = _ParticleContainer;
/// Particle type
using particle_t = typename particle_container_t::value_type;
/// Node data structure
using data_t = NodeData;
private:
friend struct scalfmm::tests::test_Node;
friend struct scalfmm::tests::test_NodeIterator;
friend struct scalfmm::tests::test_InOrderNodeIterator;
friend struct scalfmm::tests::test_PreOrderNodeIterator;
friend struct scalfmm::tests::test_PostOrderNodeIterator;
friend tree_t;
/// Children array, filled with nullptr is Node is a leaf
child_node_array_t _children{{}};
/// Node parent, nullptr if the node is a root
FNode* _parent = nullptr;
/// Node spatial bounding box
box_t _box;
/// Node depth in its tree
std::size_t _depth = 0;
/// Node index in parent child array
std::size_t _m_idx = 0;
/// Particle container
std::unique_ptr<particle_container_t> _p_container{new particle_container_t};
/// Node data
data_t _data = data_t();
/// Tree the Node belongs to
tree_t* _tree;
/// Indicates whether node is a leaf
bool _is_leaf = true;
public:
/// Near-field leaves interaction list
interaction_list_t U;
/// Mid-field node interaction list
interaction_list_t V;
/// Near-field node interaction list
interaction_list_t W;
/// Mid-field leaf interaction list
interaction_list_t X;
/** Constructor called from parent
*
* \param parent The parent node
* \param child_index The index of this node in the parent children array
*/
FNode(FNode& parent, const std::size_t& child_index) :
_parent(&parent),
_box (parent.getBox().center(), parent.getBox().corner(child_index) ),
_depth (parent.getDepth()+1),
_m_idx((parent.getIndex() << Dim) + child_index),
_tree (parent._tree )
{
if (child_index >= child_count) {
throw std::invalid_argument(std::string("Wrong child index in node contructor: got ")
+ std::to_string(child_index)
+ std::string(", expected at most ")
+ std::to_string(Dim)
+ ".");
}
}
/** Root constructor called from tree */
FNode(tree_t* tree) :
_box(tree->box()),
_tree(tree)
{
tree->leaves().insert(this);
}
/** Default constructor */
FNode() = delete;
/** Copy constructor */
FNode(const FNode& other) = delete;
/** Copy operator */
FNode& operator=(const FNode& other) = delete;
/** Move constructor */
FNode(FNode&& other) = delete;
/** Move operator */
FNode& operator=(FNode&& other) = delete;
/** Destructor */
~FNode() {
for(auto&& child : _children) {
delete child;
}
}
/// Data accessor
data_t& getData() noexcept {
return _data;
}
/// Data const accessor
const data_t& getData() const noexcept {
return _data;
}
/// Children container accessor
child_node_array_t& getChildren() noexcept {
return _children;
}
/// Children const container accessor
const child_node_array_t& getChildren() const noexcept {
return _children;
}
/** Child container accessor
*
* \param index Child index
*/
FNode* getChild(const std::size_t& index) noexcept {
return getChildren()[index];
}
/** Child container const accessor
*
* \param index Child index
*/
const FNode* getChild(const std::size_t& index) const noexcept {
return getChildren()[index];
}
/// Parent accessor
FNode* getParent() noexcept {
return _parent;
}
/// Parent const accessor
const FNode* getParent() const noexcept {
return _parent;
}
/// Depth accessor
const std::size_t& getDepth() const noexcept {
return _depth;
}
/// Morton index accessor
const std::size_t& getIndex() const noexcept {
return _m_idx;
}
/// Tree accessor
tree_t& getTree() noexcept {
return *_tree;
}
/// Tree const accessor
const tree_t& getTree() const noexcept {
return *_tree;
}
/// Box const accessor
const box_t& getBox() const noexcept {
return _box;
}
/// Particle container accessor
particle_container_t* getParticleContainer() noexcept {
return _p_container.get();
}
/// Particle container accessor
const particle_container_t* getParticleContainer() const noexcept {
return _p_container.get();
}
/// Particle count for the container
std::size_t getParticleCount() const noexcept {
if(getParticleContainer()) {
return getParticleContainer()->size();
}
return 0;
}
/** Returns true if this node and the 'other' node are adjacent
*
* The nodes are assumed to belong to the same tree.
*
* To check whether nodes are adjacent, on each axis, the distance between
* the nodes' center is compared to the sum of their half diameter. For at
* least one of the axes, the two must be equal. For the others, the
* distance must be less than or equal to the sum. This ensures that a node
* is not adjacent to one of its descendants.
*
* \param other The node to test adjacency with.
*
* \return true if this FNode and the 'other' FNode are adjacent.
*/
bool is_adjacent(const FNode& other) const noexcept {
// Sum of the half side lengh of the two nodes boxes.
// Boxes are cubes, we only need one side.
Real centers_distance = getBox().center()[0] - getBox().c1[0]
+ other.getBox().center()[0] - other.getBox().c1[0];
// Used to check that the other box isn't overlapping with this box
bool one_axis_is_at_exact_distance = false;
position_t my_center = getBox().center();
position_t other_center = other.getBox().center();
for(std::size_t i = 0; i < Dim; ++i) {
Real distance = fabs(my_center[i] - other_center[i]);
if( Ffeq(distance, centers_distance) ) {
one_axis_is_at_exact_distance = true;
} else if(distance > centers_distance) {
return false;
}
}
return one_axis_is_at_exact_distance;
}
/** Tests whether this node and the 'other' node are adjacent
*
* \return true if the nodes are adjacent.
*/
bool is_adjacent(const FNode* other) const noexcept {
if(nullptr == other) {
return false;
} else if (this == other){
return true;
}
return is_adjacent(*other);
}
/** Tests whether this node is a leaf
*
* \return true if this node is a leaf.
*/
bool is_leaf() const noexcept {
return _is_leaf;
}
/** Inserts a particle in the tree rooted at node
*
* Pushes a particle in the node particle container if it is a leaf. If it
* isn't, the particle is forwarded to the relevant child.
*
* \note The `push_back` method of #particle_container_t is used to insert in a leaf
*/
void insert(const particle_t& p) {
if(! is_leaf()) {
std::size_t child_index = box_t::space_filling_curve_t::index(p.position(), getBox().center());
getChild(child_index)->insert(p);
} else {
getParticleContainer()->push_back(p);
if(getParticleContainer()->size() > getTree().leaf_max_particle_count()) {
split();
}
}
}