};
/**
- * A class that represents a bounding box in a space with arbitrary dimension
- * <tt>spacedim</tt>.
+ * A class that represents a box of arbitrary dimension <tt>spacedim</tt> and
+ * with sides parallel to the coordinate axes. That is, a region
*
- * Objects of this class are used to represent bounding boxes. They are,
- * among other uses, useful in parallel distributed meshes to give a general
- * description of the owners of each portion of the mesh.
+ * @f[
+ * [x_0^L, x_0^U] \times ... \times [x_{spacedim-1}^L, x_{spacedim-1}^U],
+ * @f]
*
- * Bounding boxes are represented by two vertices (bottom left and top right).
- * Geometrically, a bounding box is:
- * - 1 d: a segment (represented by its vertices in the proper order)
- * - 2 d: a rectangle (represented by the vertices V at bottom left, top right)
+ * where $(x_0^L , ..., x_{spacedim-1}^L) and $(x_0^U , ..., x_{spacedim-1}^U)
+ * denote the two vertices (bottom left and top right) which are used to
+ * represent the box.
+ *
+ * Geometrically, a bounding box is thus:
+ * - 1D: a segment (represented by its vertices in the proper order)
+ * - 2D: a rectangle (represented by the vertices V at bottom left, top right)
* @code
* .--------V
* | |
* V--------.
* @endcode
*
- * - 3 d: a cuboid (in which case the two vertices V follow the convention and
+ * - 3D: a cuboid (in which case the two vertices V follow the convention and
* are not owned by the same face)
* @code
* .------V
* | |/
* V------.
* @endcode
- * Notice that the sides are always parallel to the respective axis.
*
- * @author Giovanni Alzetta, 2017.
+ * Bounding boxes are, for example, useful in parallel distributed meshes to
+ * give a general description of the owners of each portion of the mesh.
+ *
+ * Taking the cross section of a BoundingBox<spacedim> orthogonal to a given
+ * direction gives a box in one dimension lower: BoundingBox<spacedim - 1>.
+ * In 3D, the 2 coordinates of the cross section of BoundingBox<3> can be
+ * ordered in 2 different ways. That is, if we take the cross section orthogonal
+ * to the y direction we could either order a 3D-coordinate into a
+ * 2D-coordinate as $(x,z)$ or as $(z,x)$. This class uses the second
+ * convention, corresponding to the coordinates being ordered cyclicly
+ * $x \rightarrow y \rightarrow z \rightarrow x \rightarrow ... $
+ * To be precise, if we take a cross section:
+ *
+ * | Orthogonal to | Cross section coordinates ordered as |
+ * |:-------------:|:------------------------------------:|
+ * | x | (y, z) |
+ * | y | (z, x) |
+ * | z | (x, y) |
+ *
+ * This is according to the convention set by the function
+ * <code>coordinate_to_one_dim_higher</code>.
+ *
+ * @author Giovanni Alzetta, 2017, Simon Sticko 2020.
*/
template <int spacedim, typename Number = double>
class BoundingBox
double
volume() const;
+ /**
+ * Returns the point in the center of the box.
+ */
+ Point<spacedim, Number>
+ center() const;
+
+ /**
+ * Returns the side length of the box in @p direction.
+ */
+ Number
+ side_length(const unsigned int direction) const;
+
+ /**
+ * Return the lower bound of the box in @p direction.
+ */
+ Number
+ lower_bound(const unsigned int direction) const;
+
+ /**
+ * Return the upper bound of the box in @p direction.
+ */
+ Number
+ upper_bound(const unsigned int direction) const;
+
+ /**
+ * Return the bounds of the box in @p direction, as a one-dimensional box.
+ */
+ BoundingBox<1, Number>
+ bounds(const unsigned int direction) const;
+
+ /**
+ * Returns the indexth vertex of the box.
+ */
+ Point<spacedim, Number>
+ vertex(const unsigned int index) const;
+
+ /**
+ * Returns the indexth child of the box. Child is meant in the same way as for
+ * a cell.
+ */
+ BoundingBox<spacedim, Number>
+ child(const unsigned int index) const;
+
+ /**
+ * Returns the cross section of the box orthogonal to @p direction.
+ * This is a box in one dimension lower.
+ */
+ BoundingBox<spacedim - 1, Number>
+ cross_section(const unsigned int direction) const;
+
/**
* Boost serialization function
*/
std::pair<Point<spacedim, Number>, Point<spacedim, Number>> boundary_points;
};
+
+/**
+ * Returns the unit box [0,1]^dim.
+ */
+template <int dim, typename Number = double>
+BoundingBox<dim, Number>
+create_unit_box();
+
+
+/**
+ * This function defines a convention for how coordinates in dim
+ * dimensions should translate to the coordinates in dim + 1 dimensions,
+ * when one of the coordinates in dim + 1 dimensions is locked to a given
+ * value.
+ *
+ * The convention is the following: Starting from the locked coordinate we
+ * store the lower dimensional coordinates consecutively and wrapping
+ * around when going over the dimension. This relationship can be
+ * described by the following tables:
+ *
+ * 2D
+ * ------------------------------------------------
+ * | locked in 2D | 1D coordinate | 2D coordinate |
+ * |--------------------------------------------- |
+ * | x0 | (a) | (x0, a) |
+ * | x1 | (a) | (a , x1) |
+ * ------------------------------------------------
+ *
+ * 3D
+ * --------------------------------------------------
+ * | locked in 3D | 2D coordinates | 3D coordinates |
+ * |----------------------------------------------- |
+ * | x0 | (a, b) | (x0, a, b) |
+ * | x1 | (a, b) | ( b, x1, a) |
+ * | x2 | (a, b) | ( a, b, x2) |
+ * --------------------------------------------------
+ *
+ * Given a locked coordinate, this function maps a coordinate index in dim
+ * dimension to a coordinate index in dim + 1 dimensions.
+ *
+ * @param locked_coordinate should be in the range [0, dim+1).
+ * @param coordiante_in_dim should be in the range [0, dim).
+ * @return A coordinate index in the range [0, dim+1)
+ */
+template <int dim>
+inline unsigned int
+coordinate_to_one_dim_higher(const unsigned int locked_coordinate,
+ const unsigned int coordiante_in_dim)
+{
+ AssertIndexRange(locked_coordinate, dim + 1);
+ AssertIndexRange(coordiante_in_dim, dim);
+ return (locked_coordinate + coordiante_in_dim + 1) % (dim + 1);
+}
+
+
/*------------------------ Inline functions: BoundingBox --------------------*/
#ifndef DOXYGEN
//
// ---------------------------------------------------------------------
#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/geometry_info.h>
DEAL_II_NAMESPACE_OPEN
return vol;
}
+
+
+template <int spacedim, typename Number>
+Number
+BoundingBox<spacedim, Number>::lower_bound(const unsigned int direction) const
+{
+ AssertIndexRange(direction, spacedim);
+
+ return boundary_points.first[direction];
+}
+
+
+
+template <int spacedim, typename Number>
+Number
+BoundingBox<spacedim, Number>::upper_bound(const unsigned int direction) const
+{
+ AssertIndexRange(direction, spacedim);
+
+ return boundary_points.second[direction];
+}
+
+
+
+template <int spacedim, typename Number>
+Point<spacedim, Number>
+BoundingBox<spacedim, Number>::center() const
+{
+ Point<spacedim, Number> point;
+ for (unsigned int i = 0; i < spacedim; ++i)
+ point[i] = .5 * (boundary_points.first[i] + boundary_points.second[i]);
+
+ return point;
+}
+
+
+
+template <int spacedim, typename Number>
+BoundingBox<1, Number>
+BoundingBox<spacedim, Number>::bounds(const unsigned int direction) const
+{
+ AssertIndexRange(direction, spacedim);
+
+ std::pair<Point<1, Number>, Point<1, Number>> lower_upper_bounds;
+ lower_upper_bounds.first[0] = lower_bound(direction);
+ lower_upper_bounds.second[0] = upper_bound(direction);
+
+ return BoundingBox<1, Number>(lower_upper_bounds);
+}
+
+
+
+template <int spacedim, typename Number>
+Number
+BoundingBox<spacedim, Number>::side_length(const unsigned int direction) const
+{
+ AssertIndexRange(direction, spacedim);
+
+ return boundary_points.second[direction] - boundary_points.first[direction];
+}
+
+
+
+template <int spacedim, typename Number>
+Point<spacedim, Number>
+BoundingBox<spacedim, Number>::vertex(const unsigned int index) const
+{
+ AssertIndexRange(index, GeometryInfo<spacedim>::vertices_per_cell);
+
+ const Point<spacedim> unit_cell_vertex =
+ GeometryInfo<spacedim>::unit_cell_vertex(index);
+
+ Point<spacedim, Number> point;
+ for (unsigned int i = 0; i < spacedim; ++i)
+ point[i] = boundary_points.first[i] + side_length(i) * unit_cell_vertex[i];
+
+ return point;
+}
+
+
+
+template <int spacedim, typename Number>
+BoundingBox<spacedim, Number>
+BoundingBox<spacedim, Number>::child(const unsigned int index) const
+{
+ AssertIndexRange(index, GeometryInfo<spacedim>::max_children_per_cell);
+
+ // Vertex closest to child.
+ const Point<spacedim, Number> parent_vertex = vertex(index);
+ const Point<spacedim, Number> parent_center = center();
+
+ const Point<spacedim> upper_corner_unit_cell =
+ GeometryInfo<spacedim>::unit_cell_vertex(
+ GeometryInfo<spacedim>::vertices_per_cell - 1);
+
+ const Point<spacedim> lower_corner_unit_cell =
+ GeometryInfo<spacedim>::unit_cell_vertex(0);
+
+ std::pair<Point<spacedim, Number>, Point<spacedim, Number>>
+ child_lower_upper_corner;
+ for (unsigned int i = 0; i < spacedim; ++i)
+ {
+ const double child_side_length = side_length(i) / 2;
+
+ const double child_center = (parent_center[i] + parent_vertex[i]) / 2;
+
+ child_lower_upper_corner.first[i] =
+ child_center + child_side_length * (lower_corner_unit_cell[i] - .5);
+ child_lower_upper_corner.second[i] =
+ child_center + child_side_length * (upper_corner_unit_cell[i] - .5);
+ }
+
+ return BoundingBox<spacedim, Number>(child_lower_upper_corner);
+}
+
+
+
+template <int spacedim, typename Number>
+BoundingBox<spacedim - 1, Number>
+BoundingBox<spacedim, Number>::cross_section(const unsigned int direction) const
+{
+ AssertIndexRange(direction, spacedim);
+
+ std::pair<Point<spacedim - 1, Number>, Point<spacedim - 1, Number>>
+ cross_section_lower_upper_corner;
+ for (unsigned int d = 0; d < spacedim - 1; ++d)
+ {
+ const int index_to_write_from =
+ coordinate_to_one_dim_higher<spacedim - 1>(direction, d);
+
+ cross_section_lower_upper_corner.first[d] =
+ boundary_points.first[index_to_write_from];
+
+ cross_section_lower_upper_corner.second[d] =
+ boundary_points.second[index_to_write_from];
+ }
+
+ return BoundingBox<spacedim - 1, Number>(cross_section_lower_upper_corner);
+}
+
+
+
+template <int dim, typename Number>
+BoundingBox<dim, Number>
+create_unit_box()
+{
+ std::pair<Point<dim, Number>, Point<dim, Number>> lower_upper_corner;
+ for (unsigned int i = 0; i < dim; ++i)
+ {
+ lower_upper_corner.second[i] = 1;
+ }
+ return BoundingBox<dim, Number>(lower_upper_corner);
+}
+
+
#include "bounding_box.inst"
DEAL_II_NAMESPACE_CLOSE
for (deal_II_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS)
{
template class BoundingBox<deal_II_dimension, number>;
+
+ template BoundingBox<deal_II_dimension, number> create_unit_box();
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 2018 by the deal.II Authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Test the following functions of the BoundingBox class
+// cross_section
+// center
+// side_length
+// child
+// vertex
+// bounds
+// and the non-member but related function create_unit_box
+// Each function is tested in the function named test_{member_function_name}.
+
+#include <deal.II/base/bounding_box.h>
+#include <deal.II/base/geometry_info.h>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+// Print the bounds of the incoming box to deallog.
+template <int dim>
+void
+print_box(const BoundingBox<dim> &box)
+{
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ if (d > 0)
+ deallog << "x";
+ deallog << "[" << box.lower_bound(d) << ", " << box.upper_bound(d) << "]";
+ }
+ deallog << std::endl;
+}
+
+
+
+// Construct the cross section orthogonal to all axes and print them.
+template <int dim>
+void
+test_cross_section()
+{
+ deallog << "test_cross_section" << std::endl;
+
+ std::pair<Point<dim>, Point<dim>> lower_upper_corners;
+ for (int d = 0; d < dim; ++d)
+ {
+ lower_upper_corners.first[d] = -(d + 1);
+ lower_upper_corners.second[d] = d + 1;
+ }
+
+ const BoundingBox<dim> box(lower_upper_corners);
+
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ deallog << "orthogonal to " << d << std::endl;
+ const BoundingBox<dim - 1> cross_section = box.cross_section(d);
+ print_box(cross_section);
+ }
+}
+
+
+
+// Compute and print the center of the box.
+template <int dim>
+void
+test_center()
+{
+ deallog << "test_center" << std::endl;
+
+ std::pair<Point<dim>, Point<dim>> lower_upper_corners;
+ for (int d = 0; d < dim; ++d)
+ lower_upper_corners.second[d] = 1;
+
+ const BoundingBox<dim> box(lower_upper_corners);
+
+ deallog << "center " << box.center() << std::endl;
+}
+
+
+
+// Print all side lengths of a box.
+template <int dim>
+void
+test_side_length()
+{
+ deallog << "test_side_length" << std::endl;
+
+ std::pair<Point<dim>, Point<dim>> lower_upper_corners;
+ for (int d = 0; d < dim; ++d)
+ lower_upper_corners.second[d] = d + 1;
+
+ const BoundingBox<dim> box(lower_upper_corners);
+
+ for (unsigned int d = 0; d < dim; ++d)
+ deallog << box.side_length(d) << " ";
+ deallog << std::endl;
+}
+
+
+
+// Construct all possible children of a box and print them.
+template <int dim>
+void
+test_child()
+{
+ deallog << "test_child" << std::endl;
+
+ std::pair<Point<dim>, Point<dim>> lower_upper_corners;
+ for (int d = 0; d < dim; ++d)
+ {
+ lower_upper_corners.first[d] = -(d + 1);
+ lower_upper_corners.second[d] = (d + 1);
+ }
+
+ const BoundingBox<dim> box(lower_upper_corners);
+
+ for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_cell; ++i)
+ {
+ deallog << "child " << i << std::endl;
+ const BoundingBox<dim> child = box.child(i);
+ print_box(child);
+ }
+}
+
+
+
+// Print all the vertices of a box.
+template <int dim>
+void
+test_vertex()
+{
+ deallog << "test_vertex" << std::endl;
+
+ std::pair<Point<dim>, Point<dim>> lower_upper_corners;
+ for (int d = 0; d < dim; ++d)
+ {
+ lower_upper_corners.first[d] = -(d + 1);
+ lower_upper_corners.second[d] = d + 1;
+ }
+
+ const BoundingBox<dim> box(lower_upper_corners);
+
+ for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
+ {
+ deallog << "vertex " << i << " = " << box.vertex(i) << std::endl;
+ }
+}
+
+
+
+// Print the results of the bounds function for all coordinate directions.
+template <int dim>
+void
+test_bounds()
+{
+ deallog << "test_bounds" << std::endl;
+
+ std::pair<Point<dim>, Point<dim>> lower_upper_corners;
+ for (int d = 0; d < dim; ++d)
+ {
+ lower_upper_corners.first[d] = -(d + 1);
+ lower_upper_corners.second[d] = d + 1;
+ }
+
+ const BoundingBox<dim> box(lower_upper_corners);
+
+ for (unsigned int i = 0; i < dim; ++i)
+ {
+ const BoundingBox<1> bounds = box.bounds(i);
+ deallog << "Bounds in direction " << i << std::endl;
+ print_box(bounds);
+ }
+}
+
+
+
+// Test that we can call the create_unit_box function.
+template <int dim>
+void
+test_create_unit_box()
+{
+ deallog << "test_create_unit_box" << std::endl;
+ BoundingBox<dim> box = create_unit_box<dim>();
+ print_box(box);
+}
+
+
+
+template <int dim>
+void
+run_test()
+{
+ deallog << "dim = " << dim << std::endl;
+
+ test_cross_section<dim>();
+ deallog << std::endl;
+
+ test_center<dim>();
+ deallog << std::endl;
+
+ test_side_length<dim>();
+ deallog << std::endl;
+
+ test_child<dim>();
+ deallog << std::endl;
+
+ test_vertex<dim>();
+ deallog << std::endl;
+
+ test_bounds<dim>();
+ deallog << std::endl;
+
+ test_create_unit_box<dim>();
+ deallog << std::endl;
+
+ deallog << std::endl;
+}
+
+
+
+int
+main()
+{
+ initlog();
+
+ run_test<1>();
+ run_test<2>();
+ run_test<3>();
+}
--- /dev/null
+
+DEAL::dim = 1
+DEAL::test_cross_section
+DEAL::orthogonal to 0
+DEAL::
+DEAL::
+DEAL::test_center
+DEAL::center 0.500000
+DEAL::
+DEAL::test_side_length
+DEAL::1.00000
+DEAL::
+DEAL::test_child
+DEAL::child 0
+DEAL::[-1.00000, 0.00000]
+DEAL::child 1
+DEAL::[0.00000, 1.00000]
+DEAL::
+DEAL::test_vertex
+DEAL::vertex 0 = -1.00000
+DEAL::vertex 1 = 1.00000
+DEAL::
+DEAL::test_bounds
+DEAL::Bounds in direction 0
+DEAL::[-1.00000, 1.00000]
+DEAL::
+DEAL::test_create_unit_box
+DEAL::[0.00000, 1.00000]
+DEAL::
+DEAL::
+DEAL::dim = 2
+DEAL::test_cross_section
+DEAL::orthogonal to 0
+DEAL::[-2.00000, 2.00000]
+DEAL::orthogonal to 1
+DEAL::[-1.00000, 1.00000]
+DEAL::
+DEAL::test_center
+DEAL::center 0.500000 0.500000
+DEAL::
+DEAL::test_side_length
+DEAL::1.00000 2.00000
+DEAL::
+DEAL::test_child
+DEAL::child 0
+DEAL::[-1.00000, 0.00000]x[-2.00000, 0.00000]
+DEAL::child 1
+DEAL::[0.00000, 1.00000]x[-2.00000, 0.00000]
+DEAL::child 2
+DEAL::[-1.00000, 0.00000]x[0.00000, 2.00000]
+DEAL::child 3
+DEAL::[0.00000, 1.00000]x[0.00000, 2.00000]
+DEAL::
+DEAL::test_vertex
+DEAL::vertex 0 = -1.00000 -2.00000
+DEAL::vertex 1 = 1.00000 -2.00000
+DEAL::vertex 2 = -1.00000 2.00000
+DEAL::vertex 3 = 1.00000 2.00000
+DEAL::
+DEAL::test_bounds
+DEAL::Bounds in direction 0
+DEAL::[-1.00000, 1.00000]
+DEAL::Bounds in direction 1
+DEAL::[-2.00000, 2.00000]
+DEAL::
+DEAL::test_create_unit_box
+DEAL::[0.00000, 1.00000]x[0.00000, 1.00000]
+DEAL::
+DEAL::
+DEAL::dim = 3
+DEAL::test_cross_section
+DEAL::orthogonal to 0
+DEAL::[-2.00000, 2.00000]x[-3.00000, 3.00000]
+DEAL::orthogonal to 1
+DEAL::[-3.00000, 3.00000]x[-1.00000, 1.00000]
+DEAL::orthogonal to 2
+DEAL::[-1.00000, 1.00000]x[-2.00000, 2.00000]
+DEAL::
+DEAL::test_center
+DEAL::center 0.500000 0.500000 0.500000
+DEAL::
+DEAL::test_side_length
+DEAL::1.00000 2.00000 3.00000
+DEAL::
+DEAL::test_child
+DEAL::child 0
+DEAL::[-1.00000, 0.00000]x[-2.00000, 0.00000]x[-3.00000, 0.00000]
+DEAL::child 1
+DEAL::[0.00000, 1.00000]x[-2.00000, 0.00000]x[-3.00000, 0.00000]
+DEAL::child 2
+DEAL::[-1.00000, 0.00000]x[0.00000, 2.00000]x[-3.00000, 0.00000]
+DEAL::child 3
+DEAL::[0.00000, 1.00000]x[0.00000, 2.00000]x[-3.00000, 0.00000]
+DEAL::child 4
+DEAL::[-1.00000, 0.00000]x[-2.00000, 0.00000]x[0.00000, 3.00000]
+DEAL::child 5
+DEAL::[0.00000, 1.00000]x[-2.00000, 0.00000]x[0.00000, 3.00000]
+DEAL::child 6
+DEAL::[-1.00000, 0.00000]x[0.00000, 2.00000]x[0.00000, 3.00000]
+DEAL::child 7
+DEAL::[0.00000, 1.00000]x[0.00000, 2.00000]x[0.00000, 3.00000]
+DEAL::
+DEAL::test_vertex
+DEAL::vertex 0 = -1.00000 -2.00000 -3.00000
+DEAL::vertex 1 = 1.00000 -2.00000 -3.00000
+DEAL::vertex 2 = -1.00000 2.00000 -3.00000
+DEAL::vertex 3 = 1.00000 2.00000 -3.00000
+DEAL::vertex 4 = -1.00000 -2.00000 3.00000
+DEAL::vertex 5 = 1.00000 -2.00000 3.00000
+DEAL::vertex 6 = -1.00000 2.00000 3.00000
+DEAL::vertex 7 = 1.00000 2.00000 3.00000
+DEAL::
+DEAL::test_bounds
+DEAL::Bounds in direction 0
+DEAL::[-1.00000, 1.00000]
+DEAL::Bounds in direction 1
+DEAL::[-2.00000, 2.00000]
+DEAL::Bounds in direction 2
+DEAL::[-3.00000, 3.00000]
+DEAL::
+DEAL::test_create_unit_box
+DEAL::[0.00000, 1.00000]x[0.00000, 1.00000]x[0.00000, 1.00000]
+DEAL::
+DEAL::
//
// ---------------------------------------------------------------------
-#include <deal.II/non_matching/quadrature_generator.h>
+#include <deal.II/base/bounding_box.h>
#include "../tests.h"
deallog << "locked coordinate = " << locked << std::endl;
for (unsigned int i = 0; i < dim; ++i)
{
- deallog << i << " -> "
- << NonMatching::internal::QuadratureGeneratorImplementation::
- coordinate_to_one_dim_higher<dim>(locked, i)
+ deallog << i << " -> " << coordinate_to_one_dim_higher<dim>(locked, i)
<< std::endl;
}
}
projects / dealii.git / commitdiff