From e624fe10e308b4bbb0ec44494105ae64c7ccb383 Mon Sep 17 00:00:00 2001 From: wolf Date: Fri, 17 Oct 2003 19:46:59 +0000 Subject: [PATCH] Invent function neighbor_child_on_subface. Add respective test. git-svn-id: https://svn.dealii.org/trunk@8102 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_accessor.h | 240 +++++++++++------- deal.II/deal.II/include/grid/tria_accessor.h | 86 +++++-- .../include/grid/tria_accessor.templates.h | 65 +++++ deal.II/deal.II/source/dofs/dof_accessor.cc | 18 ++ deal.II/deal.II/source/grid/tria_accessor.cc | 87 +++++++ tests/bits/Makefile | 3 +- tests/bits/mesh_3d_15.cc | 146 +++++++++++ .../bits/mesh_3d_15.output | 19 ++ .../bits/mesh_3d_15.output | 19 ++ .../bits/mesh_3d_15.output | 19 ++ 10 files changed, 580 insertions(+), 122 deletions(-) create mode 100644 tests/bits/mesh_3d_15.cc create mode 100644 tests/results/i686-pc-linux-gnu+gcc2.95/bits/mesh_3d_15.output create mode 100644 tests/results/i686-pc-linux-gnu+gcc3.2/bits/mesh_3d_15.output create mode 100644 tests/results/sparc-sun-solaris2.7+gcc2.95/bits/mesh_3d_15.output diff --git a/deal.II/deal.II/include/dofs/dof_accessor.h b/deal.II/deal.II/include/dofs/dof_accessor.h index ba7f4edcc8..c008cdb963 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.h @@ -1090,9 +1090,10 @@ class DoFCellAccessor : public DoFObjectAccessor { public: /** - * Declare the data type that this accessor - * class expects to get passed from the - * iterator classes. + * Declare the data type that + * this accessor class expects to + * get passed from the iterator + * classes. */ typedef typename DoFObjectAccessor::AccessorData AccessorData; @@ -1105,21 +1106,25 @@ class DoFCellAccessor : public DoFObjectAccessor const AccessorData *local_data); /** - * Return the @p{i}th neighbor as a DoF cell - * iterator. This function is needed since - * the neighbor function of the base - * class returns a cell accessor without - * access to the DoF data. + * Return the @p{i}th neighbor as + * a DoF cell iterator. This + * function is needed since the + * neighbor function of the base + * class returns a cell accessor + * without access to the DoF + * data. */ TriaIterator > neighbor (const unsigned int) const; /** - * Return the @p{i}th child as a DoF cell - * iterator. This function is needed since - * the child function of the base - * class returns a cell accessor without - * access to the DoF data. + * Return the @p{i}th child as a + * DoF cell iterator. This + * function is needed since the + * child function of the base + * class returns a cell accessor + * without access to the DoF + * data. */ TriaIterator > child (const unsigned int) const; @@ -1135,43 +1140,70 @@ class DoFCellAccessor : public DoFObjectAccessor face (const unsigned int i) const; /** - * Return the interpolation of the given - * finite element function to the present - * cell. In the simplest case, the cell - * is a terminal one, i.e. has no children; - * then, the returned value is the vector - * of nodal values on that cell. You could - * then as well get the desired values - * through the @p{get_dof_values} function - * In the other case, when the cell has - * children, we use the restriction - * matrices provided by the finite element - * class to compute the interpolation - * from the children to the present cell. + * Return the result of the + * @p{neighbor_child_on_subface} + * function of the base class, + * but convert it so that one can + * also access the DoF data (the + * function in the base class + * only returns an iterator with + * access to the triangulation + * data). + */ + TriaIterator > + neighbor_child_on_subface (const unsigned int face_no, + const unsigned int subface_no) const; + + /** + * Return the interpolation of + * the given finite element + * function to the present + * cell. In the simplest case, + * the cell is a terminal one, + * i.e. has no children; then, + * the returned value is the + * vector of nodal values on that + * cell. You could then as well + * get the desired values through + * the @p{get_dof_values} + * function In the other case, + * when the cell has children, we + * use the restriction matrices + * provided by the finite element + * class to compute the + * interpolation from the + * children to the present cell. * - * It is assumed that both vectors already - * have the right size beforehand. This - * function assumes the existence of an - * interpolation from child cells to the - * mother cell, denoted by the restriction - * matrices of the finite element class. - * Futhermore, this interpolation should - * be possible for each child alone, i.e. - * it should be possible to compute the - * restriction by writing values obtained - * from each child directly into the output - * vector, without, for example, computing - * an average over all children. - * These properties, however, do not exist - * for all elements; an example is the - * DG(0) element, which represents - * piecewise constant elements: for these - * the restriction to mother cell could - * be the average of the values of the - * children, maybe weighted by the measure - * of each child. It is not yet decided - * what the this function does in these - * cases. + * It is assumed that both + * vectors already have the right + * size beforehand. This function + * assumes the existence of an + * interpolation from child cells + * to the mother cell, denoted by + * the restriction matrices of + * the finite element class. + * Futhermore, this interpolation + * should be possible for each + * child alone, i.e. it should + * be possible to compute the + * restriction by writing values + * obtained from each child + * directly into the output + * vector, without, for example, + * computing an average over all + * children. These properties, + * however, do not exist for all + * elements; an example is the + * DG(0) element, which + * represents piecewise constant + * elements: for these the + * restriction to mother cell + * could be the average of the + * values of the children, maybe + * weighted by the measure of + * each child. It is not yet + * decided what the this function + * does in these cases. * * Unlike the @p{get_dof_values} * function, this function is @@ -1187,57 +1219,73 @@ class DoFCellAccessor : public DoFObjectAccessor Vector &interpolated_values) const; /** - * This, again, is the counterpart to - * @p{get_interpolated_dof_values}: you - * specify the dof values on a cell and - * these are interpolated to the children - * of the present cell and set on the - * terminal cells. + * This, again, is the + * counterpart to + * @p{get_interpolated_dof_values}: + * you specify the dof values on + * a cell and these are + * interpolated to the children + * of the present cell and set on + * the terminal cells. * - * In principle, it works as follows: if - * the cell pointed to by this object is - * terminal, then the dof values are set - * in the global data vector by calling - * the @p{set_dof_values} function; - * otherwise, the values are prolonged - * to each of the children and this - * function is called for each of them. + * In principle, it works as + * follows: if the cell pointed + * to by this object is terminal, + * then the dof values are set in + * the global data vector by + * calling the @p{set_dof_values} + * function; otherwise, the + * values are prolonged to each + * of the children and this + * function is called for each of + * them. * - * Using the @p{get_interpolated_dof_values} - * and this function, you can compute the - * interpolation of a finite element - * function to a coarser grid by first - * getting the interpolated solution on a - * cell of the coarse grid and afterwards - * redistributing it using this function. + * Using the + * @p{get_interpolated_dof_values} + * and this function, you can + * compute the interpolation of a + * finite element function to a + * coarser grid by first getting + * the interpolated solution on a + * cell of the coarse grid and + * afterwards redistributing it + * using this function. * - * Note that for continuous finite - * elements, calling this function affects - * the dof values on neighboring cells as - * well. It may also violate continuity - * requirements for hanging nodes, if - * neighboring cells are less refined than - * the present one, or if their children are - * less refined than the children of this - * cell. These requirements - * are not taken care of and must be - * enforced by the user afterwards. + * Note that for continuous + * finite elements, calling this + * function affects the dof + * values on neighboring cells as + * well. It may also violate + * continuity requirements for + * hanging nodes, if neighboring + * cells are less refined than + * the present one, or if their + * children are less refined than + * the children of this + * cell. These requirements are + * not taken care of and must be + * enforced by the user + * afterwards. * - * It is assumed that both vectors already - * have the right size beforehand. This - * function relies on the existence of a - * natural interpolation property of - * finite element spaces of a cell to - * its children, denoted by the - * prolongation matrices of finite element - * classes. For some elements, the spaces - * on coarse and fine grids are not nested, - * in which case the interpolation to a - * child is not the identity; refer to the - * documentation of the respective finite - * element class for a description of what - * the prolongation matrices represent in - * this case. + * It is assumed that both + * vectors already have the right + * size beforehand. This function + * relies on the existence of a + * natural interpolation property + * of finite element spaces of a + * cell to its children, denoted + * by the prolongation matrices + * of finite element classes. For + * some elements, the spaces on + * coarse and fine grids are not + * nested, in which case the + * interpolation to a child is + * not the identity; refer to the + * documentation of the + * respective finite element + * class for a description of + * what the prolongation matrices + * represent in this case. * * Unlike the @p{set_dof_values} * function, this function is diff --git a/deal.II/deal.II/include/grid/tria_accessor.h b/deal.II/deal.II/include/grid/tria_accessor.h index 90e15c5789..664687eab9 100644 --- a/deal.II/deal.II/include/grid/tria_accessor.h +++ b/deal.II/deal.II/include/grid/tria_accessor.h @@ -84,11 +84,7 @@ class TriaAccessor TriaAccessor (const Triangulation *parent = 0, const int level = -1, const int index = -1, - const AccessorData * = 0) - : - present_level (level), - present_index (index), - tria (parent) {}; + const AccessorData * = 0); /** * Copy operator. Since this is @@ -297,10 +293,7 @@ class TriaObjectAccessor : public TriaAccessor TriaObjectAccessor (const Triangulation *parent = 0, const int level = -1, const int index = -1, - const AccessorData *local_data = 0) - : - TriaAccessor (parent, level, index, local_data) - {}; + const AccessorData *local_data = 0); /** * Copy the data of a line. Only @@ -815,10 +808,7 @@ class TriaObjectAccessor<1, dim> : public TriaAccessor TriaObjectAccessor (const Triangulation *parent = 0, const int level = -1, const int index = -1, - const AccessorData *local_data = 0) - : - TriaAccessor (parent, level, index, local_data) - {}; + const AccessorData *local_data = 0); /** * Copy the data of the given @@ -1252,10 +1242,7 @@ class TriaObjectAccessor<2, dim> : public TriaAccessor TriaObjectAccessor (const Triangulation *parent = 0, const int level = -1, const int index = -1, - const AccessorData *local_data = 0) - : - TriaAccessor (parent, level, index, local_data) - {}; + const AccessorData *local_data = 0); /** * Copy the data of the given quad. @@ -1742,10 +1729,7 @@ class TriaObjectAccessor<3, dim> : public TriaAccessor TriaObjectAccessor (const Triangulation *parent = 0, const int level = -1, const int index = -1, - const AccessorData *local_data = 0) - : - TriaAccessor (parent, level, index, local_data) - {}; + const AccessorData *local_data = 0); /** * Copy the data of the given @@ -2258,10 +2242,7 @@ class CellAccessor : public TriaObjectAccessor CellAccessor (const Triangulation *parent = 0, const int level = -1, const int index = -1, - const AccessorData *local_data = 0) - : - TriaObjectAccessor (parent, level, index, local_data) - {}; + const AccessorData *local_data = 0); /** * Return a pointer to the @@ -2440,6 +2421,48 @@ class CellAccessor : public TriaObjectAccessor */ TriaIterator > face (const unsigned int i) const; + + /** + * Return an iterator to that + * cell that neighbors the + * present cell on the given face + * and subface number. + * + * To succeed, the present cell + * must not be further refined, + * and the neighbor on the given + * face must be further refined + * exactly once; the returned + * cell is then a child of that + * neighbor. + * + * The function may not be called + * in 1d, since there we have no + * subfaces. The implementation + * of this function is rather + * straightforward in 2d, by + * first determining which face + * of the neighbor cell the + * present cell is bordering on + * (this is what the + * @p{neighbor_of_neighbor} + * function does), and then + * asking + * @p{GeometryInfo::child_cell_on_subface} + * for the index of the + * child. However, the function + * is more complicated in 3d, + * since there faces may have + * more than one orientation, and + * we have to use + * @p{face_orientation} for both + * this and the neighbor cell to + * figure out which cell we want + * to have. + */ + TriaIterator > + neighbor_child_on_subface (const unsigned int face_no, + const unsigned int subface_no) const; /** * Return the material id of this @@ -2563,6 +2586,19 @@ template <> bool CellAccessor<3>::point_inside (const Point<3> &) const; template <> bool CellAccessor<1>::has_boundary_lines () const; +template <> +TriaIterator<1,CellAccessor<1> > +CellAccessor<1>::neighbor_child_on_subface (const unsigned int, + const unsigned int) const; +template <> +TriaIterator<2,CellAccessor<2> > +CellAccessor<2>::neighbor_child_on_subface (const unsigned int, + const unsigned int) const; +template <> +TriaIterator<3,CellAccessor<3> > +CellAccessor<3>::neighbor_child_on_subface (const unsigned int, + const unsigned int) const; + template <> double TriaObjectAccessor<2, 2>::measure () const; template <> double TriaObjectAccessor<2, 3>::measure () const; template <> double TriaObjectAccessor<3, 3>::measure () const; diff --git a/deal.II/deal.II/include/grid/tria_accessor.templates.h b/deal.II/deal.II/include/grid/tria_accessor.templates.h index 7eb16fed5c..d70059e6b9 100644 --- a/deal.II/deal.II/include/grid/tria_accessor.templates.h +++ b/deal.II/deal.II/include/grid/tria_accessor.templates.h @@ -29,6 +29,19 @@ /*------------------------ Functions: TriaAccessor ---------------------------*/ +template +inline +TriaAccessor::TriaAccessor (const Triangulation *parent, + const int level, + const int index, + const AccessorData *) + : + present_level (level), + present_index (index), + tria (parent) +{} + + template inline @@ -114,6 +127,19 @@ TriaAccessor::get_triangulation () const /*------------------------ Functions: LineAccessor ---------------------------*/ +template +inline +TriaObjectAccessor<1,dim>:: +TriaObjectAccessor (const Triangulation *parent, + const int level, + const int index, + const AccessorData *local_data) + : + TriaAccessor (parent, level, index, local_data) +{} + + + template inline bool @@ -281,6 +307,19 @@ TriaObjectAccessor<1,dim>::operator -- () /*------------------------ Functions: QuadAccessor ---------------------------*/ +template +inline +TriaObjectAccessor<2,dim>:: +TriaObjectAccessor (const Triangulation *parent, + const int level, + const int index, + const AccessorData *local_data) + : + TriaAccessor (parent, level, index, local_data) +{} + + + template inline bool @@ -477,6 +516,19 @@ TriaObjectAccessor<2,dim>::operator -- () /*------------------------ Functions: HexAccessor ---------------------------*/ +template +inline +TriaObjectAccessor<3,dim>:: +TriaObjectAccessor (const Triangulation *parent, + const int level, + const int index, + const AccessorData *local_data) + : + TriaAccessor (parent, level, index, local_data) +{} + + + template inline bool @@ -1054,6 +1106,19 @@ TriaObjectAccessor::operator -- () /*------------------------ Functions: CellAccessor -----------------------*/ +template +inline +CellAccessor:: +CellAccessor (const Triangulation *parent, + const int level, + const int index, + const AccessorData *local_data) + : + TriaObjectAccessor (parent, level, index, local_data) +{} + + + template <> inline TriaIterator<1,TriaObjectAccessor<0, 1> > diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index b79f0fc196..952988c7d0 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -637,6 +637,24 @@ DoFCellAccessor<3>::face (const unsigned int i) const #endif + +template +inline +TriaIterator > +DoFCellAccessor:: +neighbor_child_on_subface (const unsigned int face, + const unsigned int subface) const +{ + const TriaIterator > q + = CellAccessor::neighbor_child_on_subface (face, subface); + return TriaIterator > (this->tria, + q->level (), + q->index (), + this->dof_handler); +} + + + template template void diff --git a/deal.II/deal.II/source/grid/tria_accessor.cc b/deal.II/deal.II/source/grid/tria_accessor.cc index 0b69a294e8..dfc23cb812 100644 --- a/deal.II/deal.II/source/grid/tria_accessor.cc +++ b/deal.II/deal.II/source/grid/tria_accessor.cc @@ -2122,6 +2122,93 @@ bool CellAccessor::has_boundary_lines () const #endif +#if deal_II_dimension == 1 + +template <> +inline +TriaIterator<1,CellAccessor<1> > +CellAccessor<1>:: +neighbor_child_on_subface (const unsigned int, + const unsigned int) const +{ + Assert (false, ExcNotImplemented()); +} + +#endif + +#if deal_II_dimension == 2 + +template <> +inline +TriaIterator<2,CellAccessor<2> > +CellAccessor<2>:: +neighbor_child_on_subface (const unsigned int face, + const unsigned int subface) const +{ + Assert (!this->has_children(), + ExcMessage ("The present cell must not have children!")); + Assert (!this->at_boundary(face), + ExcMessage ("The present cell must have a valid neighbor!")); + Assert (this->neighbor(face)->level() == this->level(), + ExcMessage ("The neighbor must be on the same level as this cell!")); + Assert (this->neighbor(face)->has_children() == true, + ExcMessage ("The neighbor must have children!")); + + const unsigned int neighbor_neighbor + = this->neighbor_of_neighbor (face); + const unsigned int neighbor_child_index + = GeometryInfo<2>::child_cell_on_face(neighbor_neighbor,subface); + + return this->neighbor(face)->child(neighbor_child_index); +} + +#endif + +#if deal_II_dimension == 3 + +template <> +inline +TriaIterator<3,CellAccessor<3> > +CellAccessor<3>:: +neighbor_child_on_subface (const unsigned int face, + const unsigned int subface) const +{ + Assert (!this->has_children(), + ExcMessage ("The present cell must not have children!")); + Assert (!this->at_boundary(face), + ExcMessage ("The present cell must have a valid neighbor!")); + Assert (this->neighbor(face)->level() == this->level(), + ExcMessage ("The neighbor must be on the same level as this cell!")); + Assert (this->neighbor(face)->has_children() == true, + ExcMessage ("The neighbor must have children!")); + + static const unsigned int subface_translation[4] + = { 0, 3, 2, 1 }; + // see whether face and + // the neighbor's + // counterface share the + // same indexing of + // children. if not so, + // translate child + // indices + const unsigned int neighbor_neighbor + = this->neighbor_of_neighbor (face); + const bool face_orientations_match + = (this->neighbor(face)->face_orientation(neighbor_neighbor) == + this->face_orientation(face)); + const unsigned int neighbor_child_index + = (GeometryInfo<3>:: + child_cell_on_face(neighbor_neighbor, + (face_orientations_match ? + subface : + subface_translation[subface]))); + + return this->neighbor(face)->child(neighbor_child_index); +} + +#endif + + // explicit instantiations template class TriaAccessor; diff --git a/tests/bits/Makefile b/tests/bits/Makefile index ec69c44d8c..e423ca1505 100644 --- a/tests/bits/Makefile +++ b/tests/bits/Makefile @@ -117,6 +117,7 @@ mesh_3d_11.exe : mesh_3d_11.g.$(OBJEXT) $(libraries) mesh_3d_12.exe : mesh_3d_12.g.$(OBJEXT) $(libraries) mesh_3d_13.exe : mesh_3d_13.g.$(OBJEXT) $(libraries) mesh_3d_14.exe : mesh_3d_14.g.$(OBJEXT) $(libraries) +mesh_3d_15.exe : mesh_3d_15.g.$(OBJEXT) $(libraries) q_points.exe : q_points.g.$(OBJEXT) $(libraries) @@ -149,7 +150,7 @@ tests = anna_1 anna_2 anna_3 anna_4 anna_5 anna_6 \ hyper_ball_3d cylinder coarsening_3d \ mesh_3d_1 mesh_3d_2 mesh_3d_3 mesh_3d_4 mesh_3d_5 \ mesh_3d_6 mesh_3d_7 mesh_3d_8 mesh_3d_9 mesh_3d_10 \ - mesh_3d_11 mesh_3d_12 mesh_3d_13 mesh_3d_14 \ + mesh_3d_11 mesh_3d_12 mesh_3d_13 mesh_3d_14 mesh_3d_15 \ q_points ############################################################ diff --git a/tests/bits/mesh_3d_15.cc b/tests/bits/mesh_3d_15.cc new file mode 100644 index 0000000000..44b0109a11 --- /dev/null +++ b/tests/bits/mesh_3d_15.cc @@ -0,0 +1,146 @@ +//---------------------------- mesh_3d_15.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003 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. +// +//---------------------------- mesh_3d_15.cc --------------------------- + + +// make sure TriaCellAccessor::neighbor_child_on_subface does what it +// is supposed to do. check it for dof accessors, since they simply +// call the tria accessors, and this way we catch both cases at the +// same time + +#include "mesh_3d.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + + +void check_this (Triangulation<3> &tria) +{ + FE_Q<3> fe(1); + + DoFHandler<3> dof_handler (tria); + dof_handler.distribute_dofs (fe); + + DoFHandler<3>::active_cell_iterator cell = dof_handler.begin_active(); + for (; cell!=dof_handler.end(); ++cell) + for (unsigned int face_no=0; face_no::faces_per_cell; + ++face_no) + if (!cell->at_boundary(face_no) + && + cell->neighbor(face_no)->has_children()) + for (unsigned int subface_no=0; + subface_no::subfaces_per_face; + ++subface_no) + { + // get an iterator + // pointing to the cell + // behind the present + // subface + const DoFHandler<3>::cell_iterator + neighbor = cell->neighbor(face_no); + const unsigned int neighbor_neighbor + = cell->neighbor_of_neighbor (face_no); + + // see whether face and + // the neighbor's + // counterface share the + // same indexing of + // children. if not so, + // translate child + // indices + const bool face_orientations_match + = (neighbor->face_orientation(neighbor_neighbor) == + cell->face_orientation(face_no)); + static const unsigned int subface_translation[4] + = { 0, 3, 2, 1 }; + const unsigned int neighbor_child_index + = (GeometryInfo<3>:: + child_cell_on_face(neighbor_neighbor, + (face_orientations_match ? + subface_no : + subface_translation[subface_no]))); + const DoFHandler<3>::active_cell_iterator neighbor_child + = neighbor->child(neighbor_child_index); + + Assert (neighbor_child == + cell->neighbor_child_on_subface (face_no, + subface_no), + ExcInternalError()); + } +} + + + +void check (Triangulation<3> &tria) +{ + (++tria.begin_active())->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + + deallog << "Initial check" << std::endl; + check_this (tria); + + for (unsigned int r=0; r<3; ++r) + { + tria.refine_global (1); + deallog << "Check " << r << std::endl; + check_this (tria); + } + + coarsen_global (tria); + deallog << "Check " << 1 << std::endl; + check_this (tria); + + tria.refine_global (1); + deallog << "Check " << 2 << std::endl; + check_this (tria); +} + + +int main () +{ + std::ofstream logfile("mesh_3d_15.output"); + deallog.attach(logfile); + deallog.depth_console(0); + + { + Triangulation<3> coarse_grid; + create_two_cubes (coarse_grid); + check (coarse_grid); + } + + { + Triangulation<3> coarse_grid; + create_L_shape (coarse_grid); + check (coarse_grid); + } + + { + Triangulation<3> coarse_grid; + GridGenerator::hyper_ball (coarse_grid); + check (coarse_grid); + } + +} + + + diff --git a/tests/results/i686-pc-linux-gnu+gcc2.95/bits/mesh_3d_15.output b/tests/results/i686-pc-linux-gnu+gcc2.95/bits/mesh_3d_15.output new file mode 100644 index 0000000000..0b20bdf25f --- /dev/null +++ b/tests/results/i686-pc-linux-gnu+gcc2.95/bits/mesh_3d_15.output @@ -0,0 +1,19 @@ + +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 diff --git a/tests/results/i686-pc-linux-gnu+gcc3.2/bits/mesh_3d_15.output b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/mesh_3d_15.output new file mode 100644 index 0000000000..0b20bdf25f --- /dev/null +++ b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/mesh_3d_15.output @@ -0,0 +1,19 @@ + +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 diff --git a/tests/results/sparc-sun-solaris2.7+gcc2.95/bits/mesh_3d_15.output b/tests/results/sparc-sun-solaris2.7+gcc2.95/bits/mesh_3d_15.output new file mode 100644 index 0000000000..0b20bdf25f --- /dev/null +++ b/tests/results/sparc-sun-solaris2.7+gcc2.95/bits/mesh_3d_15.output @@ -0,0 +1,19 @@ + +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Initial check +DEAL::Check 0 +DEAL::Check 1 +DEAL::Check 2 +DEAL::Check 1 +DEAL::Check 2 -- 2.39.5