From 90ca9698756dd3beada0dd85437ee5b9ef83b83e Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 7 Jan 2002 09:40:29 +0000 Subject: [PATCH] FilteredIterator class git-svn-id: https://svn.dealii.org/trunk@5346 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/grid/filtered_iterator.h | 1042 +++++++++++++++++ deal.II/doc/news/2001/c-3-2.html | 9 +- tests/deal.II/filtered_iterator.cc | 222 ++++ tests/deal.II/filtered_iterator.checked | 6 + 4 files changed, 1278 insertions(+), 1 deletion(-) create mode 100644 deal.II/deal.II/include/grid/filtered_iterator.h create mode 100644 tests/deal.II/filtered_iterator.cc create mode 100644 tests/deal.II/filtered_iterator.checked diff --git a/deal.II/deal.II/include/grid/filtered_iterator.h b/deal.II/deal.II/include/grid/filtered_iterator.h new file mode 100644 index 0000000000..3242582367 --- /dev/null +++ b/deal.II/deal.II/include/grid/filtered_iterator.h @@ -0,0 +1,1042 @@ +//------------------------ filtered_iterator.h ----------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2002 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. +// +//------------------------ filtered_iterator.h ----------------------- +#ifndef __deal2__filtered_iterator_h +#define __deal2__filtered_iterator_h + + + +/** + * In this namespace a number of classes is declared that may be used + * as filters in the @ref{FilteredIterator} class. The filters either + * check for binary information (for example, the @p{Active} filter + * class checks whether the object pointed to is active), or for + * valued information by comparison with prescribed values (for + * example, the @ref{LevelEqualTo} filter class checks whether the + * level of the object pointed to by the iterator under consideration + * is equal to a value that was given to the filter upon construction. + * + * For examples of use of these classes as well as requirements on + * filters see the general description of the @ref{FilteredIterator} + * class. + * + * @author Wolfgang Bangerth, 2002 + */ +namespace IteratorFilters +{ + /** + * Filter that evaluates to true if + * either the iterator points to an + * active object or an iterator + * past the end. + */ + class Active + { + public: + /** + * Evaluate the iterator and + * return true if the object is + * active or past the end. + */ + template + bool operator () (const Iterator &i) const; + }; + + /** + * Filter that evaluates to true if + * either the iterator points to an + * object for which the user flag + * is set or an iterator past the + * end. + */ + class UserFlagSet + { + public: + /** + * Evaluate the iterator and + * return true if the object + * has a set user flag or past + * the end. + */ + template + bool operator () (const Iterator &i) const; + }; + + + /** + * Filter that evaluates to true if + * either the iterator points to an + * object for which the user flag + * is not set or an iterator past + * the end. Inverse filter to the + * previous class. + */ + class UserFlagNotSet + { + public: + /** + * Evaluate the iterator and + * return true if the object + * has an unset user flag or + * past the end. + */ + template + bool operator () (const Iterator &i) const; + }; + + + /** + * Filter for iterators that + * evaluates to true if either the + * iterator is past the end or the + * level of the object pointed to + * is equal to a value given to the + * constructor. + */ + class LevelEqualTo + { + public: + /** + * Constructor. Store the level + * which iterators shall have + * to be evaluated to true. + */ + LevelEqualTo (const unsigned int level); + + /** + * Evaluation operator. Returns + * true if either the level of + * the object pointed to is + * equal to the stored value or + * the iterator is past the + * end. + */ + template + bool operator () (const Iterator &i) const; + + protected: + /** + * Stored value to compare the + * level with. + */ + const unsigned int level; + }; + + + + /** + * Filter for iterators that + * evaluates to true if either the + * iterator is past the end or the + * subdomain id of the object + * pointed to is equal to a value + * given to the constructor, + * assuming that the iterator + * allows querying for a subdomain + * id). + */ + class SubdomainEqualTo + { + public: + /** + * Constructor. Store the + * subdomain which iterators + * shall have to be evaluated + * to true. + */ + SubdomainEqualTo (const unsigned int subdomain_id); + + /** + * Evaluation operator. Returns + * true if either the subdomain + * of the object pointed to is + * equal to the stored value or + * the iterator is past the + * end. + */ + template + bool operator () (const Iterator &i) const; + + protected: + /** + * Stored value to compare the + * subdomain with. + */ + const unsigned int subdomain_id; + }; +}; + + +/** + * This class provides a certain view on a range of triangulation or + * DoFHandler iterators by only iterating over elements that satisfy a + * given filter (called a @em{predicate}, following the notation of + * the C++ standard library). Once initialized with a predicate and a + * value for the iterator, a filtered iterator hops to the next or + * previous element that satisfies the predicate if operators ++ or -- + * are invoked. Intermediate iterator values that lie in between but + * do not satisfy the predicate are skipped. It is thus very simple to + * write loops over a certain class of objects without the need to + * explicitely write down the condition they have to satisfy in each + * loop iteration. This in particular is helpful if functions are + * called with a pair of iterators denoting a range on which they + * shall act, by choosing a filtered iterator instead of usual ones. + * + * + * @sect3{Predicates} + * + * The object that represent the condition an iterator has to satisfy + * only have to provide an interface that allows to call the + * evaluation operator, i.e. @p{operator()}. This includes function + * pointers as well as classes that implement an @p{operator ()}. + * + * An example of a simple valid predicate is the following: given the function + * @begin{verbatim} + * template + * bool level_equal_to_3 (const Iterator c) + * { + * return (static_cast(c->level()) == 3); + * }; + * @end{verbatim} + * then + * @begin{verbatim} + * &level_equal_to_3::active_cell_iterator> + * @end{verbatim} + * is a valid predicate. + * + * Likewise, given the following binary function + * @begin{verbatim} + * template + * bool level_equal_to (const Iterator c, + * const unsigned int level) + * { + * return (static_cast(c->level()) == level); + * }; + * @end{verbatim} + * then + * @begin{verbatim} + * std::bind2nd (std::ptr_fun(&level_equal_to), 3) + * @end{verbatim} + * is another valid predicate (here: a function that returns true if + * either the iterator is past the end or the level is equal to the + * second argument; this second argument is bound to a fixed value + * using the @p{std::bind2nd} function). + * + * Finally, classes can be predicates. The following class is one: + * @begin{verbatim} + * class Active + * { + * public: + * template + * bool operator () (const Iterator &i) const { + * return (i->active()); + * } + * }; + * @end{verbatim} + * and objects of this type can be used as predicates. Likewise, this + * more complicated one can also be used: + * @begin{verbatim} + * class SubdomainEqualTo + * { + * public: + * SubdomainEqualTo (const unsigned int subdomain_id) + * : subdomain_id (subdomain_id) {}; + * + * template + * bool operator () (const Iterator &i) const { + * return (i->subdomain_id() == subdomain_id); + * } + * + * private: + * const unsigned int subdomain_id; + * }; + * @end{verbatim} + * Objects like @p{SubdomainEqualTo(3)} can then be used as predicates. + * + * Since whenever a predicate is evaluated it is checked that the + * iterator checked is actually valid (i.e. not past the end), no + * checks for this case have to be performed inside predicates. + * + * A number of filter classes are already implemented in the + * @ref{IteratorFilters} namespace, but writing different ones is + * simple following the examples above. + * + * + * @sect3{Initialization of filtered iterators} + * + * Filtered iterators are given a predicate at construction time which + * cannot be changed any more. This behaviour would be expected if the + * predicate would have been given as a template parameter to the + * class, but since that would make the declaration of filtered + * iterators a nightmare, we rather give the predicate as an + * unchangeable entity to the constructor. Note that one can assign a + * filtered iterator with one predicate to another filtered iterator + * with another type; yet, this does @em{not} change the predicate of + * the assigned-to iterator, only the pointer indicating the iterator + * is changed. + * + * If a filtered iterator is not assigned a value of the underlying + * (unfiltered) iterator type, the default value is taken. If, + * however, a value is given to the constructor, that value has either + * to be past the end, or has to satisfy the predicate. For example, + * if the predicate only evaluates to true if the level of an object + * is equal to three, then @p{tria.begin_active(3)} would be a valid + * choice while @p{tria.begin()} would not since the latter also + * returns iterators to non-active cells which always start at level + * 0. + * + * Since one often only has some iterator and wants to set a filtered + * iterator to the first one that satisfies a predicate (for example, + * the first one for which the user flag is set, or the first one with + * a given subdomain id), there are assignement functions + * @p{set_to_next_positive} and @p{set_to_previous_positive} that + * assign the next or last previous iterator that satisfies the + * predicate, i.e. they follow the list of iterators in either + * direction until they find a matching one (or the past-the-end + * iterator). Like the @p{operator=} they return the resulting value + * of the filtered iterator. + * + * + * @sect3{Examples} + * + * The following call counts the number of active cells that + * have a set user flag: + * @begin{verbatim} + * FilteredIterator::active_cell_iterator> + * begin (IteratorFilters::UserFlagSet()), + * end (IteratorFilters::UserFlagSet()); + * begin.set_to_next_positive(tria.begin_active()); + * end = tria.end(); + * n_flagged_cells = std::distance (begin, end); + * @begin{verbatim} + * Note that by the @p{set_to_next_positive} call the first cell with + * a set user flag was assigned to the @p{begin} iterator. For the + * @{end} iterator, no such call was necessary, since the past-the-end + * iterator always satisfies all predicates. + * + * The same can be achieved by the following snippet, though harder to read: + * @begin{verbatim} + * typedef FilteredIterator::active_cell_iterator> FI; + * n_flagged_cells = + * std::distance (FI(IteratorFilters::UserFlagSet()) + * .set_to_next_positive(tria.begin_active()), + * FI(IteratorFilters::UserFlagSet(), tria.end())); + * @begin{verbatim} + * It relies on the fact that if we create an unnamed filtered + * iterator with a given predicate but no iterator value and assign it + * the next positive value with respect to this predicate, it returns + * itself which is then used as the first parameter to the + * @p{std::distance} function. This procedure is not necessary for the + * end element to this function here, since the past-the-end iterator + * always satisfies the predicate so that we can assign this value to + * the filtered iterator directly in the constructor. + * + * Finally, the following loop only assembles the matrix on cells with + * subdomain id equal to three: + * @begin{verbatim} + * FilteredIterator::active_cell_iterator> + * cell (FilteredIterator::SubdomainEqualTo(3)), + * endc (FilteredIterator::SubdomainEqualTo(3), tria.end()); + * cell.set_to_next_positive (tria.begin_active()); + * for (; cell!=endc; ++cell) + * assemble_local_matrix (cell); + * @end{verbatim} + * + * Since comparison between filtered and unfiltered iterators is + * defined, we could as well have let the @p{endc} variable in the + * last example be of type + * @p{Triangulation::active_cell_iterator} since it is unchanged + * and its value does not depend on the filter. + * + * @author Wolfgang Bangerth, 2002 + */ +template +class FilteredIterator : public BaseIterator +{ + public: + /** + * Typedef to the accessor type + * of the underlying iterator. + */ + typedef typename BaseIterator::AccessorType AccessorType; + + /** + * Constructor. Set the iterator + * to the default state and use + * the given predicate for + * filtering subsequent + * assignement and iteration. + */ + template + FilteredIterator (Predicate p); + + /** + * Constructor. Use the given + * predicate for filtering and + * initialize the iterator with + * the given value. This initial + * value has to satisfy the + * predicate. + */ + template + FilteredIterator (Predicate p, + const BaseIterator &bi); + + /** + * Copy constructor. Copy the + * predicate and iterator value + * of the given argument. + */ + FilteredIterator (const FilteredIterator &fi); + + /** + * Destructor. + */ + ~FilteredIterator (); + + /** + * Assignment operator. Copy the + * iterator value of the + * argument, but as discussed in + * the class documentation, the + * predicate of the argument is + * not copied. The iterator value + * underlying the argument has to + * satisfy the predicate of the + * object assigned to, as given + * at its construction time. + */ + FilteredIterator & operator = (const FilteredIterator &fi); + + /** + * Assignment operator. Copy the + * iterator value of the + * argument, and keep the + * predicate of this object. The + * given iterator value has to + * satisfy the predicate of the + * object assigned to, as given + * at its construction time. + */ + FilteredIterator & operator = (const BaseIterator &fi); + + /** + * Search for the next iterator + * from @p{bi} onwards that + * satisfies the predicate of + * this object and assign it to + * this object. + * + * Since filtered iterators are + * automatically converted to the + * underlying base iterator type, + * you can also give a filtered + * iterator as argument to this + * function. + */ + FilteredIterator & + set_to_next_positive (const BaseIterator &bi); + + /** + * As above, but search for the + * previous iterator from @p{bi} + * backwards that satisfies the + * predicate of this object and + * assign it to this object. + * + * Since filtered iterators are + * automatically converted to the + * underlying base iterator type, + * you can also give a filtered + * iterator as argument to this + * function. + */ + FilteredIterator & + set_to_previous_positive (const BaseIterator &bi); + + /** + * Compare for equality of the + * underlying iterator values of + * this and the given object. + * + * We do not compare for equality + * of the predicates. + */ + bool operator == (const FilteredIterator &fi) const; + + /** + * Compare for equality of the + * underlying iterator value of + * this object with the given + * object. + * + * The predicate of this object + * is irrelevant for this + * operation. + */ + bool operator == (const BaseIterator &fi) const; + + /** + * Compare for inequality of the + * underlying iterator values of + * this and the given object. + * + * We do not compare for equality + * of the predicates. + */ + bool operator != (const FilteredIterator &fi) const; + + /** + * Compare for inequality of the + * underlying iterator value of + * this object with the given + * object. + * + * The predicate of this object + * is irrelevant for this + * operation. + */ + bool operator != (const BaseIterator &fi) const; + + /** + * Compare for ordering of the + * underlying iterator values of + * this and the given object. + * + * We do not compare the + * predicates. + */ + bool operator < (const FilteredIterator &fi) const; + + /** + * Compare for ordering of the + * underlying iterator value of + * this object with the given + * object. + * + * The predicate of this object + * is irrelevant for this + * operation. + */ + bool operator < (const BaseIterator &fi) const; + + /** + * Prefix advancement operator: + * move to the next iterator + * value satisfying the predicate + * and return the new iterator + * value. + */ + FilteredIterator & operator ++ (); + + /** + * Postfix advancement operator: + * move to the next iterator + * value satisfying the predicate + * and return the old iterator + * value. + */ + FilteredIterator operator ++ (int); + + /** + * Prefix decrement operator: + * move to the previous iterator + * value satisfying the predicate + * and return the new iterator + * value. + */ + FilteredIterator & operator -- (); + + /** + * Postfix advancement operator: + * move to the previous iterator + * value satisfying the predicate + * and return the old iterator + * value. + */ + FilteredIterator operator -- (int); + + /** + * Exception. + */ + DeclException1 (ExcInvalidElement, + BaseIterator, + << "The element " << arg1 + << " with which you want to compare or which you want to" + << " assign is invalid since it does not satisfy the predicate."); + + private: + + /** + * Base class to encapsulate a + * predicate object. Since + * predicates can be of different + * types and we do not want to + * code these types into the + * template parameter list of the + * filtered iterator class, we + * use a base class with an + * abstract function and + * templatized derived classes + * that implement the use of + * actual predicate types through + * the virtual function. + */ + class PredicateBase + { + public: + /** + * Mark the destructor + * virtual to allow + * destruction through + * pointers to the base + * class. + */ + virtual ~PredicateBase () {}; + + /** + * Abstract function which in + * derived classes denotes + * the evaluation of the + * predicate on the give + * iterator. + */ + virtual bool operator () (const BaseIterator &bi) const = 0; + + /** + * Generate a copy of this + * object, i.e. of the actual + * type of this pointer. + */ + virtual PredicateBase * clone () const = 0; + }; + + + /** + * Actual implementation of the + * above abstract base class. Use + * a template parameter to denote + * the actual type of the + * predicate and store a copy of + * it. When the virtual function + * is called evaluate the given + * iterator with the stored copy + * of the predicate. + */ + template + class PredicateTemplate : public PredicateBase + { + public: + /** + * Constructor. Take a + * predicate and store a copy + * of it. + */ + PredicateTemplate (const Predicate &predicate); + + /** + * Evaluate the iterator with + * the stored copy of the + * predicate. + */ + virtual bool operator () (const BaseIterator &bi) const; + + /** + * Generate a copy of this + * object, i.e. of the actual + * type of this pointer. + */ + virtual PredicateBase * clone () const; + + private: + /** + * Copy of the predicate. + */ + const Predicate predicate; + }; + + /** + * Pointer to an object that + * encapsulated the actual data + * type of the predicate given to + * the constructor. + */ + const PredicateBase * predicate; +}; + + + +/* ------------------ Inline functions and templates ------------ */ + + +template +template +inline +FilteredIterator:: +FilteredIterator (Predicate p) + : + predicate (new PredicateTemplate(p)) +{}; + + + +template +template +inline +FilteredIterator:: +FilteredIterator (Predicate p, + const BaseIterator &bi) + : + BaseIterator (bi), + predicate (new PredicateTemplate(p)) +{ + Assert ((state() != IteratorState::valid) || (*predicate) (*this), + ExcInvalidElement(bi)); +}; + + + +template +inline +FilteredIterator:: +FilteredIterator (const FilteredIterator &fi) + : + BaseIterator (static_cast(fi)), + predicate (fi.predicate->clone ()) +{}; + + + +template +inline +FilteredIterator:: +~FilteredIterator () +{ + delete predicate; + predicate = 0; +}; + + + +template +inline +FilteredIterator & +FilteredIterator:: +operator = (const FilteredIterator &fi) +{ + Assert ((fi.state() != IteratorState::valid) || (*predicate)(fi), + ExcInvalidElement(fi)); + BaseIterator::operator = (fi); + return *this; +}; + + + +template +inline +FilteredIterator & +FilteredIterator:: +operator = (const BaseIterator &bi) +{ + Assert ((bi.state() != IteratorState::valid) || (*predicate)(bi), + ExcInvalidElement(bi)); + BaseIterator::operator = (bi); + return *this; +}; + + + +template +inline +FilteredIterator & +FilteredIterator:: +set_to_next_positive (const BaseIterator &bi) +{ + BaseIterator::operator = (bi); + while ((state() == IteratorState::valid) && + ( ! (*predicate)(*this))) + BaseIterator::operator++ (); + + return *this; +}; + + + +template +inline +FilteredIterator & +FilteredIterator:: +set_to_previous_positive (const BaseIterator &bi) +{ + BaseIterator::operator = (bi); + while ((state() == IteratorState::valid) && + ( ! predicate(*this))) + BaseIterator::operator-- (); + + return *this; +}; + + + +template +inline +bool +FilteredIterator:: +operator == (const FilteredIterator &fi) const +{ + return (static_cast(*this) + == + static_cast(fi)); +}; + + + +template +inline +bool +FilteredIterator:: +operator != (const FilteredIterator &fi) const +{ + return (static_cast(*this) + != + static_cast(fi)); +}; + + + +template +inline +bool +FilteredIterator:: +operator < (const FilteredIterator &fi) const +{ + return (static_cast(*this) + < + static_cast(fi)); +}; + + + + +template +inline +bool +FilteredIterator:: +operator == (const BaseIterator &bi) const +{ + return (static_cast(*this) == bi); +}; + + + +template +inline +bool +FilteredIterator:: +operator != (const BaseIterator &bi) const +{ + return (static_cast(*this) != bi); +}; + + + +template +inline +bool +FilteredIterator:: +operator < (const BaseIterator &bi) const +{ + return (static_cast(*this) < bi); +}; + + +template +inline +FilteredIterator & +FilteredIterator:: +operator ++ () +{ + if (state() == IteratorState::valid) + do + BaseIterator::operator++ (); + while ((state() == IteratorState::valid) && + !(*predicate) (*this)); + return *this; +}; + + + +template +inline +FilteredIterator +FilteredIterator:: +operator ++ (int) +{ + const FilteredIterator old_state = *this; + + if (state() == IteratorState::valid) + do + BaseIterator::operator++ (); + while ((state() == IteratorState::valid) && + !(*predicate) (*this)); + return old_state; +}; + + + + +template +inline +FilteredIterator & +FilteredIterator:: +operator -- () +{ + if (state() == IteratorState::valid) + do + BaseIterator::operator-- (); + while ((state() == IteratorState::valid) && + !(*predicate) (*this)); + return *this; +}; + + + +template +inline +FilteredIterator +FilteredIterator:: +operator -- (int) +{ + const FilteredIterator old_state = *this; + + if (state() == IteratorState::valid) + do + BaseIterator::operator-- (); + while ((state() == IteratorState::valid) && + !(*predicate) (*this)); + return old_state; +}; + + + +template +template +inline +FilteredIterator::PredicateTemplate:: +PredicateTemplate (const Predicate &predicate) + : + predicate (predicate) +{}; + + + +template +template +bool +FilteredIterator::PredicateTemplate:: +operator () (const BaseIterator &bi) const +{ + return predicate(bi); +}; + + + +template +template +FilteredIterator::PredicateBase * +FilteredIterator::PredicateTemplate:: +clone () const +{ + return new PredicateTemplate (predicate); +}; + + + +namespace IteratorFilters +{ + +// ---------------- IteratorFilters::Active --------- + + template + inline + bool + Active::operator () (const Iterator &i) const + { + return (i->active()); + }; + + +// ---------------- IteratorFilters::UserFlagSet --------- + + template + inline + bool + UserFlagSet::operator () (const Iterator &i) const + { + return (i->user_flag_set()); + }; + + +// ---------------- IteratorFilters::UserFlagNotSet --------- + + template + inline + bool + UserFlagNotSet::operator () (const Iterator &i) const + { + return (! i->user_flag_set()); + }; + + +// ---------------- IteratorFilters::LevelEqualTo --------- + inline + LevelEqualTo::LevelEqualTo (const unsigned int level) + : + level (level) + {}; + + + + template + inline + bool + LevelEqualTo::operator () (const Iterator &i) const + { + return (static_cast(i->level()) == level); + }; + + + +// ---------------- IteratorFilters::SubdomainEqualTo --------- + inline + SubdomainEqualTo::SubdomainEqualTo (const unsigned int subdomain_id) + : + subdomain_id (subdomain_id) + {}; + + + + template + inline + bool + SubdomainEqualTo::operator () (const Iterator &i) const + { + return (static_cast(i->subdomain_id()) == subdomain_id); + }; +}; + + +/*------------------------- filtered_iterator.h ------------------------*/ +#endif +/*------------------------- filtered_iterator.h ------------------------*/ + + diff --git a/deal.II/doc/news/2001/c-3-2.html b/deal.II/doc/news/2001/c-3-2.html index 2645479070..9c668991a6 100644 --- a/deal.II/doc/news/2001/c-3-2.html +++ b/deal.II/doc/news/2001/c-3-2.html @@ -259,6 +259,14 @@ documentation, etc.

