// $Id$
// Version: $Name$
//
-// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
+// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*/
template <int dim, int spacedim>
static
- void partition_triangulation (const unsigned int n_partitions,
- Triangulation<dim, spacedim> &triangulation);
+ void
+ partition_triangulation (const unsigned int n_partitions,
+ Triangulation<dim, spacedim> &triangulation);
/**
* This function does the same as the
*/
template <int dim, int spacedim>
static
- void partition_triangulation (const unsigned int n_partitions,
- const SparsityPattern &cell_connection_graph,
- Triangulation<dim,spacedim> &triangulation);
+ void
+ partition_triangulation (const unsigned int n_partitions,
+ const SparsityPattern &cell_connection_graph,
+ Triangulation<dim,spacedim> &triangulation);
/**
* For each active cell, return in the
create_union_triangulation (const Triangulation<dim, spacedim> &triangulation_1,
const Triangulation<dim, spacedim> &triangulation_2,
Triangulation<dim, spacedim> &result);
+
+ /**
+ * Given a triangulation and a
+ * list of cells whose children
+ * have become distorted as a
+ * result of mesh refinement, try
+ * to fix these cells up by
+ * moving the center node around.
+ *
+ * The function returns a list of
+ * cells with distorted children
+ * that couldn't be fixed up for
+ * whatever reason. The returned
+ * list is therefore a subset of
+ * the input argument.
+ *
+ * For a definition of the
+ * concept of distorted cells,
+ * see the
+ * @ref GlossDistorted "glossary entry".
+ * The first argument passed to the
+ * current function is typically
+ * the exception thrown by the
+ * Triangulation::execute_coarsening_and_refinement
+ * function.
+ */
+ template <int dim, int spacedim>
+ static
+ typename Triangulation<dim,spacedim>::DistortedCellList
+ fix_up_distorted_child_cells (const typename Triangulation<dim,spacedim>::DistortedCellList &distorted_cells,
+ Triangulation<dim,spacedim> &triangulation);
/**
* Exception
#include <multigrid/mg_dof_accessor.h>
#include <cmath>
+#include <numeric>
+
DEAL_II_NAMESPACE_OPEN
}
+namespace internal
+{
+ namespace GridTools
+ {
+ namespace FixUpDistortedChildCells
+ {
+ template <int dim, int spacedim>
+ double
+ objective_function (const typename dealii::Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Point<spacedim> &cell_mid_point)
+ {
+ // everything below is wrong
+ // if not for the following
+ // condition
+ Assert (cell->refinement_case() == RefinementCase<dim>::isotropic_refinement,
+ ExcNotImplemented());
+ // first calculate the
+ // average jacobian
+ // determinant for the parent
+ // cell
+ Point<spacedim> parent_vertices
+ [GeometryInfo<dim>::vertices_per_cell];
+ double parent_determinants
+ [GeometryInfo<dim>::vertices_per_cell];
+
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ parent_vertices[i] = cell->vertex(i);
+
+ GeometryInfo<dim>::jacobian_determinants_at_vertices (parent_vertices,
+ parent_determinants);
+
+ const double average_parent_jacobian
+ = std::accumulate (&parent_determinants[0],
+ &parent_determinants[GeometryInfo<dim>::vertices_per_cell],
+ 0.);
+
+ // now do the same
+ // computation for the
+ // children where we use the
+ // given location for the
+ // cell mid point instead of
+ // the one the triangulation
+ // currently reports
+ Point<spacedim> child_vertices
+ [GeometryInfo<dim>::max_children_per_cell]
+ [GeometryInfo<dim>::vertices_per_cell];
+ double child_determinants
+ [GeometryInfo<dim>::max_children_per_cell]
+ [GeometryInfo<dim>::vertices_per_cell];
+
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ child_vertices[c][i] = cell->child(c)->vertex(i);
+
+ // replace mid-cell
+ // vertex. note that for
+ // child i, the mid-cell
+ // vertex happens to have the
+ // number
+ // max_children_per_cell-i
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ child_vertices[c][GeometryInfo<dim>::max_children_per_cell-c-1]
+ = cell_mid_point;
+
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ GeometryInfo<dim>::jacobian_determinants_at_vertices (child_vertices[c],
+ child_determinants[c]);
+
+ // on a uniformly refined
+ // hypercube cell, the child
+ // jacobians should all be
+ // smaller by a factor of
+ // 2^dim than the ones of the
+ // parent. as a consequence,
+ // we'll use the squared
+ // deviation from this ideal
+ // value as an objective
+ // function
+ double objective = 0;
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ objective += std::pow (child_determinants[c][i] -
+ average_parent_jacobian/std::pow(2.,1.*dim),
+ 2);
+
+ return objective;
+ }
+
+ }
+ }
+}
+
+
+
+template <int dim, int spacedim>
+typename Triangulation<dim,spacedim>::DistortedCellList
+GridTools::
+fix_up_distorted_child_cells (const typename Triangulation<dim,spacedim>::DistortedCellList &distorted_cells,
+ Triangulation<dim,spacedim> &/*triangulation*/)
+{
+ typename Triangulation<dim,spacedim>::DistortedCellList unfixable_subset;
+
+ // loop over all cells that we have
+ // to fix up
+ for (typename std::list<typename Triangulation<dim,spacedim>::active_cell_iterator>::const_iterator
+ cell_ptr = distorted_cells.distorted_cells.begin();
+ cell_ptr != distorted_cells.distorted_cells.end(); ++cell_ptr)
+ {
+ const typename Triangulation<dim,spacedim>::cell_iterator
+ cell = *cell_ptr;
+
+ // right now we can only deal
+ // with cells that have been
+ // refined isotropically
+ // because that is the only
+ // case where we have a cell
+ // mid-point that can be moved
+ // around without having to
+ // consider boundary
+ // information
+ Assert (cell->has_children(), ExcInternalError());
+ Assert (cell->refinement_case() == RefinementCase<dim>::isotropic_refinement,
+ ExcNotImplemented());
+
+ // get the current location of
+ // the cell mid-vertex:
+ Point<spacedim> cell_mid_point
+ = cell->child(0)->vertex (GeometryInfo<dim>::max_children_per_cell-1);
+
+ // now do a few steepest
+ // descent steps to reduce the
+ // objective function
+ unsigned int iteration = 0;
+ do
+ {
+ // choose a step length
+ // that is initially 1/10
+ // of the child cells'
+ // diameter, and a sequence
+ // whose sum does not
+ // converge (to avoid
+ // premature termination of
+ // the iteration)
+ const double step_length = cell->diameter() / 10 / (iteration + 1);
+
+ // compute the objective
+ // function and its derivative
+ const double val = internal::GridTools::FixUpDistortedChildCells::
+ objective_function<dim,spacedim> (cell, cell_mid_point);
+
+ Tensor<1,dim> gradient;
+ for (unsigned int d=0; d<dim; ++d)
+ {
+ Point<dim> h;
+ h[d] = step_length/2;
+ gradient[d] = (internal::GridTools::FixUpDistortedChildCells::
+ objective_function<dim,spacedim> (cell,
+ cell_mid_point + h)
+ -
+ internal::GridTools::FixUpDistortedChildCells::
+ objective_function<dim,spacedim> (cell,
+ cell_mid_point - h))
+ /
+ step_length;
+ }
+
+ // so we need to go in
+ // direction -gradient. the
+ // optimal value of the
+ // objective function is
+ // zero, so assuming that
+ // the model is quadratic
+ // we would have to go
+ // -2*val/||gradient|| in
+ // this direction, make
+ // sure we go at most
+ // step_length into this
+ // direction
+ cell_mid_point -= std::min(2*val / (gradient*gradient),
+ step_length / gradient.norm()) *
+ gradient;
+
+ ++iteration;
+ }
+ while (iteration < 10);
+
+
+ // verify that the new location
+ // is indeed better than the
+ // one before in terms of the
+ // minimal jacobian determinant
+ double old_min_jacobian, new_min_jacobian;
+
+ for (unsigned int test=0; test<2; ++test)
+ {
+ Point<spacedim> child_vertices
+ [GeometryInfo<dim>::max_children_per_cell]
+ [GeometryInfo<dim>::vertices_per_cell];
+ double child_determinants
+ [GeometryInfo<dim>::max_children_per_cell]
+ [GeometryInfo<dim>::vertices_per_cell];
+
+ if (test == 1)
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ child_vertices[c][i] = cell->child(c)->vertex(i);
+
+ // replace mid-cell
+ // vertex. note that for
+ // child i, the mid-cell
+ // vertex happens to have the
+ // number
+ // max_children_per_cell-i
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ child_vertices[c][GeometryInfo<dim>::max_children_per_cell-c-1]
+ = cell_mid_point;
+
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ GeometryInfo<dim>::jacobian_determinants_at_vertices (child_vertices[c],
+ child_determinants[c]);
+
+ double min = child_determinants[0][0];
+ for (unsigned int c=0; c<cell->n_children(); ++c)
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ min = std::min (min,
+ child_determinants[c][i]);
+
+ if (test == 0)
+ old_min_jacobian = min;
+ else
+ new_min_jacobian = min;
+ }
+
+ // if new minimum jacobian is
+ // better than before, and if
+ // it is positive, then set the
+ // new mid point. otherwise
+ // return this cell as one of
+ // those that can't apparently
+ // be fixed
+ if ((new_min_jacobian > old_min_jacobian)
+ &&
+ (new_min_jacobian > 0))
+ cell->child(0)->vertex (GeometryInfo<dim>::max_children_per_cell-1)
+ = cell_mid_point;
+ else
+ unfixable_subset.distorted_cells.push_back (cell);
+ }
+
+ return unfixable_subset;
+}
+
// explicit instantiations
const Triangulation<deal_II_dimension> &triangulation_2,
Triangulation<deal_II_dimension> &result);
+template
+Triangulation<deal_II_dimension,deal_II_dimension>::DistortedCellList
+GridTools::
+fix_up_distorted_child_cells (const Triangulation<deal_II_dimension,deal_II_dimension>::DistortedCellList &distorted_cells,
+ Triangulation<deal_II_dimension,deal_II_dimension> &triangulation);
+
+
#if deal_II_dimension != 3
Triangulation<deal_II_dimension, deal_II_dimension+1> &);
-
#endif
* Triangulation::execute_coarsening_and_refinement can then decide
* what to do with this situation.
*
+ * One way to deal with this problem is to use the
+ * GridTools::fix_up_distorted_child_cells function that attempts to
+ * fix up exactly these cells if possible by moving around the node at
+ * the center of the cell.
+ *
*
* <dt class="glossary">@anchor GlossFaceOrientation Face orientation</dt>
* <dd>In a triangulation, the normal vector to a face
--- /dev/null
+//---------------------------- distorted_cells_04.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2003, 2004, 2005, 2009 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.
+//
+//---------------------------- distorted_cells_04.cc ---------------------------
+
+
+// like _03, but catch the exception and pass it to GridTools::fix_up_distorted_child_cells
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <base/quadrature_lib.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_reordering.h>
+#include <grid/grid_tools.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_boundary.h>
+#include <grid/grid_out.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+#include <fstream>
+
+
+template <int dim>
+class MyBoundary : public Boundary<dim>
+{
+ virtual Point<dim>
+ get_new_point_on_line (const typename Triangulation<dim>::line_iterator &line) const
+ {
+ deallog << "Finding point between "
+ << line->vertex(0) << " and "
+ << line->vertex(1) << std::endl;
+
+ // in 2d, find a point that
+ // lies on the opposite side
+ // of the quad. in 3d, choose
+ // the midpoint of the edge
+ if (dim == 2)
+ return Point<dim>(0,0.75);
+ else
+ return (line->vertex(0) + line->vertex(1)) / 2;
+ }
+
+ virtual Point<dim>
+ get_new_point_on_quad (const typename Triangulation<dim>::quad_iterator &quad) const
+ {
+ deallog << "Finding point between "
+ << quad->vertex(0) << " and "
+ << quad->vertex(1) << " and "
+ << quad->vertex(2) << " and "
+ << quad->vertex(3) << std::endl;
+
+ return Point<dim>(0,0,.75);
+ }
+};
+
+
+
+template <int dim>
+void check ()
+{
+ MyBoundary<dim> my_boundary;
+
+ // create a single square/cube
+ Triangulation<dim> coarse_grid;
+ GridGenerator::hyper_cube (coarse_grid, -1, 1);
+
+ // set bottom face to use MyBoundary
+ for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+ if (coarse_grid.begin_active()->face(f)->center()[dim-1] == -1)
+ coarse_grid.begin_active()->face(f)->set_boundary_indicator (1);
+ coarse_grid.set_boundary (1, my_boundary);
+
+ // now try to refine this one
+ // cell. we should get an exception
+ try
+ {
+ coarse_grid.begin_active()->set_refine_flag ();
+ coarse_grid.execute_coarsening_and_refinement ();
+ }
+ catch (typename Triangulation<dim>::DistortedCellList &dcv)
+ {
+ typename Triangulation<dim>::DistortedCellList
+ subset = GridTools::fix_up_distorted_child_cells (dcv,
+ coarse_grid);
+ Assert (subset.distorted_cells.size() == 0,
+ ExcInternalError());
+ }
+
+ Assert (coarse_grid.n_levels() == 2, ExcInternalError());
+ Assert (coarse_grid.n_active_cells() == 1<<dim, ExcInternalError());
+
+ // output the coordinates of the
+ // child cells
+ GridOut().write_gnuplot (coarse_grid, deallog.get_file_stream());
+}
+
+
+int main ()
+{
+ std::ofstream logfile("distorted_cells_04/output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+ deallog.threshold_double(1.e-10);
+
+ check<2> ();
+ check<3> ();
+}
+
+
+
--- /dev/null
+4 1 0 0 0
+1 0 0 0
+2 1 0 0
+3 0 1 0
+4 1 1 0
+1 0 quad 1 2 3 4
--- /dev/null
+8 1 0 0 0
+1 0 0 0
+2 1 0 0
+3 0 1 0
+4 1 1 0
+5 0 0 1
+6 1 0 1
+7 0 1 1
+8 1 1 1
+1 0 hex 1 2 3 4 5 6 7 8
--- /dev/null
+
+DEAL::Finding point between -1.00000 -1.00000 and 1.00000 -1.00000
+-1.00000 -1.00000 1 0
+0.00000 0.750000 1 0
+-3.50619e-10 0.846953 1 0
+-1.00000 0.00000 1 0
+-1.00000 -1.00000 1 0
+
+
+0.00000 0.750000 1 0
+1.00000 -1.00000 1 0
+1.00000 0.00000 1 0
+-3.50619e-10 0.846953 1 0
+0.00000 0.750000 1 0
+
+
+-1.00000 0.00000 1 0
+-3.50619e-10 0.846953 1 0
+0.00000 1.00000 1 0
+-1.00000 1.00000 1 0
+-1.00000 0.00000 1 0
+
+
+-3.50619e-10 0.846953 1 0
+1.00000 0.00000 1 0
+1.00000 1.00000 1 0
+0.00000 1.00000 1 0
+-3.50619e-10 0.846953 1 0
+
+
+DEAL::Finding point between -1.00000 -1.00000 -1.00000 and 1.00000 -1.00000 -1.00000 and -1.00000 1.00000 -1.00000 and 1.00000 1.00000 -1.00000
+-1.00000 -1.00000 -1.00000 1 0
+0.00000 -1.00000 -1.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+-1.00000 -1.00000 0.00000 1 0
+-1.00000 -1.00000 -1.00000 1 0
+
+-1.00000 0.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+-1.00000 0.00000 0.00000 1 0
+-1.00000 0.00000 -1.00000 1 0
+
+-1.00000 -1.00000 -1.00000 1 0
+-1.00000 0.00000 -1.00000 1 0
+
+0.00000 -1.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+
+-1.00000 -1.00000 0.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+0.00000 -1.00000 -1.00000 1 0
+1.00000 -1.00000 -1.00000 1 0
+1.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+1.00000 0.00000 -1.00000 1 0
+1.00000 0.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 0.00000 0.750000 1 0
+
+0.00000 -1.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+
+1.00000 -1.00000 -1.00000 1 0
+1.00000 0.00000 -1.00000 1 0
+
+1.00000 -1.00000 0.00000 1 0
+1.00000 0.00000 0.00000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+
+-1.00000 0.00000 -1.00000 1 0
+0.00000 0.00000 0.750000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+-1.00000 0.00000 0.00000 1 0
+-1.00000 0.00000 -1.00000 1 0
+
+-1.00000 1.00000 -1.00000 1 0
+0.00000 1.00000 -1.00000 1 0
+0.00000 1.00000 0.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+-1.00000 1.00000 -1.00000 1 0
+
+-1.00000 0.00000 -1.00000 1 0
+-1.00000 1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+0.00000 1.00000 -1.00000 1 0
+
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 1.00000 0.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+1.00000 0.00000 -1.00000 1 0
+1.00000 0.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 0.00000 0.750000 1 0
+
+0.00000 1.00000 -1.00000 1 0
+1.00000 1.00000 -1.00000 1 0
+1.00000 1.00000 0.00000 1 0
+0.00000 1.00000 0.00000 1 0
+0.00000 1.00000 -1.00000 1 0
+
+0.00000 0.00000 0.750000 1 0
+0.00000 1.00000 -1.00000 1 0
+
+1.00000 0.00000 -1.00000 1 0
+1.00000 1.00000 -1.00000 1 0
+
+1.00000 0.00000 0.00000 1 0
+1.00000 1.00000 0.00000 1 0
+
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 1.00000 0.00000 1 0
+
+-1.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+0.00000 -1.00000 1.00000 1 0
+-1.00000 -1.00000 1.00000 1 0
+-1.00000 -1.00000 0.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+-1.00000 -1.00000 0.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+
+0.00000 -1.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+
+-1.00000 -1.00000 1.00000 1 0
+-1.00000 0.00000 1.00000 1 0
+
+0.00000 -1.00000 0.00000 1 0
+1.00000 -1.00000 0.00000 1 0
+1.00000 -1.00000 1.00000 1 0
+0.00000 -1.00000 1.00000 1 0
+0.00000 -1.00000 0.00000 1 0
+
+1.07919e-09 -5.43281e-10 0.876731 1 0
+1.00000 0.00000 0.00000 1 0
+1.00000 0.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+
+0.00000 -1.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+
+1.00000 -1.00000 0.00000 1 0
+1.00000 0.00000 0.00000 1 0
+
+1.00000 -1.00000 1.00000 1 0
+1.00000 0.00000 1.00000 1 0
+
+0.00000 -1.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 1.00000 1 0
+-1.00000 0.00000 0.00000 1 0
+
+-1.00000 1.00000 0.00000 1 0
+0.00000 1.00000 0.00000 1 0
+0.00000 1.00000 1.00000 1 0
+-1.00000 1.00000 1.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+
+-1.00000 0.00000 0.00000 1 0
+-1.00000 1.00000 0.00000 1 0
+
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 1.00000 0.00000 1 0
+
+0.00000 0.00000 1.00000 1 0
+0.00000 1.00000 1.00000 1 0
+
+-1.00000 0.00000 1.00000 1 0
+-1.00000 1.00000 1.00000 1 0
+
+1.07919e-09 -5.43281e-10 0.876731 1 0
+1.00000 0.00000 0.00000 1 0
+1.00000 0.00000 1.00000 1 0
+0.00000 0.00000 1.00000 1 0
+1.07919e-09 -5.43281e-10 0.876731 1 0
+
+0.00000 1.00000 0.00000 1 0
+1.00000 1.00000 0.00000 1 0
+1.00000 1.00000 1.00000 1 0
+0.00000 1.00000 1.00000 1 0
+0.00000 1.00000 0.00000 1 0
+
+1.07919e-09 -5.43281e-10 0.876731 1 0
+0.00000 1.00000 0.00000 1 0
+
+1.00000 0.00000 0.00000 1 0
+1.00000 1.00000 0.00000 1 0
+
+1.00000 0.00000 1.00000 1 0
+1.00000 1.00000 1.00000 1 0
+
+0.00000 0.00000 1.00000 1 0
+0.00000 1.00000 1.00000 1 0
+