From: rschulz Date: Wed, 10 May 2006 14:17:29 +0000 (+0000) Subject: added functions to handle all kinds of vertex-finding related stuff + testcases... X-Git-Tag: v8.0.0~11749 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4b03e8f02a877dc410dd2e3e0bde5398aa9e89ea;p=dealii.git added functions to handle all kinds of vertex-finding related stuff + testcases... git-svn-id: https://svn.dealii.org/trunk@13080 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/geometry_info.h b/deal.II/base/include/base/geometry_info.h index ee54ef4152..aef09b7d06 100644 --- a/deal.II/base/include/base/geometry_info.h +++ b/deal.II/base/include/base/geometry_info.h @@ -656,8 +656,21 @@ struct GeometryInfo * \end{verbatim} */ static const unsigned int dx_to_deal[vertices_per_cell]; - - /** + + + /** + * This field stores for each vertex + * to which faces it belongs. In any + * given dimension, the number of + * faces is equal to the dimension. + * The first index in this 2D-array + * runs over all vertices, the second + * index over @p dim faces to which + * the vertex belongs + */ + static const unsigned int vertex_to_face[vertices_per_cell][dim]; + + /** * This field stores which child * cells are adjacent to a * certain face of the mother @@ -855,6 +868,37 @@ struct GeometryInfo static bool is_inside_unit_cell (const Point &p); /** + * Return true if the given point + * is inside the unit cell of the + * present space dimension. This + * function accepts an additional + * parameter which specifies how + * much the point position may + * actually deviate from a true + * unit cell (may be less than zero) + */ + static bool is_inside_unit_cell (const Point &p, + const double eps); + + /** + * Projects a given point onto the + * unit cell, i.e. each coordinate + * outside [0..1] is modified + * to lie within that interval + */ + static void project_to_unit_cell (Point &p); + + /** + * Returns the infinity norm of + * the vector between a given point @p p + * outside the unit cell to the closest + * unit cell boundary. + * For points inside the cell, this is + * defined as zero. + */ + static double distance_to_unit_cell(Point &p); + + /** * For each face of the reference * cell, this field stores the * coordinate direction in which @@ -1131,6 +1175,39 @@ GeometryInfo<3>::is_inside_unit_cell (const Point<3> &p) (p[2] >= 0.) && (p[2] <= 1.); } +template <> +inline +bool +GeometryInfo<1>::is_inside_unit_cell (const Point<1> &p, const double eps) +{ + return (p[0] >= -eps) && (p[0] <= 1.+eps); +} + + + +template <> +inline +bool +GeometryInfo<2>::is_inside_unit_cell (const Point<2> &p, const double eps) +{ + const double l = -eps, u = 1+eps; + return (p[0] >= l) && (p[0] <= u) && + (p[1] >= l) && (p[1] <= u); +} + + + +template <> +inline +bool +GeometryInfo<3>::is_inside_unit_cell (const Point<3> &p, const double eps) +{ + const double l = -eps, u = 1.0+eps; + return (p[0] >= l) && (p[0] <= u) && + (p[1] >= l) && (p[1] <= u) && + (p[2] >= l) && (p[2] <= u); +} + #endif // DOXYGEN #endif diff --git a/deal.II/base/source/geometry_info.cc b/deal.II/base/source/geometry_info.cc index d261495211..647e0426b7 100644 --- a/deal.II/base/source/geometry_info.cc +++ b/deal.II/base/source/geometry_info.cc @@ -178,6 +178,51 @@ const unsigned int GeometryInfo<4>::dx_to_deal[GeometryInfo<4>::vertices_per_cel invalid_unsigned_int, invalid_unsigned_int}; +template <> +const unsigned int GeometryInfo<1>::vertex_to_face + [GeometryInfo<1>::vertices_per_cell][1] += { { 0 }, + { 1 } }; + +template <> +const unsigned int GeometryInfo<2>::vertex_to_face + [GeometryInfo<2>::vertices_per_cell][2] += { { 0, 2 }, + { 1, 2 }, + { 0, 3 }, + { 1, 3 } }; + +template <> +const unsigned int GeometryInfo<3>::vertex_to_face + [GeometryInfo<3>::vertices_per_cell][3] += { { 0, 2, 4 }, + { 1, 2, 4 }, + { 0, 3, 4 }, + { 1, 3, 4 }, + { 0, 2, 5 }, + { 1, 2, 5 }, + { 0, 3, 5 }, + { 1, 3, 5 } }; + +template <> +const unsigned int GeometryInfo<4>::vertex_to_face + [GeometryInfo<4>::vertices_per_cell][4] += { { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }, + { invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int, invalid_unsigned_int }}; template <> unsigned int @@ -381,8 +426,29 @@ GeometryInfo::face_to_cell_vertices (const unsigned int face, return child_cell_on_face(face, vertex, face_orientation); } +template +void +GeometryInfo::project_to_unit_cell (Point &p) +{ + for(unsigned i=0; i 1.) p[i] = 1.; +} - +template +double +GeometryInfo::distance_to_unit_cell (Point &p) +{ + double result = 0.0; + + for(unsigned i=0; i result) + result = -p[i]; + else if ((p[i]-1.) > result) + result = (p[i] - 1.); + + return result; +} template class GeometryInfo<1>; template class GeometryInfo<2>; diff --git a/deal.II/deal.II/include/grid/grid_tools.h b/deal.II/deal.II/include/grid/grid_tools.h index deed8080be..d3a18b2c88 100644 --- a/deal.II/deal.II/include/grid/grid_tools.h +++ b/deal.II/deal.II/include/grid/grid_tools.h @@ -15,6 +15,7 @@ #include +#include #include #include #include @@ -193,10 +194,52 @@ class GridTools void scale (const double scaling_factor, Triangulation &triangulation); + /** + * Find and return the number of + * the used vertex in a given + * Container that is located closest + * to a given point @p p. The + * type of the first parameter + * may be either Triangulation, + * DoFHandler, hp::DoFHandler, or + * MGDoFHandler. + * + * @author Ralf B. Schulz, 2006 + */ + template class Container> + static + unsigned int + find_closest_vertex (const Container &container, + const Point &p); + + /** + * Find and return a vector of + * iterators to active cells that + * surround a given vertex @p vertex. + * The type of the first parameter + * may be either Triangulation, + * DoFHandler, hp::DoFHandler, or + * MGDoFHandler. + * + * For locally refined grids, the + * vertex itself might not be a vertex + * of all adjacent cells, but will + * always be located on a face or an + * edge of the adjacent cells returned. + * + * @author Ralf B. Schulz, + * Wolfgang Bangert, 2006 + */ + template class Container> + static + std::vector::active_cell_iterator> + find_cells_adjacent_to_vertex(const Container &container, + const unsigned int vertex); + /** * Find and return an iterator to * the active cell that surrounds - * a given point @p ref. The + * a given point @p p. The * type of the first parameter * may be either * Triangulation, @@ -283,6 +326,69 @@ class GridTools find_active_cell_around_point (const Container &container, const Point &p); + /** + * Find and return an iterator to + * the active cell that surrounds + * a given point @p p. The + * type of the first parameter + * may be either + * Triangulation, + * DoFHandler, hp::DoFHandler, or + * MGDoFHandler, i.e. we + * can find the cell around a + * point for iterators into each + * of these classes. + * + * This function works with + * arbitrary boundary mappings, + * using a different algorithm than + * the version of this function above. + * The algorithm used in this + * function proceeds by first + * looking for the vertex that is + * closest to the given point, + * using find_closest_vertex(). + * Then, only in adjacent cells + * to this vertex it is checked + * whether or not the point is + * inside a given cell. + * + * The function returns an iterator + * to the cell, as well as the local + * position of the point inside + * the unit cell. This local position + * might be located slightly outside + * an actual unit cell. + * + * If a point lies on the + * boundary of two or more cells, + * then the algorithm returns + * the cell (A) in which the local + * coordinate is exactly within the + * unit cell (however, for most + * cases, on the boundary the unit + * cell position will be located + * slightly outside the unit cell) + * or (B) the cell of highest + * refinement level; and if there + * are several cells of the same + * refinement level, then it returns + * (C) the one with the lowest distance + * to the actual unit cell. + * + * However, if you are trying + * to locate a vertex, and if the vertex + * can be matched exactly to a cell, + * it is not guaranteed that the most + * refined cell will be returned. + */ + template class Container> + static + std::pair::active_cell_iterator, Point > + find_active_cell_around_point (const Mapping &mapping, + const Container &container, + const Point &p); + /** * Use the METIS partitioner to generate * a partitioning of the active cells @@ -494,6 +600,11 @@ class GridTools << "The point <" << arg1 << "> could not be found inside any of the " << "subcells of a coarse grid cell."); + + DeclException1 (ExcVertexNotUsed, + unsigned int, + << "The given vertex " << arg1 + << " is not used in the given triangulation"); }; diff --git a/deal.II/deal.II/source/grid/grid_tools.cc b/deal.II/deal.II/source/grid/grid_tools.cc index 0241820e18..f7ceeeb8bf 100644 --- a/deal.II/deal.II/source/grid/grid_tools.cc +++ b/deal.II/deal.II/source/grid/grid_tools.cc @@ -29,6 +29,24 @@ #include +// This anonymous namespace contains utility functions to extract the +// triangulation from any container such as DoFHandler, MGDoFHandler, +// and the like +namespace +{ + template + const Triangulation &get_tria(const Triangulation &tria) + { + return tria; + } + + template class Container> + const Triangulation &get_tria(const Container &container) + { + return container.get_tria(); + } +}; + #if deal_II_dimension != 1 template @@ -400,6 +418,122 @@ GridTools::scale (const double scaling_factor, } +template class Container> +unsigned int +GridTools::find_closest_vertex (const Container &container, + const Point &p) +{ + // first get the underlying + // triangulation from the + // container and determine vertices + // and used vertices + const Triangulation &tria = get_tria(container); + + const std::vector< Point > &vertices = tria.get_vertices(); + const std::vector< bool > &used = tria.get_used_vertices(); + + // At the beginning, the first + // used vertex is the closest one + std::vector::const_iterator first = + std::find(used.begin(), used.end(), true); + + // Assert that at least one vertex + // is actually used + Assert(first != used.end(), ExcInternalError()); + + unsigned int best_vertex = std::distance(used.begin(), first); + double best_dist = (p - vertices[best_vertex]).square(); + + // For all remaining vertices, test + // whether they are any closer + for(unsigned int j = best_vertex+1; j < vertices.size(); j++) + if(used[j]) + { + double dist = (p - vertices[j]).square(); + if(dist < best_dist) + { + best_vertex = j; + best_dist = dist; + } + } + + return best_vertex; +} + + +template class Container> +std::vector::active_cell_iterator> +GridTools::find_cells_adjacent_to_vertex(const Container &container, + const unsigned int vertex) +{ + // make sure that the given vertex is + // an active vertex of the underlying + // triangulation + Assert(vertex < get_tria(container).n_vertices(), + ExcIndexRange(0,get_tria(container).n_vertices(),vertex)); + Assert(get_tria(container).get_used_vertices()[vertex], + ExcVertexNotUsed(vertex)); + + std::vector::active_cell_iterator> adjacent_cells; + + typename Container::active_cell_iterator + cell = container.begin_active(), + endc = container.end(); + + // go through all active cells and look + // if the vertex is part of that cell + for(; cell != endc; ++cell) + for(unsigned v = 0; v < GeometryInfo::vertices_per_cell; v++) + if(cell->vertex_index(v) == (int)vertex) + { + // OK, we found a cell that contains + // the particular vertex. We add it + // to the list. + adjacent_cells.push_back(cell); + + // Now we need to make sure that the + // vertex is not a locally refined + // vertex not being part of the + // neighboring cells. So we loop over + // all faces to which this vortex + // belongs, check the level of + // the neighbor, and if it is coarser, + // then check whether the vortex is + // part of that neighbor or not. + for(unsigned vface = 0; vface < dim; vface++) + { + const unsigned face = + GeometryInfo::vertex_to_face[v][vface]; + if(!cell->at_boundary(face)) { + typename Container::cell_iterator + nb = cell->neighbor(face); + + // Here we check whether the neighbor + // is coarser. If it is, we search for + // the vertex in this coarser cell and + // only if not found we will add the + // coarser cell itself + if(nb->level() < cell->level()) { + bool found = false; + for(unsigned v=0; v::vertices_per_cell; v++) + if(cell->vertex_index(v) == (int)vertex) { + found = true; + break; + } + if(!found) + adjacent_cells.push_back(nb); + } + } + } + + break; + } + + Assert(adjacent_cells.size() > 0, ExcInternalError()); + + return adjacent_cells; +} + template typename Container::active_cell_iterator @@ -448,6 +582,66 @@ GridTools::find_active_cell_around_point (const Container &container, } +template class Container> +std::pair::active_cell_iterator, Point > +GridTools::find_active_cell_around_point (const Mapping &mapping, + const Container &container, + const Point &p) +{ + typedef typename Container::active_cell_iterator cell_iterator; + + double best_distance = -1.; + std::pair > best_cell; + + // Find closest vertex and determine + // all adjacent cells + unsigned int vertex = find_closest_vertex(container, p); + + std::vector adjacent_cells = + find_cells_adjacent_to_vertex(container, vertex); + + typename std::vector::const_iterator + cell = adjacent_cells.begin(), + endc = adjacent_cells.end(); + + for(; cell != endc; ++cell) + { + Point p_cell = mapping.transform_real_to_unit_cell(*cell, p); + + // calculate the infinity norm of + // the distance vector to the unit cell. + double dist = GeometryInfo::distance_to_unit_cell(p_cell); + + // If the point is strictly inside the + // unit cell, the cell can be directly + // returned. + if(dist == 0.0) + return std::make_pair(*cell, p_cell); + + // Otherwise, we will compare the new + // cell with a previously found one; + // If there is not a previous one, take + // the current one directly; if there is, + // only replace it if we have found a more + // refined cell or if we have found a cell + // of the same level, but the distance to + // the unit cell is smaller + else if( best_cell.first.state() != IteratorState::valid || + best_cell.first->level() < (*cell)->level() || + ( best_cell.first->level() == (*cell)->level() && + best_distance > dist ) ) + { + best_distance = dist; + best_cell = std::make_pair(*cell, p_cell); + } + } + + Assert(best_cell.first.state() == IteratorState::valid, ExcInternalError()); + Assert(best_distance > 0.0, ExcInternalError()); + + return best_cell; +} + template void @@ -723,6 +917,46 @@ template void GridTools::scale (const double, Triangulation &); +template +unsigned int +GridTools::find_closest_vertex (const Triangulation &, + const Point &); + +template +unsigned int +GridTools::find_closest_vertex (const DoFHandler &, + const Point &); + +template +unsigned int +GridTools::find_closest_vertex (const hp::DoFHandler &, + const Point &); + +template +unsigned int +GridTools::find_closest_vertex (const MGDoFHandler &, + const Point &); + +template +std::vector::active_cell_iterator> +GridTools::find_cells_adjacent_to_vertex(const Triangulation &, + const unsigned int); + +template +std::vector::active_cell_iterator> +GridTools::find_cells_adjacent_to_vertex(const DoFHandler &, + const unsigned int); + +template +std::vector::active_cell_iterator> +GridTools::find_cells_adjacent_to_vertex(const hp::DoFHandler &, + const unsigned int); + +template +std::vector::active_cell_iterator> +GridTools::find_cells_adjacent_to_vertex(const MGDoFHandler &, + const unsigned int); + template Triangulation::active_cell_iterator GridTools::find_active_cell_around_point (const Triangulation &, @@ -743,6 +977,30 @@ MGDoFHandler::active_cell_iterator GridTools::find_active_cell_around_point (const MGDoFHandler &, const Point &p); +template +std::pair::active_cell_iterator, Point > +GridTools::find_active_cell_around_point (const Mapping &, + const Triangulation &, + const Point &); + +template +std::pair::active_cell_iterator, Point > +GridTools::find_active_cell_around_point (const Mapping &, + const DoFHandler &, + const Point &); + +template +std::pair::active_cell_iterator, Point > +GridTools::find_active_cell_around_point (const Mapping &, + const MGDoFHandler &, + const Point &); + +template +std::pair::active_cell_iterator, Point > +GridTools::find_active_cell_around_point (const Mapping &, + const hp::DoFHandler &, + const Point &); + template void GridTools::partition_triangulation (const unsigned int, diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 847fbf5657..32346a41a5 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -290,6 +290,17 @@ inconvenience this causes.
    +
  1. New: GeometryInfo offers several new functions, + is_inside_unit_cell with an epsilon parameter to specify allowable + offsets from the actual unit cell, distance_to_unit_cell returning the + infinity norm of the distance of a given point to the unit cell, and + project_to_unit_cell returning the projection of a point onto the unit + cell. Also, a new member vertex_to_face allow to determine to which + faces of a cell a vertex belongs. +
    + (Ralf B. Schulz 2006/05/10) +

    +
  2. Improved: DataOutBase::OutputFormat has a new value none, writing no output at all. This way, the writing of output files can be controlled more @@ -635,6 +646,21 @@ inconvenience this causes.

    deal.II

      +
    1. + New: In GridTools, several functions have been added. + GridTools::find_closest_vertex searches for the vertex + located at closest distance to a given point. + GridTools::find_cells_adjacent_to_vertex allows to determine + all cells adjacent to a given vertex. And finally, a new version of + find_active_cell_around_point, which takes additionally + a mapping as parameter, implements a new and faster algorithm to + determine the active cell in which a given point is located. For + points located on boundaries and edges, it is in most cases also able + to give the finest cell. +
      + (Ralf B. Schulz 2006/05/10) +

      +
    2. Changed: The DoFObjectAccessor::get_dof_values and DoFObjectAccessor::set_dof_values were part of the diff --git a/tests/bits/Makefile b/tests/bits/Makefile index 53e7381099..e0522f2222 100644 --- a/tests/bits/Makefile +++ b/tests/bits/Makefile @@ -51,6 +51,8 @@ tests_x = geometry_info_* \ mapping_cartesian_1 \ mapping_q4_3d \ q_points \ + find_closest_vertex_* \ + find_cells_adjacent_to_vertex_* \ find_cell_* \ sparsity_pattern_* \ sparse_matrix_* \ diff --git a/tests/bits/find_cell_alt_1.cc b/tests/bits/find_cell_alt_1.cc new file mode 100644 index 0000000000..6455feccbf --- /dev/null +++ b/tests/bits/find_cell_alt_1.cc @@ -0,0 +1,80 @@ +//---------------------------- find_cell_alt_1.cc ----------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- find_cell_alt_1.cc ------------------------ + + +// same as find_cell_1, but for the alternative algorithm +// take a 2d mesh and check that we can find an arbitrary point's cell +// in it + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + + + +void check (Triangulation<2> &tria) +{ + MappingQ<2> map(3); // Let's take a higher order mapping + + Point<2> p (1./3., 1./2.); + + std::pair::active_cell_iterator, Point<2> > + cell = GridTools::find_active_cell_around_point (map, tria, p); + + deallog << cell.first << std::endl; + for (unsigned int v=0; v::vertices_per_cell; ++v) + deallog << "<" << cell.first->vertex(v) << "> "; + deallog << "[ " << cell.second << "] "; + deallog << std::endl; + + Assert (p.distance (cell.first->center()) < cell.first->diameter()/2, + ExcInternalError()); +} + + +int main () +{ + std::ofstream logfile("find_cell_alt_1/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + { + Triangulation<2> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (2); + check (coarse_grid); + } + + { + Triangulation<2> coarse_grid; + GridGenerator::hyper_ball (coarse_grid); + static const HyperBallBoundary<2> boundary; + coarse_grid.set_boundary (0, boundary); + coarse_grid.refine_global (2); + check (coarse_grid); + } +} + + + diff --git a/tests/bits/find_cell_alt_1/cmp/generic b/tests/bits/find_cell_alt_1/cmp/generic new file mode 100644 index 0000000000..164678470c --- /dev/null +++ b/tests/bits/find_cell_alt_1/cmp/generic @@ -0,0 +1,5 @@ + +DEAL::2.3 +DEAL::<0.250000 0.250000> <0.500000 0.250000> <0.250000 0.500000> <0.500000 0.500000> [ 0.333333 1.00000] +DEAL::2.78 +DEAL::<0.250000 0.518306> <0.198223 0.405600> <0.500000 0.500000> <0.396447 0.396447> [ 0.106539 0.363417] diff --git a/tests/bits/find_cell_alt_2.cc b/tests/bits/find_cell_alt_2.cc new file mode 100644 index 0000000000..44d2c5cee2 --- /dev/null +++ b/tests/bits/find_cell_alt_2.cc @@ -0,0 +1,79 @@ +//---------------------------- find_cell_alt_2.cc ------------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- find_cell_alt_2.cc ------------------------ + + +// same as find_cell_alt_1, but in 3d + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + + + +void check (Triangulation<3> &tria) +{ + MappingQ<3> map(3); + + Point<3> p(1./3.,1./2.,1./5.); + + std::pair::active_cell_iterator, Point<3> > cell + = GridTools::find_active_cell_around_point (map, tria, p); + + deallog << cell.first << std::endl; + for (unsigned int v=0; v::vertices_per_cell; ++v) + deallog << "<" << cell.first->vertex(v) << "> "; + deallog << "[ " << cell.second << "] "; + + deallog << std::endl; + + Assert (p.distance (cell.first->center()) < cell.first->diameter()/2, + ExcInternalError()); +} + + +int main () +{ + std::ofstream logfile("find_cell_alt_2/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (2); + check (coarse_grid); + } + + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_ball (coarse_grid); + static const HyperBallBoundary<3> boundary; + coarse_grid.set_boundary (0, boundary); + coarse_grid.refine_global (2); + check (coarse_grid); + } +} + + + diff --git a/tests/bits/find_cell_alt_2/cmp/generic b/tests/bits/find_cell_alt_2/cmp/generic new file mode 100644 index 0000000000..14e42574c4 --- /dev/null +++ b/tests/bits/find_cell_alt_2/cmp/generic @@ -0,0 +1,5 @@ + +DEAL::2.3 +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> [ 0.333333 1.00000 0.800000] +DEAL::2.435 +DEAL::<0.286905 0.658777 -1.73472e-18> <0.209333 0.437629 -3.46945e-18> <0.562887 0.562887 0.00000> <0.418667 0.418667 0.00000> <0.269808 0.604598 0.269808> <0.203251 0.422066 0.203251> <0.530301 0.530301 0.245590> <0.406502 0.406502 0.197169> [ 0.494208 0.406025 0.866848] diff --git a/tests/bits/find_cell_alt_3.cc b/tests/bits/find_cell_alt_3.cc new file mode 100644 index 0000000000..9a02f11132 --- /dev/null +++ b/tests/bits/find_cell_alt_3.cc @@ -0,0 +1,84 @@ +//---------------------------- find_cell_alt_3.cc ------------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- find_cell_alt_3.cc ------------------------ + + +// like find_cell_2, but with the strange meshes from the mesh_3d_* tests + +#include "../tests.h" +#include "mesh_3d.h" + +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + +void check (Triangulation<3> &tria) +{ + MappingQ<3> map(3); + + Point<3> p(1./3.,1./2.,-1./5.); + + std::pair::active_cell_iterator, Point<3> > cell + = GridTools::find_active_cell_around_point (map, tria, p); + + deallog << cell.first << std::endl; + for (unsigned int v=0; v::vertices_per_cell; ++v) + deallog << "<" << cell.first->vertex(v) << "> "; + deallog << "[ " << cell.second << "] "; + deallog << std::endl; + + Assert (p.distance (cell.first->center()) < cell.first->diameter()/2, + ExcInternalError()); +} + + +int main () +{ + std::ofstream logfile("find_cell_alt_3/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + { + Triangulation<3> coarse_grid; + create_two_cubes (coarse_grid); + coarse_grid.refine_global (1); + check (coarse_grid); + } + + { + Triangulation<3> coarse_grid; + create_L_shape (coarse_grid); + coarse_grid.refine_global (1); + check (coarse_grid); + } + + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_ball (coarse_grid); + coarse_grid.refine_global (1); + check (coarse_grid); + } + +} + + + diff --git a/tests/bits/find_cell_alt_3/cmp/generic b/tests/bits/find_cell_alt_3/cmp/generic new file mode 100644 index 0000000000..cbd5329fc6 --- /dev/null +++ b/tests/bits/find_cell_alt_3/cmp/generic @@ -0,0 +1,7 @@ + +DEAL::1.0 +DEAL::<0.00000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.00000 0.00000 -0.500000> <0.500000 0.00000 -0.500000> <0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 0.500000 -0.500000> <0.500000 0.500000 -0.500000> [ 0.666667 0.400000 1.00000] +DEAL::1.16 +DEAL::<0.00000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.00000 0.00000 -0.500000> <0.500000 0.00000 -0.500000> <0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 0.500000 -0.500000> <0.500000 0.500000 -0.500000> [ 0.666667 0.400000 1.00000] +DEAL::1.50 +DEAL::<0.00000 0.577350 -0.577350> <0.00000 0.394338 -0.394338> <0.577350 0.577350 -0.577350> <0.394338 0.394338 -0.394338> <0.00000 0.577350 0.00000> <0.00000 0.394338 6.93889e-18> <0.577350 0.577350 0.00000> <0.394338 0.394338 0.00000> [ 0.408625 0.668369 0.613047] diff --git a/tests/bits/find_cell_alt_4.cc b/tests/bits/find_cell_alt_4.cc new file mode 100644 index 0000000000..3f5ed64519 --- /dev/null +++ b/tests/bits/find_cell_alt_4.cc @@ -0,0 +1,75 @@ +//---------------------------- find_cell_alt_4.cc ------------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005, 2006 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- find_cell_alt_4.cc ------------------------ + + +// take a 3d mesh and check that we can find an arbitrary point's cell +// in it. + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + +void check (Triangulation<3> &tria) +{ + MappingQ1<3> map; + + Point<3> p (0.75,0,0); + + std::pair::active_cell_iterator, Point<3> > cell + = GridTools::find_active_cell_around_point (map, tria, p); + + deallog << cell.first << std::endl; + for (unsigned int v=0; v::vertices_per_cell; ++v) + deallog << "<" << cell.first->vertex(v) << "> "; + deallog << "[ " << cell.second << "] "; + deallog << std::endl; + + Assert (GeometryInfo<3>::distance_to_unit_cell(cell.second) < 1e-10, + ExcInternalError()); +} + + +int main () +{ + std::ofstream logfile("find_cell_alt_4/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + try + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (3); + check (coarse_grid); + } + catch (const std::exception &exc) + { + // we shouldn't get here... + deallog << "Caught an error..." << std::endl; + deallog << exc.what() << std::endl; + } +} + + + diff --git a/tests/bits/find_cell_alt_4/cmp/generic b/tests/bits/find_cell_alt_4/cmp/generic new file mode 100644 index 0000000000..943ccd8e97 --- /dev/null +++ b/tests/bits/find_cell_alt_4/cmp/generic @@ -0,0 +1,3 @@ + +DEAL::3.72 +DEAL::<0.750000 0.00000 0.00000> <0.875000 0.00000 0.00000> <0.750000 0.125000 0.00000> <0.875000 0.125000 0.00000> [ -3.33067e-16 0.00000 5.55112e-17] diff --git a/tests/bits/find_cell_alt_5.cc b/tests/bits/find_cell_alt_5.cc new file mode 100644 index 0000000000..58734539cf --- /dev/null +++ b/tests/bits/find_cell_alt_5.cc @@ -0,0 +1,74 @@ +//---------------------------- find_cell_alt_5.cc ------------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005, 2006 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- find_cell_alt_5.cc ------------------------ + + +// take a 3d mesh and check that we can find an arbitrary point's cell +// in it. + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + +void check (Triangulation<3> &tria) +{ + MappingQ1<3> map; + Point<3> p (0.75,0.75,0.75); + + std::pair::active_cell_iterator, Point<3> > cell + = GridTools::find_active_cell_around_point (map, tria, p); + + deallog << cell.first << std::endl; + for (unsigned int v=0; v::vertices_per_cell; ++v) + deallog << "<" << cell.first->vertex(v) << "> "; + deallog << "[ " << cell.second << "] "; + deallog << std::endl; + + Assert (GeometryInfo<3>::distance_to_unit_cell(cell.second) < 1e-10, + ExcInternalError()); +} + + +int main () +{ + std::ofstream logfile("find_cell_alt_5/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + try + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (3); + check (coarse_grid); + } + catch (const std::exception &exc) + { + // we shouldn't get here... + deallog << "Caught an error..." << std::endl; + deallog << exc.what() << std::endl; + } +} + + + diff --git a/tests/bits/find_cell_alt_5/cmp/generic b/tests/bits/find_cell_alt_5/cmp/generic new file mode 100644 index 0000000000..d16f6f931b --- /dev/null +++ b/tests/bits/find_cell_alt_5/cmp/generic @@ -0,0 +1,3 @@ + +DEAL::3.455 +DEAL::<0.625000 0.625000 0.625000> <0.750000 0.625000 0.625000> <0.625000 0.750000 0.625000> <0.750000 0.750000 0.625000> [ 1.00000 1.00000 1.000000] diff --git a/tests/bits/find_cell_alt_6.cc b/tests/bits/find_cell_alt_6.cc new file mode 100644 index 0000000000..07a3741d95 --- /dev/null +++ b/tests/bits/find_cell_alt_6.cc @@ -0,0 +1,87 @@ +// ----------------------- find_cell_alt_6.cc -------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2006 by the deal.II authors and Ralf B. Schulz +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +// ----------------------- find_cell_alt_6.cc -------------------------- + + +// On a 2D mesh of the following structure look for the cells surrounding +// each vertex, using the find_active_cell_around_point with Mapping: +// +// x-----x-----x +// | | | +// | | | +// | | | +// x--x--x-----x +// | | | | +// x--x--x x +// | | | | +// x--x--x-----x + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + + + + +void check (Triangulation<2> &tria) +{ + const std::vector > &v = tria.get_vertices(); + MappingQ1<2> map; + + for(unsigned i=0; i::active_cell_iterator, Point<2> > + cell = GridTools::find_active_cell_around_point(map, tria, v[i]); + + deallog << "Vertex <" << v[i] << "> found in cell "; + for (unsigned int v=0; v::vertices_per_cell; ++v) + deallog << "<" << cell.first->vertex(v) << "> "; + deallog << " [local: " << cell.second << "]" << std::endl; + } +} + + +int main () +{ + std::ofstream logfile("find_cell_alt_6/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + try + { + Triangulation<2> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (1); + coarse_grid.begin_active()->set_refine_flag(); + coarse_grid.execute_coarsening_and_refinement(); + check (coarse_grid); + } + catch (const std::exception &exc) + { + // we shouldn't get here... + deallog << "Caught an error..." << std::endl; + deallog << exc.what() << std::endl; + } +} + + + diff --git a/tests/bits/find_cell_alt_6/cmp/generic b/tests/bits/find_cell_alt_6/cmp/generic new file mode 100644 index 0000000000..dab7ccfba1 --- /dev/null +++ b/tests/bits/find_cell_alt_6/cmp/generic @@ -0,0 +1,15 @@ + +DEAL::Vertex <0.00000 0.00000> found in cell <0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> [local: 0.00000 0.00000] +DEAL::Vertex <1.00000 0.00000> found in cell <0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> [local: 1.00000 0.00000] +DEAL::Vertex <0.00000 1.00000> found in cell <0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> [local: 0.00000 1.00000] +DEAL::Vertex <1.00000 1.00000> found in cell <0.500000 0.500000> <1.00000 0.500000> <0.500000 1.00000> <1.00000 1.00000> [local: 1.00000 1.00000] +DEAL::Vertex <0.500000 0.00000> found in cell <0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> [local: 0.00000 0.00000] +DEAL::Vertex <0.00000 0.500000> found in cell <0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> [local: 0.00000 0.00000] +DEAL::Vertex <1.00000 0.500000> found in cell <0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> [local: 1.00000 1.00000] +DEAL::Vertex <0.500000 1.00000> found in cell <0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> [local: 1.00000 1.00000] +DEAL::Vertex <0.500000 0.500000> found in cell <0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> [local: 0.00000 1.00000] +DEAL::Vertex <0.250000 0.00000> found in cell <0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> [local: 1.00000 0.00000] +DEAL::Vertex <0.00000 0.250000> found in cell <0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> [local: 0.00000 1.00000] +DEAL::Vertex <0.500000 0.250000> found in cell <0.250000 0.00000> <0.500000 0.00000> <0.250000 0.250000> <0.500000 0.250000> [local: 1.00000 1.00000] +DEAL::Vertex <0.250000 0.500000> found in cell <0.00000 0.250000> <0.250000 0.250000> <0.00000 0.500000> <0.250000 0.500000> [local: 1.00000 1.00000] +DEAL::Vertex <0.250000 0.250000> found in cell <0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> [local: 1.00000 1.00000] diff --git a/tests/bits/find_cells_adjacent_to_vertex_1.cc b/tests/bits/find_cells_adjacent_to_vertex_1.cc new file mode 100644 index 0000000000..91b91fee9d --- /dev/null +++ b/tests/bits/find_cells_adjacent_to_vertex_1.cc @@ -0,0 +1,85 @@ +// ---------------- find_cells_adjacent_to_vertex_1.cc ------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2006 by the deal.II authors and Ralf B. Schulz +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +// ---------------- find_cells_adjacent_to_vertex_1.cc ------------------ + + +// On a 2D mesh of the following structure look for the cells adjacent to +// each vertex: +// +// x-----x-----x +// | | | +// | | | +// | | | +// x--x--x-----x +// | | | | +// x--x--x x +// | | | | +// x--x--x-----x + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + + + + +void check (Triangulation<2> &tria) +{ + for(unsigned i=0; i::active_cell_iterator> + cells = GridTools::find_cells_adjacent_to_vertex(tria, i); + + deallog << "Vertex " << i << " at " << tria.get_vertices()[i] << ": " << cells.size() << " cells" << std::endl; + + for(unsigned c=0; c::vertices_per_cell; ++v) + deallog << "<" << cells[c]->vertex(v) << "> "; + deallog << std::endl; + } + } +} + + +int main () +{ + std::ofstream logfile("find_cells_adjacent_to_vertex_1/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + try + { + Triangulation<2> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (1); + coarse_grid.begin_active()->set_refine_flag(); + coarse_grid.execute_coarsening_and_refinement(); + check (coarse_grid); + } + catch (const std::exception &exc) + { + // we shouldn't get here... + deallog << "Caught an error..." << std::endl; + deallog << exc.what() << std::endl; + } +} + + + diff --git a/tests/bits/find_cells_adjacent_to_vertex_1/cmp/generic b/tests/bits/find_cells_adjacent_to_vertex_1/cmp/generic new file mode 100644 index 0000000000..cdeff8a936 --- /dev/null +++ b/tests/bits/find_cells_adjacent_to_vertex_1/cmp/generic @@ -0,0 +1,43 @@ + +DEAL::Vertex 0 at 0.00000 0.00000: 1 cells +DEAL::<0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> +DEAL::Vertex 1 at 1.00000 0.00000: 1 cells +DEAL::<0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> +DEAL::Vertex 2 at 0.00000 1.00000: 1 cells +DEAL::<0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> +DEAL::Vertex 3 at 1.00000 1.00000: 1 cells +DEAL::<0.500000 0.500000> <1.00000 0.500000> <0.500000 1.00000> <1.00000 1.00000> +DEAL::Vertex 4 at 0.500000 0.00000: 2 cells +DEAL::<0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> +DEAL::<0.250000 0.00000> <0.500000 0.00000> <0.250000 0.250000> <0.500000 0.250000> +DEAL::Vertex 5 at 0.00000 0.500000: 2 cells +DEAL::<0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> +DEAL::<0.00000 0.250000> <0.250000 0.250000> <0.00000 0.500000> <0.250000 0.500000> +DEAL::Vertex 6 at 1.00000 0.500000: 2 cells +DEAL::<0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> +DEAL::<0.500000 0.500000> <1.00000 0.500000> <0.500000 1.00000> <1.00000 1.00000> +DEAL::Vertex 7 at 0.500000 1.00000: 2 cells +DEAL::<0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> +DEAL::<0.500000 0.500000> <1.00000 0.500000> <0.500000 1.00000> <1.00000 1.00000> +DEAL::Vertex 8 at 0.500000 0.500000: 4 cells +DEAL::<0.500000 0.00000> <1.00000 0.00000> <0.500000 0.500000> <1.00000 0.500000> +DEAL::<0.00000 0.500000> <0.500000 0.500000> <0.00000 1.00000> <0.500000 1.00000> +DEAL::<0.500000 0.500000> <1.00000 0.500000> <0.500000 1.00000> <1.00000 1.00000> +DEAL::<0.250000 0.250000> <0.500000 0.250000> <0.250000 0.500000> <0.500000 0.500000> +DEAL::Vertex 9 at 0.250000 0.00000: 2 cells +DEAL::<0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> +DEAL::<0.250000 0.00000> <0.500000 0.00000> <0.250000 0.250000> <0.500000 0.250000> +DEAL::Vertex 10 at 0.00000 0.250000: 2 cells +DEAL::<0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> +DEAL::<0.00000 0.250000> <0.250000 0.250000> <0.00000 0.500000> <0.250000 0.500000> +DEAL::Vertex 11 at 0.500000 0.250000: 2 cells +DEAL::<0.250000 0.00000> <0.500000 0.00000> <0.250000 0.250000> <0.500000 0.250000> +DEAL::<0.250000 0.250000> <0.500000 0.250000> <0.250000 0.500000> <0.500000 0.500000> +DEAL::Vertex 12 at 0.250000 0.500000: 2 cells +DEAL::<0.00000 0.250000> <0.250000 0.250000> <0.00000 0.500000> <0.250000 0.500000> +DEAL::<0.250000 0.250000> <0.500000 0.250000> <0.250000 0.500000> <0.500000 0.500000> +DEAL::Vertex 13 at 0.250000 0.250000: 4 cells +DEAL::<0.00000 0.00000> <0.250000 0.00000> <0.00000 0.250000> <0.250000 0.250000> +DEAL::<0.250000 0.00000> <0.500000 0.00000> <0.250000 0.250000> <0.500000 0.250000> +DEAL::<0.00000 0.250000> <0.250000 0.250000> <0.00000 0.500000> <0.250000 0.500000> +DEAL::<0.250000 0.250000> <0.500000 0.250000> <0.250000 0.500000> <0.500000 0.500000> diff --git a/tests/bits/find_cells_adjacent_to_vertex_2.cc b/tests/bits/find_cells_adjacent_to_vertex_2.cc new file mode 100644 index 0000000000..d8e3e0080e --- /dev/null +++ b/tests/bits/find_cells_adjacent_to_vertex_2.cc @@ -0,0 +1,74 @@ +// ---------------- find_cells_adjacent_to_vertex_2.cc ------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2006 by the deal.II authors and Ralf B. Schulz +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +// ---------------- find_cells_adjacent_to_vertex_2.cc ------------------ + + +// Same as the first test, but on a 3D grid of the same structure + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + + + + +void check (Triangulation<3> &tria) +{ + for(unsigned i=0; i::active_cell_iterator> + cells = GridTools::find_cells_adjacent_to_vertex(tria, i); + + deallog << "Vertex " << i << " at " << tria.get_vertices()[i] << ": " << cells.size() << " cells" << std::endl; + + for(unsigned c=0; c::vertices_per_cell; ++v) + deallog << "<" << cells[c]->vertex(v) << "> "; + deallog << std::endl; + } + } +} + + +int main () +{ + std::ofstream logfile("find_cells_adjacent_to_vertex_2/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + try + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (1); + coarse_grid.begin_active()->set_refine_flag(); + coarse_grid.execute_coarsening_and_refinement(); + check (coarse_grid); + } + catch (const std::exception &exc) + { + // we shouldn't get here... + deallog << "Caught an error..." << std::endl; + deallog << exc.what() << std::endl; + } +} + + + diff --git a/tests/bits/find_cells_adjacent_to_vertex_2/cmp/generic b/tests/bits/find_cells_adjacent_to_vertex_2/cmp/generic new file mode 100644 index 0000000000..b63e40d23c --- /dev/null +++ b/tests/bits/find_cells_adjacent_to_vertex_2/cmp/generic @@ -0,0 +1,167 @@ + +DEAL::Vertex 0 at 0.00000 0.00000 0.00000: 1 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::Vertex 1 at 1.00000 0.00000 0.00000: 1 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::Vertex 2 at 0.00000 1.00000 0.00000: 1 cells +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::Vertex 3 at 1.00000 1.00000 0.00000: 1 cells +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::Vertex 4 at 0.00000 0.00000 1.00000: 1 cells +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::Vertex 5 at 1.00000 0.00000 1.00000: 1 cells +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::Vertex 6 at 0.00000 1.00000 1.00000: 1 cells +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::Vertex 7 at 1.00000 1.00000 1.00000: 1 cells +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 8 at 0.500000 0.00000 0.00000: 2 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::Vertex 9 at 0.00000 0.500000 0.00000: 2 cells +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::Vertex 10 at 0.00000 0.00000 0.500000: 2 cells +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::Vertex 11 at 1.00000 0.500000 0.00000: 2 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::Vertex 12 at 1.00000 0.00000 0.500000: 2 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::Vertex 13 at 0.500000 1.00000 0.00000: 2 cells +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::Vertex 14 at 0.00000 1.00000 0.500000: 2 cells +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::Vertex 15 at 1.00000 1.00000 0.500000: 2 cells +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 16 at 0.500000 0.00000 1.00000: 2 cells +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::Vertex 17 at 0.00000 0.500000 1.00000: 2 cells +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::Vertex 18 at 1.00000 0.500000 1.00000: 2 cells +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 19 at 0.500000 1.00000 1.00000: 2 cells +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 20 at 0.500000 0.00000 0.500000: 4 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::Vertex 21 at 0.500000 0.500000 0.00000: 4 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::Vertex 22 at 0.00000 0.500000 0.500000: 4 cells +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::Vertex 23 at 1.00000 0.500000 0.500000: 4 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 24 at 0.500000 1.00000 0.500000: 4 cells +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 25 at 0.500000 0.500000 1.00000: 4 cells +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::Vertex 26 at 0.500000 0.500000 0.500000: 8 cells +DEAL::<0.500000 0.00000 0.00000> <1.00000 0.00000 0.00000> <0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> +DEAL::<0.00000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.00000 1.00000 0.00000> <0.500000 1.00000 0.00000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> +DEAL::<0.500000 0.500000 0.00000> <1.00000 0.500000 0.00000> <0.500000 1.00000 0.00000> <1.00000 1.00000 0.00000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> +DEAL::<0.00000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 0.00000 1.00000> <0.500000 0.00000 1.00000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> +DEAL::<0.500000 0.00000 0.500000> <1.00000 0.00000 0.500000> <0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 0.00000 1.00000> <1.00000 0.00000 1.00000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> +DEAL::<0.00000 0.500000 0.500000> <0.500000 0.500000 0.500000> <0.00000 1.00000 0.500000> <0.500000 1.00000 0.500000> <0.00000 0.500000 1.00000> <0.500000 0.500000 1.00000> <0.00000 1.00000 1.00000> <0.500000 1.00000 1.00000> +DEAL::<0.500000 0.500000 0.500000> <1.00000 0.500000 0.500000> <0.500000 1.00000 0.500000> <1.00000 1.00000 0.500000> <0.500000 0.500000 1.00000> <1.00000 0.500000 1.00000> <0.500000 1.00000 1.00000> <1.00000 1.00000 1.00000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 27 at 0.250000 0.00000 0.00000: 2 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::Vertex 28 at 0.00000 0.250000 0.00000: 2 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::Vertex 29 at 0.00000 0.00000 0.250000: 2 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::Vertex 30 at 0.250000 0.00000 0.500000: 2 cells +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::Vertex 31 at 0.500000 0.00000 0.250000: 2 cells +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::Vertex 32 at 0.500000 0.250000 0.00000: 2 cells +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::Vertex 33 at 0.250000 0.500000 0.00000: 2 cells +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::Vertex 34 at 0.00000 0.500000 0.250000: 2 cells +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::Vertex 35 at 0.00000 0.250000 0.500000: 2 cells +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::Vertex 36 at 0.500000 0.250000 0.500000: 2 cells +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 37 at 0.250000 0.500000 0.500000: 2 cells +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 38 at 0.500000 0.500000 0.250000: 2 cells +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 39 at 0.250000 0.00000 0.250000: 4 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::Vertex 40 at 0.250000 0.250000 0.00000: 4 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::Vertex 41 at 0.00000 0.250000 0.250000: 4 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::Vertex 42 at 0.500000 0.250000 0.250000: 4 cells +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 43 at 0.250000 0.500000 0.250000: 4 cells +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 44 at 0.250000 0.250000 0.500000: 4 cells +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> +DEAL::Vertex 45 at 0.250000 0.250000 0.250000: 8 cells +DEAL::<0.00000 0.00000 0.00000> <0.250000 0.00000 0.00000> <0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> +DEAL::<0.250000 0.00000 0.00000> <0.500000 0.00000 0.00000> <0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> +DEAL::<0.00000 0.250000 0.00000> <0.250000 0.250000 0.00000> <0.00000 0.500000 0.00000> <0.250000 0.500000 0.00000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> +DEAL::<0.250000 0.250000 0.00000> <0.500000 0.250000 0.00000> <0.250000 0.500000 0.00000> <0.500000 0.500000 0.00000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> +DEAL::<0.00000 0.00000 0.250000> <0.250000 0.00000 0.250000> <0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.00000 0.500000> <0.250000 0.00000 0.500000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> +DEAL::<0.250000 0.00000 0.250000> <0.500000 0.00000 0.250000> <0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.00000 0.500000> <0.500000 0.00000 0.500000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> +DEAL::<0.00000 0.250000 0.250000> <0.250000 0.250000 0.250000> <0.00000 0.500000 0.250000> <0.250000 0.500000 0.250000> <0.00000 0.250000 0.500000> <0.250000 0.250000 0.500000> <0.00000 0.500000 0.500000> <0.250000 0.500000 0.500000> +DEAL::<0.250000 0.250000 0.250000> <0.500000 0.250000 0.250000> <0.250000 0.500000 0.250000> <0.500000 0.500000 0.250000> <0.250000 0.250000 0.500000> <0.500000 0.250000 0.500000> <0.250000 0.500000 0.500000> <0.500000 0.500000 0.500000> diff --git a/tests/bits/find_closest_vertex_1.cc b/tests/bits/find_closest_vertex_1.cc new file mode 100644 index 0000000000..e587ffbfec --- /dev/null +++ b/tests/bits/find_closest_vertex_1.cc @@ -0,0 +1,66 @@ +//----------------------- find_closest_vertex_1.cc ---------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2006 by the deal.II authors and Ralf B. Schulz +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//----------------------- find_closest_vertex_1.cc ---------------------- + + +// take a 3d mesh, take all vertices, shift them a little bit and check that +// we correctly identify the closest vertex position +// The result should be an increasing sequence of numbers + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include + +#include + + + + +void check (Triangulation<3> &tria) +{ + const std::vector > &v = tria.get_vertices(); + for(unsigned i=0; i(0.01, -0.01, 0.01)) << "] "; + + deallog << std::endl; +} + + +int main () +{ + std::ofstream logfile("find_closest_vertex_1/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + try + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_cube (coarse_grid); + coarse_grid.refine_global (3); + check (coarse_grid); + } + catch (const std::exception &exc) + { + // we shouldn't get here... + deallog << "Caught an error..." << std::endl; + deallog << exc.what() << std::endl; + } +} + + + diff --git a/tests/bits/find_closest_vertex_1/cmp/generic b/tests/bits/find_closest_vertex_1/cmp/generic new file mode 100644 index 0000000000..0aaa5b23eb --- /dev/null +++ b/tests/bits/find_closest_vertex_1/cmp/generic @@ -0,0 +1,2 @@ + +DEAL::[0] [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75] [76] [77] [78] [79] [80] [81] [82] [83] [84] [85] [86] [87] [88] [89] [90] [91] [92] [93] [94] [95] [96] [97] [98] [99] [100] [101] [102] [103] [104] [105] [106] [107] [108] [109] [110] [111] [112] [113] [114] [115] [116] [117] [118] [119] [120] [121] [122] [123] [124] [125] [126] [127] [128] [129] [130] [131] [132] [133] [134] [135] [136] [137] [138] [139] [140] [141] [142] [143] [144] [145] [146] [147] [148] [149] [150] [151] [152] [153] [154] [155] [156] [157] [158] [159] [160] [161] [162] [163] [164] [165] [166] [167] [168] [169] [170] [171] [172] [173] [174] [175] [176] [177] [178] [179] [180] [181] [182] [183] [184] [185] [186] [187] [188] [189] [190] [191] [192] [193] [194] [195] [196] [197] [198] [199] [200] [201] [202] [203] [204] [205] [206] [207] [208] [209] [210] [211] [212] [213] [214] [215] [216] [217] [218] [219] [220] [221] [222] [223] [224] [225] [226] [227] [228] [229] [230] [231] [232] [233] [234] [235] [236] [237] [238] [239] [240] [241] [242] [243] [244] [245] [246] [247] [248] [249] [250] [251] [252] [253] [254] [255] [256] [257] [258] [259] [260] [261] [262] [263] [264] [265] [266] [267] [268] [269] [270] [271] [272] [273] [274] [275] [276] [277] [278] [279] [280] [281] [282] [283] [284] [285] [286] [287] [288] [289] [290] [291] [292] [293] [294] [295] [296] [297] [298] [299] [300] [301] [302] [303] [304] [305] [306] [307] [308] [309] [310] [311] [312] [313] [314] [315] [316] [317] [318] [319] [320] [321] [322] [323] [324] [325] [326] [327] [328] [329] [330] [331] [332] [333] [334] [335] [336] [337] [338] [339] [340] [341] [342] [343] [344] [345] [346] [347] [348] [349] [350] [351] [352] [353] [354] [355] [356] [357] [358] [359] [360] [361] [362] [363] [364] [365] [366] [367] [368] [369] [370] [371] [372] [373] [374] [375] [376] [377] [378] [379] [380] [381] [382] [383] [384] [385] [386] [387] [388] [389] [390] [391] [392] [393] [394] [395] [396] [397] [398] [399] [400] [401] [402] [403] [404] [405] [406] [407] [408] [409] [410] [411] [412] [413] [414] [415] [416] [417] [418] [419] [420] [421] [422] [423] [424] [425] [426] [427] [428] [429] [430] [431] [432] [433] [434] [435] [436] [437] [438] [439] [440] [441] [442] [443] [444] [445] [446] [447] [448] [449] [450] [451] [452] [453] [454] [455] [456] [457] [458] [459] [460] [461] [462] [463] [464] [465] [466] [467] [468] [469] [470] [471] [472] [473] [474] [475] [476] [477] [478] [479] [480] [481] [482] [483] [484] [485] [486] [487] [488] [489] [490] [491] [492] [493] [494] [495] [496] [497] [498] [499] [500] [501] [502] [503] [504] [505] [506] [507] [508] [509] [510] [511] [512] [513] [514] [515] [516] [517] [518] [519] [520] [521] [522] [523] [524] [525] [526] [527] [528] [529] [530] [531] [532] [533] [534] [535] [536] [537] [538] [539] [540] [541] [542] [543] [544] [545] [546] [547] [548] [549] [550] [551] [552] [553] [554] [555] [556] [557] [558] [559] [560] [561] [562] [563] [564] [565] [566] [567] [568] [569] [570] [571] [572] [573] [574] [575] [576] [577] [578] [579] [580] [581] [582] [583] [584] [585] [586] [587] [588] [589] [590] [591] [592] [593] [594] [595] [596] [597] [598] [599] [600] [601] [602] [603] [604] [605] [606] [607] [608] [609] [610] [611] [612] [613] [614] [615] [616] [617] [618] [619] [620] [621] [622] [623] [624] [625] [626] [627] [628] [629] [630] [631] [632] [633] [634] [635] [636] [637] [638] [639] [640] [641] [642] [643] [644] [645] [646] [647] [648] [649] [650] [651] [652] [653] [654] [655] [656] [657] [658] [659] [660] [661] [662] [663] [664] [665] [666] [667] [668] [669] [670] [671] [672] [673] [674] [675] [676] [677] [678] [679] [680] [681] [682] [683] [684] [685] [686] [687] [688] [689] [690] [691] [692] [693] [694] [695] [696] [697] [698] [699] [700] [701] [702] [703] [704] [705] [706] [707] [708] [709] [710] [711] [712] [713] [714] [715] [716] [717] [718] [719] [720] [721] [722] [723] [724] [725] [726] [727] [728]