From 98ea1c8f2213dbd7e7eee6637bdf3e0296e4dee8 Mon Sep 17 00:00:00 2001 From: kanschat Date: Wed, 12 Sep 2012 15:26:51 +0000 Subject: [PATCH] remove specializations in TriaObjects git-svn-id: https://svn.dealii.org/trunk@26313 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/grid/tria_object.h | 2 +- deal.II/include/deal.II/grid/tria_objects.h | 104 ++++--- deal.II/source/grid/tria.cc | 21 +- deal.II/source/grid/tria_objects.cc | 302 +++----------------- deal.II/source/grid/tria_objects.inst.in | 26 +- 5 files changed, 111 insertions(+), 344 deletions(-) diff --git a/deal.II/include/deal.II/grid/tria_object.h b/deal.II/include/deal.II/grid/tria_object.h index 8fa11c17ae..2cd849306e 100644 --- a/deal.II/include/deal.II/grid/tria_object.h +++ b/deal.II/include/deal.II/grid/tria_object.h @@ -26,7 +26,7 @@ namespace internal /** * Class template for the structdim-dimensional cells constituting - * a Triangulation of dimension structdim or lower dimensional + * a dealii::Triangulation of dimension structdim or lower dimensional * objects of higher dimensions. They are characterized by the * (global) indices of their faces, which are cells of dimension * structdim-1 or vertices if structdim=1. diff --git a/deal.II/include/deal.II/grid/tria_objects.h b/deal.II/include/deal.II/grid/tria_objects.h index f64e7fd589..b9c60b5567 100644 --- a/deal.II/include/deal.II/grid/tria_objects.h +++ b/deal.II/include/deal.II/grid/tria_objects.h @@ -21,13 +21,12 @@ DEAL_II_NAMESPACE_OPEN -//TODO: See if we can unify the class hierarchy a bit using partial -// specialization of classes here, e.g. declare a general class -// TriaObjects. -// To consider for this: in that case we would have to duplicate quite a few -// things, e.g. TriaObjectsHex is derived from TriaObjects and -// declares mainly additional data. This would have to be changed in case of a -// partial specialization. +//TODO: This should all be cleaned up. Currently, only a single +//function in the library makes use of the odd specializations, and +//this function is Triangulation::execute_refinement() in 3D. I +//assume, that the other refinement functions would profit from using +//next_free_single_object() and next_free_pair_object, but they seem +//to get around it. //TODO: The TriaObjects class contains a std::vector. This is only an //efficient storage scheme if G is relatively well packed, i.e. it's not a @@ -35,7 +34,8 @@ DEAL_II_NAMESPACE_OPEN //actually the case. template class Triangulation; - +template class TriaRawIterator; +template class TriaAccessor; namespace internal { @@ -51,7 +51,7 @@ namespace internal * Objects of these classes are included in the TriaLevel and TriaFaces * classes. * - * @author Tobias Leicht, Guido Kanschat, 2006, 2007 + * @author Tobias Leicht, Guido Kanschat, 2006, 2007, 2012 */ template @@ -124,7 +124,7 @@ namespace internal * already been processed. * * You can clear all used flags using - * Triangulation::clear_user_flags(). + * dealii::Triangulation::clear_user_flags(). */ std::vector user_flags; @@ -221,48 +221,52 @@ namespace internal /** * Return an iterator to the * next free slot for a - * single line. Only - * implemented for - * G=TriaObject<1> - * . + * single object. This + * function is only used by + * dealii::Triangulation::execute_refinement() + * in 3D. + * + * @warning Interestingly, + * this function is not used + * for 1D or 2D + * triangulations, where it + * seems the authors of the + * refinement function insist + * on reimplementing its + * contents. + * + * @todo This function is + * not instantiated for the + * codim-one case */ template - typename dealii::Triangulation::raw_line_iterator - next_free_single_line (const dealii::Triangulation &tria); + dealii::TriaRawIterator > + next_free_single_object (const dealii::Triangulation &tria); /** * Return an iterator to the * next free slot for a pair - * of lines. Only implemented - * for G=TriaObject<1> - * . - */ - template - typename dealii::Triangulation::raw_line_iterator - next_free_pair_line (const dealii::Triangulation &tria); - - /** - * Return an iterator to the - * next free slot for a - * single quad. Only - * implemented for - * G=TriaObject@<2@> - * . + * of objects. This + * function is only used by + * dealii::Triangulation::execute_refinement() + * in 3D. + * + * @warning Interestingly, + * this function is not used + * for 1D or 2D + * triangulations, where it + * seems the authors of the + * refinement function insist + * on reimplementing its + * contents. + * + * @todo This function is + * not instantiated for the + * codim-one case */ template - typename dealii::Triangulation::raw_quad_iterator - next_free_single_quad (const dealii::Triangulation &tria); - - /** - * Return an iterator to the - * next free slot for a pair - * of quads. Only implemented - * for G=TriaObject@<2@> - * . - */ - template - typename dealii::Triangulation::raw_quad_iterator - next_free_pair_quad (const dealii::Triangulation &tria); + dealii::TriaRawIterator > + next_free_pair_object (const dealii::Triangulation &tria); /** * Return an iterator to the @@ -401,7 +405,7 @@ namespace internal "but you can only ask for " << arg2 <<"_iterators."); /** - * Triangulation objects can + * dealii::Triangulation objects can * either access a user * pointer or a user * index. What you tried to @@ -930,16 +934,6 @@ namespace internal // declaration of explicit specializations - template<> - void - TriaObjects >::reserve_space (const unsigned int new_lines_in_pairs, - const unsigned int new_lines_single); - - template<> - void - TriaObjects >::reserve_space (const unsigned int new_quads_in_pairs, - const unsigned int new_quads_single); - template<> void TriaObjects >::monitor_memory (const unsigned int) const; diff --git a/deal.II/source/grid/tria.cc b/deal.II/source/grid/tria.cc index 3743949c8d..0e0f18ad22 100644 --- a/deal.II/source/grid/tria.cc +++ b/deal.II/source/grid/tria.cc @@ -5374,7 +5374,7 @@ namespace internal // up the two child lines // (++ takes care of the // end of the vector) - next_unused_line=triangulation.faces->lines.next_free_pair_line(triangulation); + next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation); Assert(next_unused_line.state() == IteratorState::valid, ExcInternalError()); @@ -5506,7 +5506,7 @@ namespace internal // the quad typename Triangulation::raw_line_iterator new_line; - new_line=triangulation.faces->lines.next_free_single_line(triangulation); + new_line=triangulation.faces->lines.next_free_single_object(triangulation); Assert (new_line->used() == false, ExcCellShouldBeUnused()); @@ -5561,7 +5561,7 @@ namespace internal // created quads. typename Triangulation::raw_quad_iterator new_quads[2]; - next_unused_quad=triangulation.faces->quads.next_free_pair_quad(triangulation); + next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation); new_quads[0] = next_unused_quad; Assert (new_quads[0]->used() == false, ExcCellShouldBeUnused()); @@ -5691,7 +5691,7 @@ namespace internal // all info typename Triangulation::raw_line_iterator new_child[2]; - new_child[0]=new_child[1]=triangulation.faces->lines.next_free_pair_line(triangulation); + new_child[0]=new_child[1]=triangulation.faces->lines.next_free_pair_object(triangulation); ++new_child[1]; new_child[0]->set_used_flag(); @@ -5960,7 +5960,7 @@ namespace internal // now search a slot for the two // child lines - next_unused_line=triangulation.faces->lines.next_free_pair_line(triangulation); + next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation); // set the child // pointer of the present @@ -6096,7 +6096,7 @@ namespace internal // anisotropically and the // two lines end up as // children of new line - next_unused_line=triangulation.faces->lines.next_free_pair_line(triangulation); + next_unused_line=triangulation.faces->lines.next_free_pair_object(triangulation); new_lines[i] = next_unused_line; ++next_unused_line; @@ -6193,7 +6193,7 @@ namespace internal // created quads. typename Triangulation::raw_quad_iterator new_quads[4]; - next_unused_quad=triangulation.faces->quads.next_free_pair_quad(triangulation); + next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation); new_quads[0] = next_unused_quad; Assert (new_quads[0]->used() == false, ExcCellShouldBeUnused()); @@ -6202,7 +6202,7 @@ namespace internal new_quads[1] = next_unused_quad; Assert (new_quads[1]->used() == false, ExcCellShouldBeUnused()); - next_unused_quad=triangulation.faces->quads.next_free_pair_quad(triangulation); + next_unused_quad=triangulation.faces->quads.next_free_pair_object(triangulation); new_quads[2] = next_unused_quad; Assert (new_quads[2]->used() == false, ExcCellShouldBeUnused()); @@ -6366,7 +6366,7 @@ namespace internal new_lines(n_new_lines); for (unsigned int i=0; ilines.next_free_single_line(triangulation); + new_lines[i] = triangulation.faces->lines.next_free_single_object(triangulation); Assert (new_lines[i]->used() == false, ExcCellShouldBeUnused()); @@ -6385,7 +6385,7 @@ namespace internal new_quads(n_new_quads); for (unsigned int i=0; iquads.next_free_single_quad(triangulation); + new_quads[i] = triangulation.faces->quads.next_free_single_object(triangulation); Assert (new_quads[i]->used() == false, ExcCellShouldBeUnused()); @@ -6412,7 +6412,6 @@ namespace internal { if (i%2==0) next_unused_hex=triangulation.levels[level+1]->cells.next_free_hex(triangulation,level+1); - else ++next_unused_hex; diff --git a/deal.II/source/grid/tria_objects.cc b/deal.II/source/grid/tria_objects.cc index 7ce3ccfd69..ff30938e7e 100644 --- a/deal.II/source/grid/tria_objects.cc +++ b/deal.II/source/grid/tria_objects.cc @@ -29,105 +29,10 @@ namespace internal { namespace Triangulation { - template<> - void - TriaObjects >::reserve_space (const unsigned int new_lines_in_pairs, - const unsigned int new_lines_single) - { - Assert(new_lines_in_pairs%2==0, ExcInternalError()); - - next_free_single=0; - next_free_pair=0; - reverse_order_next_free_single=false; - - // count the number of lines, of - // unused single lines and of - // unused pairs of lines - unsigned int n_lines=0; - unsigned int n_unused_pairs=0; - unsigned int n_unused_singles=0; - for (unsigned int i=0; i0) - new_size+=additional_single_lines; - - // only allocate space if necessary - if (new_size>cells.size()) - { - cells.reserve (new_size); - cells.insert (cells.end(), - new_size-cells.size(), - TriaObject<1> ()); - - used.reserve (new_size); - used.insert (used.end(), - new_size-used.size(), - false); - - user_flags.reserve (new_size); - user_flags.insert (user_flags.end(), - new_size-user_flags.size(), - false); - - children.reserve (new_size); - children.insert (children.end(), - new_size-children.size(), - -1); - - boundary_or_material_id.reserve (new_size); - boundary_or_material_id.insert (boundary_or_material_id.end(), - new_size-boundary_or_material_id.size(), - BoundaryOrMaterialId()); - - user_data.reserve (new_size); - user_data.insert (user_data.end(), - new_size-user_data.size(), - UserData()); - } - - if (n_unused_singles==0) - { - next_free_single=new_size-1; - reverse_order_next_free_single=true; - } - } - - - template <> + template template - typename dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_single_line (const dealii::Triangulation &tria) + dealii::TriaRawIterator > + TriaObjects::next_free_single_object (const dealii::Triangulation &tria) { // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this. @@ -166,19 +71,19 @@ namespace internal if (pos>0) next_free_single=pos-1; else - // no valid single line anymore - return tria.end_line(); + // no valid single object anymore + return dealii::TriaRawIterator >(&tria, -1, -1); } - return typename dealii::Triangulation::raw_line_iterator(&tria,0,pos); + return dealii::TriaRawIterator >(&tria, 0, pos); } - template <> + template template - typename dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_pair_line (const dealii::Triangulation &tria) + dealii::TriaRawIterator > + TriaObjects::next_free_pair_object (const dealii::Triangulation &tria) { // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this. @@ -194,58 +99,35 @@ namespace internal } if (pos>=last) // no free slot - return tria.end_line(); + return dealii::TriaRawIterator >(&tria, -1, -1); else next_free_pair=pos+2; - return typename dealii::Triangulation::raw_line_iterator(&tria,0,pos); - } - - - - template <> - template - typename dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_single_quad (const dealii::Triangulation &tria) - { - Assert(false, ExcWrongIterator("quad","line")); - return tria.end_quad(); - } - - - - template <> - template - typename dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_pair_quad (const dealii::Triangulation &tria) - { - Assert(false, ExcWrongIterator("quad","line")); - return tria.end_quad(); + return dealii::TriaRawIterator >(&tria, 0, pos); } - - template<> + template void - TriaObjects >::reserve_space (const unsigned int new_quads_in_pairs, - const unsigned int new_quads_single) + TriaObjects::reserve_space (const unsigned int new_objects_in_pairs, + const unsigned int new_objects_single) { - Assert(new_quads_in_pairs%2==0, ExcInternalError()); + Assert(new_objects_in_pairs%2==0, ExcInternalError()); next_free_single=0; next_free_pair=0; reverse_order_next_free_single=false; - // count the number of lines, of - // unused single lines and of - // unused pairs of lines - unsigned int n_quads=0; + // count the number of objects, of + // unused single objects and of + // unused pairs of objects + unsigned int n_objects=0; unsigned int n_unused_pairs=0; unsigned int n_unused_singles=0; for (unsigned int i=0; i0) - new_size+=additional_single_quads; + used.size() + new_objects_in_pairs - 2*n_unused_pairs; + if (additional_single_objects>0) + new_size+=additional_single_objects; // only allocate space if necessary if (new_size>cells.size()) @@ -284,7 +166,7 @@ namespace internal cells.reserve (new_size); cells.insert (cells.end(), new_size-cells.size(), - TriaObject<2> ()); + G ()); used.reserve (new_size); used.insert (used.end(), @@ -296,16 +178,19 @@ namespace internal new_size-user_flags.size(), false); - children.reserve (2*new_size); + const unsigned int factor = GeometryInfo::max_children_per_cell / 2; + children.reserve (factor*new_size); children.insert (children.end(), - 2*new_size-children.size(), + factor*new_size-children.size(), -1); - refinement_cases.reserve (new_size); - refinement_cases.insert (refinement_cases.end(), - new_size - refinement_cases.size(), - RefinementCase<2>::no_refinement); - + if (G::dimension > 1) + { + refinement_cases.reserve (new_size); + refinement_cases.insert (refinement_cases.end(), + new_size - refinement_cases.size(), + RefinementCase::no_refinement); + } boundary_or_material_id.reserve (new_size); boundary_or_material_id.insert (boundary_or_material_id.end(), @@ -326,109 +211,6 @@ namespace internal } - - - template <> - template - typename dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_single_quad (const dealii::Triangulation &tria) - { - // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this. - - int pos=next_free_single, - last=used.size()-1; - if (!reverse_order_next_free_single) - { - // first sweep forward, only use - // really single slots, do not use - // pair slots - for (; pos=last) - { - reverse_order_next_free_single=true; - next_free_single=used.size()-1; - pos=used.size()-1; - } - else - next_free_single=pos+1; - } - - if(reverse_order_next_free_single) - { - // second sweep, use all slots, even - // in pairs - for(;pos>=0;--pos) - if (!used[pos]) - break; - if (pos>0) - next_free_single=pos-1; - else - // no valid single quad anymore - return tria.end_quad(); - } - - return typename dealii::Triangulation::raw_quad_iterator(&tria,0,pos); - } - - - - template <> - template - typename dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_pair_quad (const dealii::Triangulation &tria) - { - // TODO: Think of a way to ensure that we are using the correct triangulation, i.e. the one containing *this. - - int pos=next_free_pair, - last=used.size()-1; - for (; pos=last) - // no free slot - return tria.end_quad(); - else - next_free_pair=pos+2; - - return typename dealii::Triangulation::raw_quad_iterator(&tria,0,pos); - } - - - template <> - template - typename dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_single_line (const dealii::Triangulation &tria) - { - Assert(false, ExcWrongIterator("line","quad")); - return tria.end_line(); - } - - - - template <> - template - typename dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_pair_line (const dealii::Triangulation &tria) - { - Assert(false, ExcWrongIterator("line","quad")); - return tria.end_line(); - } - - - - template <> template typename dealii::Triangulation::raw_hex_iterator @@ -536,9 +318,9 @@ namespace internal next_free_pair=0; reverse_order_next_free_single=false; - // count the number of lines, of unused - // single lines and of unused pairs of - // lines + // count the number of objects, of unused + // single objects and of unused pairs of + // objects unsigned int n_quads=0; unsigned int n_unused_pairs=0; unsigned int n_unused_singles=0; diff --git a/deal.II/source/grid/tria_objects.inst.in b/deal.II/source/grid/tria_objects.inst.in index 76fc1259ad..777e0ba6cb 100644 --- a/deal.II/source/grid/tria_objects.inst.in +++ b/deal.II/source/grid/tria_objects.inst.in @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2006, 2007, 2010 by the deal.II authors +// Copyright (C) 2006, 2007, 2010, 2012 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -14,23 +14,15 @@ for (deal_II_dimension : DIMENSIONS) { #if deal_II_dimension > 1 - template dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_single_line(const dealii::Triangulation &); - template dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_pair_line(const dealii::Triangulation &); - template dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_single_quad(const dealii::Triangulation &); - template dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_pair_quad(const dealii::Triangulation &); - template dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_single_line(const dealii::Triangulation &); - template dealii::Triangulation::raw_line_iterator - TriaObjects >::next_free_pair_line(const dealii::Triangulation &); - template dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_single_quad(const dealii::Triangulation &); - template dealii::Triangulation::raw_quad_iterator - TriaObjects >::next_free_pair_quad(const dealii::Triangulation &); + template dealii::TriaRawIterator > + TriaObjects >::next_free_single_object (const dealii::Triangulation &tria); + template dealii::TriaRawIterator > + TriaObjects >::next_free_pair_object (const dealii::Triangulation &tria); + template dealii::TriaRawIterator > + TriaObjects >::next_free_single_object (const dealii::Triangulation &tria); + template dealii::TriaRawIterator > + TriaObjects >::next_free_pair_object (const dealii::Triangulation &tria); #endif #if deal_II_dimension == 3 template dealii::Triangulation::raw_hex_iterator -- 2.39.5