From: Denis Davydov Date: Mon, 31 Aug 2015 14:05:08 +0000 (+0200) Subject: introduced FECollection::find_least_face_dominating_fe(const std::set&) const --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 2365ee1d5f..ac681fb290 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -240,6 +240,13 @@ inconvenience this causes.
    +
  1. New: introduced hp::FECollection::find_least_face_dominating_fe(const std::set &fes) + which aims to find the least dominating finite element w.r.t. those provided + as fe_indices in @p fes. +
    + (Denis Davydov, Wolfgang Bangerth, 2015/08/31) +
    2. +
  2. New: Introduce an option for FE_Nothing to dominate any other FE. Therefore at interfaces where, for example, a Q1 meets an FE_Nothing, we will force the traces of the two functions to be the same. Because the diff --git a/include/deal.II/hp/fe_collection.h b/include/deal.II/hp/fe_collection.h index 40f85fb6d3..8433cb308e 100644 --- a/include/deal.II/hp/fe_collection.h +++ b/include/deal.II/hp/fe_collection.h @@ -187,6 +187,31 @@ namespace hp */ bool hp_constraints_are_implemented () const; + /** + * Try to find a least dominant finite element inside this FECollection + * which dominates other finite elements provided as fe_indices in @p fes . + * For example, if FECollection consists of {Q1,Q2,Q3,Q4} and we are looking + * for the least dominant FE for Q3 and Q4 (@p fes is {2,3}), then the + * answer is Q3 and therefore this function will return its index in + * FECollection, namely 2. + * + * If we were not able to find a finite element, the function returns + * numbers::invalid_unsigned_int . + * + * Note that for the cases like when FECollection consists of + * {FE_Nothing x FE_Nothing, Q1xQ2, Q2xQ1} with @p fes = {1}, the function + * will not find the most dominating element as the default behavior of + * FE_Nothing is to return FiniteElementDomination::no_requirements when + * comparing for face domination. This, therefore, can't be considered as a + * dominating element in the sense defined in FiniteElementDomination . + * + * Finally, for the purpose of this function we consider that an element + * dominates itself. Thus if FECollection contains {Q1,Q2,Q4,Q3} and @p fes + * = {3}, the function returns 3. + */ + unsigned int + find_least_face_dominating_fe (const std::set &fes) const; + /** * Return a component mask with as many elements as this object has vector * components and of which exactly the one component is true that diff --git a/source/hp/fe_collection.cc b/source/hp/fe_collection.cc index 6ddd321596..b4161d04eb 100644 --- a/source/hp/fe_collection.cc +++ b/source/hp/fe_collection.cc @@ -21,6 +21,70 @@ DEAL_II_NAMESPACE_OPEN namespace hp { + template + unsigned int + FECollection::find_least_face_dominating_fe (const std::set &fes) const + { + // If the set of elements to be dominated contains only a single element X, + // then by definition the dominating set contains this single element X + // (because each element can dominate itself). There may also be others, + // say Y1...YN. Next you have to find one or more elements in the dominating + // set {X,Y1...YN} that is the weakest. Well, you can't find one that is + // weaker than X because if it were, it would not dominate X. In other words, + // X is guaranteed to be in the subset of {X,Y1...YN} of weakest dominating + // elements. Since we only guarantee that the function returns one of them, + // we may as well return X right away. + if (fes.size()==1) + return *fes.begin(); + + const hp::FECollection &fe_collection = *this; + std::set candidate_fes; + + // first loop over all FEs and check which can dominate those given in FEs: + for (unsigned int cur_fe = 0; cur_fe < fe_collection.size(); cur_fe++) + { + FiniteElementDomination::Domination domination = FiniteElementDomination::no_requirements; + // check if cur_fe can dominate all FEs in @p fes: + for (std::set::const_iterator it = fes.begin(); + it!=fes.end(); ++it) + { + Assert (*it < fe_collection.size(), + ExcIndexRangeType (*it, 0, fe_collection.size())); + domination = domination & + fe_collection[cur_fe].compare_for_face_domination + (fe_collection[*it]); + } + + // if we found dominating element, keep them in a set. + if (domination == FiniteElementDomination::this_element_dominates) + candidate_fes.insert(cur_fe); + } + + // among the ones we found, pick one that is dominated by all others and + // thus should represent the largest FE space. + if (candidate_fes.size() == 1) + { + return *candidate_fes.begin(); + } + else + for (std::set::const_iterator it = candidate_fes.begin(); it!=candidate_fes.end(); ++it) + { + FiniteElementDomination::Domination domination = FiniteElementDomination::no_requirements; + for (std::set::const_iterator ito = candidate_fes.begin(); ito!=candidate_fes.end(); ++ito) + if (it != ito) + { + domination = domination & + fe_collection[*it].compare_for_face_domination(fe_collection[*ito]); + } + + if (domination == FiniteElementDomination::other_element_dominates || + domination == FiniteElementDomination::either_element_can_dominate /*covers cases like candidate_fes={Q1,Q1}*/) + return *it; + } + // We couldn't find the FE, return invalid_unsigned_int : + return numbers::invalid_unsigned_int; + } + template FECollection::FECollection () {} diff --git a/tests/hp/fe_collection_05.cc b/tests/hp/fe_collection_05.cc new file mode 100644 index 0000000000..278b36f1ba --- /dev/null +++ b/tests/hp/fe_collection_05.cc @@ -0,0 +1,179 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2005 - 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + + +// test the results of FECollection::find_least_face_dominating_fe(), namely for: +// {Q1,Q2,Q3,Q4} with {2,3} => Q3 2 +// {Q1xQ1, Q2xQ2, Q3xQ4, Q4xQ3} with {2,3} => Q2xQ2 1 +// {Q1xQ1, Q3xQ4, Q4xQ3} with {1,2} => Q1xQ1 0 +// {0x0, 0x0, Q1x0, 0xQ1} with {2,3} => none invalid_unsigned_int +// {0x0, 0x0, Q1x0, 0xQ1} with {2,3} => 0x0 0 (with dominating FE_Nothing) +// {Q1xQ1,Q1xQ1,Q2xQ1,Q1,Q2} with {2,3} => Q1 0 +// {Q4xQ4, Q5xQ5, Q3xQ4, Q4xQ3} with {2,3} => none invalid_unsigned_int +// {Q1,Q2,Q4,Q3} with {3} => Q3 3 + + +#include "../tests.h" +#include +#include +#include +#include +#include + +#include + + +template +void test () +{ + std::set fes; + fes.insert(2); + fes.insert(3); + + // {Q1,Q2,Q3,Q4} + { + hp::FECollection fe_collection; + fe_collection.push_back (FE_Q(1)); + fe_collection.push_back (FE_Q(2)); + fe_collection.push_back (FE_Q(3)); + fe_collection.push_back (FE_Q(4)); + deallog << fe_collection.find_least_face_dominating_fe(fes) << std::endl; + } + + // {Q1xQ1, Q2xQ2, Q3xQ4, Q4xQ3} + { + hp::FECollection fe_collection; + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Q(1),1)); + fe_collection.push_back (FESystem(FE_Q(2),1, + FE_Q(2),1)); + fe_collection.push_back (FESystem(FE_Q(3),1, + FE_Q(4),1)); + fe_collection.push_back (FESystem(FE_Q(4),1, + FE_Q(3),1)); + deallog << fe_collection.find_least_face_dominating_fe(fes) << std::endl; + } + + // {Q1xQ1, Q3xQ4, Q4xQ3} + { + hp::FECollection fe_collection; + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Q(1),1)); + fe_collection.push_back (FESystem(FE_Q(3),1, + FE_Q(4),1)); + fe_collection.push_back (FESystem(FE_Q(4),1, + FE_Q(3),1)); + std::set fes; + fes.insert(1); + fes.insert(2); + deallog << fe_collection.find_least_face_dominating_fe(fes) << std::endl; + } + + // {0x0, 0x0, Q1x0, 0xQ1} + { + hp::FECollection fe_collection; + fe_collection.push_back (FESystem(FE_Nothing(),1, + FE_Nothing(),1)); + fe_collection.push_back (FESystem(FE_Nothing(),1, + FE_Nothing(),1)); + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Nothing(),1)); + fe_collection.push_back (FESystem(FE_Nothing(),1, + FE_Q(1),1)); + const unsigned int ind = fe_collection.find_least_face_dominating_fe(fes); + if (ind == numbers::invalid_unsigned_int) + deallog << "numbers::invalid_unsigned_int" << std::endl; + else + deallog << ind << std::endl; + } + + // dominating FE_Nothing + // {0x0, 0x0, Q1x0, 0xQ1} + { + hp::FECollection fe_collection; + fe_collection.push_back (FESystem(FE_Nothing(1,true),1, + FE_Nothing(1,true),1)); + fe_collection.push_back (FESystem(FE_Nothing(1,true),1, + FE_Nothing(1,true),1)); + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Nothing(1,true),1)); + fe_collection.push_back (FESystem(FE_Nothing(1,true),1, + FE_Q(1),1)); + deallog << fe_collection.find_least_face_dominating_fe(fes) << std::endl; + } + + + // {Q1xQ1,Q1xQ1,Q2xQ1,Q1,Q2} + { + hp::FECollection fe_collection; + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Q(1),1)); + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Q(1),1)); + fe_collection.push_back (FESystem(FE_Q(2),1, + FE_Q(1),1)); + fe_collection.push_back (FESystem(FE_Q(1),1, + FE_Q(2),1)); + deallog << fe_collection.find_least_face_dominating_fe(fes) << std::endl; + } + + // {Q4xQ4, Q5xQ5, Q3xQ4, Q4xQ3 + { + hp::FECollection fe_collection; + fe_collection.push_back (FESystem(FE_Q(4),1, + FE_Q(4),1)); + fe_collection.push_back (FESystem(FE_Q(5),1, + FE_Q(5),1)); + fe_collection.push_back (FESystem(FE_Q(3),1, + FE_Q(4),1)); + fe_collection.push_back (FESystem(FE_Q(4),1, + FE_Q(3),1)); + const unsigned int ind = fe_collection.find_least_face_dominating_fe(fes); + if (ind == numbers::invalid_unsigned_int) + deallog << "numbers::invalid_unsigned_int" << std::endl; + else + deallog << ind << std::endl; + } + + // {Q1,Q2,Q4,Q3} + { + hp::FECollection fe_collection; + fe_collection.push_back (FE_Q(1)); + fe_collection.push_back (FE_Q(2)); + fe_collection.push_back (FE_Q(4)); + fe_collection.push_back (FE_Q(3)); + std::set fes; + fes.insert(3); + deallog << fe_collection.find_least_face_dominating_fe(fes) << std::endl; + } + + +} + +int main () +{ + std::ofstream logfile("output"); + logfile.precision(2); + + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test<2> (); + test<3> (); + + deallog << "OK" << std::endl; +} diff --git a/tests/hp/fe_collection_05.output b/tests/hp/fe_collection_05.output new file mode 100644 index 0000000000..66c6056a96 --- /dev/null +++ b/tests/hp/fe_collection_05.output @@ -0,0 +1,18 @@ + +DEAL::2 +DEAL::1 +DEAL::0 +DEAL::numbers::invalid_unsigned_int +DEAL::0 +DEAL::0 +DEAL::numbers::invalid_unsigned_int +DEAL::3 +DEAL::2 +DEAL::1 +DEAL::0 +DEAL::numbers::invalid_unsigned_int +DEAL::0 +DEAL::0 +DEAL::numbers::invalid_unsigned_int +DEAL::3 +DEAL::OK