deal.II

    +
  1. + New: The FilteredIterator class + provides a view on ranges of iterators by iterating over only + those objects that satisfy a certain predicate. +
    + (GK 2002/01/07) +

    +
  2. Improved: It is now possible to read in unconnected domains through the GridIn class, since @@ -286,7 +294,6 @@ documentation, etc. (GK 2001/12/07)

    -
  3. New: FiniteElement::has_support_on_face allows to check diff --git a/tests/deal.II/filtered_iterator.cc b/tests/deal.II/filtered_iterator.cc new file mode 100644 index 0000000000..5f7ff76bc7 --- /dev/null +++ b/tests/deal.II/filtered_iterator.cc @@ -0,0 +1,222 @@ +//------------------------- filtered_iterator.cc ------------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2001, 2002 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. +// +//------------------------- filtered_iterator.cc ------------------------ + +// check filtered iterators + + +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + + +std::ofstream logfile("filtered_iterator.output"); + + +DeclException2 (ExcNumberMismatch, + int, int, + << "The numbers " << arg1 << " and " << arg2 + << " should be equation, but are not."); + + + +typedef Triangulation<2>::active_cell_iterator active_cell_iterator; + +template +bool level_equal_to_3 (const Iterator c) +{ + return (static_cast(c->level()) == 3); +}; + + + +template +bool level_equal_to (const Iterator c, + const unsigned int level) +{ + return (static_cast(c->level()) == level); +}; + + +void test () +{ + Triangulation<2> tria; + GridGenerator::hyper_cube(tria, -1, 1); + tria.refine_global (1); + tria.begin_active()->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + tria.refine_global (2); + + // we now have a number of cells, + // flag them with some subdomain + // ids based on their position, in + // particular we take the quadrant + // (octant) + active_cell_iterator cell = tria.begin_active (), + endc = tria.end (); + for (; cell!=endc; ++cell) + { + unsigned int subdomain = 0; + for (unsigned int d=0; d<2; ++d) + if (cell->center()(d) > 0) + subdomain |= (1<set_subdomain_id (subdomain); + }; + + + // check 1: count number of cells + // on some level + if (true) + { + const IteratorFilters::LevelEqualTo predicate(3); + FilteredIterator + begin (predicate), + end (predicate, tria.end()); + begin.set_to_next_positive (tria.begin_active ()); + + Assert (std::distance (begin, end) == + static_cast(tria.n_active_cells (3)), + ExcInternalError()); + logfile << "Check 1: " + << (std::distance (begin, end) == + static_cast(tria.n_active_cells (3)) + ? + "OK" : "Failed") + << std::endl; + }; + + + // check 2: count number of cells + // on some level in a different way + if (true) + { + FilteredIterator + begin (&level_equal_to_3, + tria.begin_active (3)), + end (&level_equal_to_3, + tria.end()); + + Assert (std::distance (begin, end) == + static_cast(tria.n_active_cells (3)), + ExcInternalError()); + logfile << "Check 2: " + << (std::distance (begin, end) == + static_cast(tria.n_active_cells (3)) + ? + "OK" : "Failed") + << std::endl; + }; + + + // check 3: count number of cells + // on some level in yet a different + // way + if (true) + { + FilteredIterator + begin (std::bind2nd (std::ptr_fun(&level_equal_to), 3), + tria.begin_active (3)), + end (std::bind2nd (std::ptr_fun(&level_equal_to), 3), + tria.end()); + + Assert (std::distance (begin, end) == + static_cast(tria.n_active_cells (3)), + ExcInternalError()); + logfile << "Check 3: " + << (std::distance (begin, end) == + static_cast(tria.n_active_cells (3)) + ? + "OK" : "Failed") + << std::endl; + }; + + + // check 4: and yet another possibility + if (true) + { + typedef FilteredIterator FI; + + Assert (std::distance (FI(std::bind2nd (std::ptr_fun(&level_equal_to), 3)).set_to_next_positive(tria.begin_active()), + FI(std::bind2nd (std::ptr_fun(&level_equal_to), 3), tria.end())) == + static_cast(tria.n_active_cells (3)), + ExcInternalError()); + logfile << "Check 4: " + << (std::distance (FI(std::bind2nd (std::ptr_fun(&level_equal_to), 3)).set_to_next_positive(tria.begin_active()), + FI(std::bind2nd (std::ptr_fun(&level_equal_to), 3), tria.end())) == + static_cast(tria.n_active_cells (3)) + ? + "OK" : "Failed") + << std::endl; + }; + + + // check 5: check that we loop over + // all cells with a given subdomain + // id + if (true) + { + typedef FilteredIterator FI; + const IteratorFilters::SubdomainEqualTo predicate(1); + FI cell (predicate); + cell.set_to_next_positive (tria.begin_active()); + active_cell_iterator endc (tria.end()); + active_cell_iterator cell1 = tria.begin_active (); + + while (cell1->subdomain_id () != 1) + ++cell1; + + while (true) + { + // move filtered iterator ahead + ++cell; + // move unfiltered iterator + // ahead + ++cell1; + while ((cell1 != endc) && + (cell1->subdomain_id () != 1)) + ++cell1; + + Assert (cell == cell1, ExcInternalError()); + Assert (cell1 == cell, ExcInternalError()); + + if (cell.state() != IteratorState::valid) + break; + }; + Assert (cell == endc, ExcInternalError()); + Assert (cell1 == endc, ExcInternalError()); + + logfile << "Check 5: OK" << std::endl; + }; +}; + + +int main () +{ + logfile.precision(4); + deallog.attach(logfile); + deallog.depth_console(0); + + test (); + + return 0; +}; + diff --git a/tests/deal.II/filtered_iterator.checked b/tests/deal.II/filtered_iterator.checked new file mode 100644 index 0000000000..81c0478bbd --- /dev/null +++ b/tests/deal.II/filtered_iterator.checked @@ -0,0 +1,6 @@ + +Check 1: OK +Check 2: OK +Check 3: OK +Check 4: OK +Check 5: OK -- 2.39.5