From 57f36131cb4d0dee36c120958c164a65a68e26e9 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Wed, 30 Dec 2020 21:49:02 -0700 Subject: [PATCH] Added tests for FENothing constraints. --- tests/hp/fe_nothing_dominates.h | 92 +++++++++++++++++++++++++ tests/hp/fe_nothing_dominates_01.cc | 84 ++++++++++++++++++++++ tests/hp/fe_nothing_dominates_01.output | 75 ++++++++++++++++++++ tests/hp/fe_nothing_dominates_02.cc | 80 +++++++++++++++++++++ tests/hp/fe_nothing_dominates_02.output | 37 ++++++++++ tests/hp/fe_nothing_dominates_03.cc | 79 +++++++++++++++++++++ tests/hp/fe_nothing_dominates_03.output | 49 +++++++++++++ 7 files changed, 496 insertions(+) create mode 100644 tests/hp/fe_nothing_dominates.h create mode 100644 tests/hp/fe_nothing_dominates_01.cc create mode 100644 tests/hp/fe_nothing_dominates_01.output create mode 100644 tests/hp/fe_nothing_dominates_02.cc create mode 100644 tests/hp/fe_nothing_dominates_02.output create mode 100644 tests/hp/fe_nothing_dominates_03.cc create mode 100644 tests/hp/fe_nothing_dominates_03.output diff --git a/tests/hp/fe_nothing_dominates.h b/tests/hp/fe_nothing_dominates.h new file mode 100644 index 0000000000..1836c0bac3 --- /dev/null +++ b/tests/hp/fe_nothing_dominates.h @@ -0,0 +1,92 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// Check for continuity requirements of two neighboring finite elements +// and project a given function subject to these constraints. + + +#include + +#include +#include + +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +#include "../test_grids.h" + + +template +void +project(const hp::FECollection &fe_collection, + const hp::QCollection & q_collection, + const Function & function) +{ +#ifdef DEBUG + Assert(fe_collection.size() == 2, ExcInternalError()); + Assert(q_collection.size() == 2, ExcInternalError()); + for (unsigned int f = 0; f < fe_collection.size(); ++f) + Assert(fe_collection[f].n_components() == function.n_components, + ExcInternalError()); +#endif + + // setup + // +---+---+ + // | 0 | 1 | + // +---+---+ + Triangulation tria; + TestGrids::hyper_line(tria, 2); + + DoFHandler dofh(tria); + (++(dofh.begin_active()))->set_active_fe_index(1); + dofh.distribute_dofs(fe_collection); + + // make constraints + AffineConstraints constraints; + DoFTools::make_hanging_node_constraints(dofh, constraints); + constraints.close(); + deallog << "constraints:" << std::endl; + constraints.print(deallog.get_file_stream()); + + // project function with constraints + Vector solution(dofh.n_dofs()); + VectorTools::project(dofh, constraints, q_collection, function, solution); + + // verify output + deallog << "dof values:" << std::endl; + Vector cell_values; + for (const auto &cell : dofh.active_cell_iterators()) + { + cell_values.reinit(cell->get_fe().n_dofs_per_cell()); + cell->get_dof_values(solution, cell_values); + + deallog << " cell " << cell->active_cell_index() << ":"; + for (const auto &value : cell_values) + deallog << " " << value; + deallog << std::endl; + } + + deallog << "OK" << std::endl; +} diff --git a/tests/hp/fe_nothing_dominates_01.cc b/tests/hp/fe_nothing_dominates_01.cc new file mode 100644 index 0000000000..8ddb354370 --- /dev/null +++ b/tests/hp/fe_nothing_dominates_01.cc @@ -0,0 +1,84 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// Check for continuity requirements for neighboring FE_Nothing and FE_Q +// elements. +// +// Note: output corresponds to current results, which are not correct. + + +#include +#include +#include + +#include +#include + +#include +#include + +#include "../tests.h" + +#include "fe_nothing_dominates.h" + + +template +void +test(const bool fe_nothing_dominates) +{ + Functions::ConstantFunction function(1.); + + { + deallog << "(FE_Nothing, FE_Q)" << std::endl; + hp::FECollection fe_collection(FE_Nothing(/*n_components=*/1, + fe_nothing_dominates), + FE_Q(1)); + hp::QCollection q_collection(Quadrature(1), QGauss(2)); + project(fe_collection, q_collection, function); + } + + { + deallog << "(FE_Q, FE_Nothing)" << std::endl; + hp::FECollection fe_collection(FE_Q(1), + FE_Nothing(/*n_components=*/1, + fe_nothing_dominates)); + hp::QCollection q_collection(QGauss(2), Quadrature(1)); + project(fe_collection, q_collection, function); + } +} + + +int +main() +{ + initlog(); + + deallog << std::boolalpha; + for (const bool fe_nothing_dominates : {false, true}) + { + deallog << "FE_Nothing dominates: " << fe_nothing_dominates << std::endl; + deallog.push("1d"); + test<1>(fe_nothing_dominates); + deallog.pop(); + deallog.push("2d"); + test<2>(fe_nothing_dominates); + deallog.pop(); + deallog.push("3d"); + test<3>(fe_nothing_dominates); + deallog.pop(); + } +} diff --git a/tests/hp/fe_nothing_dominates_01.output b/tests/hp/fe_nothing_dominates_01.output new file mode 100644 index 0000000000..8428e72878 --- /dev/null +++ b/tests/hp/fe_nothing_dominates_01.output @@ -0,0 +1,75 @@ + +DEAL::FE_Nothing dominates: false +DEAL:1d::(FE_Nothing, FE_Q) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: +DEAL:1d:: cell 1: 1 1 +DEAL:1d::OK +DEAL:1d::(FE_Q, FE_Nothing) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: 1 1 +DEAL:1d:: cell 1: +DEAL:1d::OK +DEAL:2d::(FE_Nothing, FE_Q) +DEAL:2d::constraints: +DEAL:2d::dof values: +DEAL:2d:: cell 0: +DEAL:2d:: cell 1: 1 1 1 1 +DEAL:2d::OK +DEAL:2d::(FE_Q, FE_Nothing) +DEAL:2d::constraints: +DEAL:2d::dof values: +DEAL:2d:: cell 0: 1 1 1 1 +DEAL:2d:: cell 1: +DEAL:2d::OK +DEAL:3d::(FE_Nothing, FE_Q) +DEAL:3d::constraints: +DEAL:3d::dof values: +DEAL:3d:: cell 0: +DEAL:3d:: cell 1: 1 1 1 1 1 1 1 1 +DEAL:3d::OK +DEAL:3d::(FE_Q, FE_Nothing) +DEAL:3d::constraints: +DEAL:3d::dof values: +DEAL:3d:: cell 0: 1 1 1 1 1 1 1 1 +DEAL:3d:: cell 1: +DEAL:3d::OK +DEAL::FE_Nothing dominates: true +DEAL:1d::(FE_Nothing, FE_Q) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: +DEAL:1d:: cell 1: 1 1 +DEAL:1d::OK +DEAL:1d::(FE_Q, FE_Nothing) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: 1 1 +DEAL:1d:: cell 1: +DEAL:1d::OK +DEAL:2d::(FE_Nothing, FE_Q) +DEAL:2d::constraints: +DEAL:2d::dof values: +DEAL:2d:: cell 0: +DEAL:2d:: cell 1: 1 1 1 1 +DEAL:2d::OK +DEAL:2d::(FE_Q, FE_Nothing) +DEAL:2d::constraints: +DEAL:2d::dof values: +DEAL:2d:: cell 0: 1 1 1 1 +DEAL:2d:: cell 1: +DEAL:2d::OK +DEAL:3d::(FE_Nothing, FE_Q) +DEAL:3d::constraints: +DEAL:3d::dof values: +DEAL:3d:: cell 0: +DEAL:3d:: cell 1: 1 1 1 1 1 1 1 1 +DEAL:3d::OK +DEAL:3d::(FE_Q, FE_Nothing) +DEAL:3d::constraints: +DEAL:3d::dof values: +DEAL:3d:: cell 0: 1 1 1 1 1 1 1 1 +DEAL:3d:: cell 1: +DEAL:3d::OK diff --git a/tests/hp/fe_nothing_dominates_02.cc b/tests/hp/fe_nothing_dominates_02.cc new file mode 100644 index 0000000000..843b20c302 --- /dev/null +++ b/tests/hp/fe_nothing_dominates_02.cc @@ -0,0 +1,80 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// Check for continuity requirements for neighboring +// (FE_NothingxFE_Nothing) and (FE_QxFE_Q) elements. +// The twist: only one FE_Nothing elements dominates. +// +// Note: output corresponds to current results, which are not correct. + + +#include +#include +#include + +#include +#include +#include + +#include +#include + +#include "../tests.h" + +#include "fe_nothing_dominates.h" + + +template +void +test() +{ + Functions::ConstantFunction function(1., 2); + + { + deallog << "(FE_QxFE_Q, FE_Nothing(true)xFE_Nothing(false))" << std::endl; + hp::FECollection fe_collection( + FESystem(FE_Q(1), FE_Q(1)), + FESystem(FE_Nothing(1, true), FE_Nothing(1, false))); + hp::QCollection q_collection(QGauss(2), Quadrature(1)); + project(fe_collection, q_collection, function); + } + { + deallog << "(FE_Nothing(true)xFE_Nothing(false), FE_QxFE_Q)" << std::endl; + hp::FECollection fe_collection( + FESystem(FE_Nothing(1, true), FE_Nothing(1, false)), + FESystem(FE_Q(1), FE_Q(1))); + hp::QCollection q_collection(Quadrature(1), QGauss(2)); + project(fe_collection, q_collection, function); + } +} + + +int +main() +{ + initlog(); + + deallog.push("1d"); + test<1>(); + deallog.pop(); + deallog.push("2d"); + test<2>(); + deallog.pop(); + deallog.push("3d"); + test<3>(); + deallog.pop(); +} diff --git a/tests/hp/fe_nothing_dominates_02.output b/tests/hp/fe_nothing_dominates_02.output new file mode 100644 index 0000000000..7ee5d38428 --- /dev/null +++ b/tests/hp/fe_nothing_dominates_02.output @@ -0,0 +1,37 @@ + +DEAL:1d::(FE_QxFE_Q, FE_Nothing(true)xFE_Nothing(false)) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: 1 1 1 1 +DEAL:1d:: cell 1: +DEAL:1d::OK +DEAL:1d::(FE_Nothing(true)xFE_Nothing(false), FE_QxFE_Q) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: +DEAL:1d:: cell 1: 1 1 1 1 +DEAL:1d::OK +DEAL:2d::(FE_QxFE_Q, FE_Nothing(true)xFE_Nothing(false)) +DEAL:2d::constraints: +DEAL:2d::dof values: +DEAL:2d:: cell 0: 1 1 1 1 1 1 1 1 +DEAL:2d:: cell 1: +DEAL:2d::OK +DEAL:2d::(FE_Nothing(true)xFE_Nothing(false), FE_QxFE_Q) +DEAL:2d::constraints: +DEAL:2d::dof values: +DEAL:2d:: cell 0: +DEAL:2d:: cell 1: 1 1 1 1 1 1 1 1 +DEAL:2d::OK +DEAL:3d::(FE_QxFE_Q, FE_Nothing(true)xFE_Nothing(false)) +DEAL:3d::constraints: +DEAL:3d::dof values: +DEAL:3d:: cell 0: 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 +DEAL:3d:: cell 1: +DEAL:3d::OK +DEAL:3d::(FE_Nothing(true)xFE_Nothing(false), FE_QxFE_Q) +DEAL:3d::constraints: +DEAL:3d::dof values: +DEAL:3d:: cell 0: +DEAL:3d:: cell 1: 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 +DEAL:3d::OK diff --git a/tests/hp/fe_nothing_dominates_03.cc b/tests/hp/fe_nothing_dominates_03.cc new file mode 100644 index 0000000000..6e37f67671 --- /dev/null +++ b/tests/hp/fe_nothing_dominates_03.cc @@ -0,0 +1,79 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// Check for continuity requirements for neighboring +// (FE_QxFE_Nothing) and (FE_NothingxFE_Q) elements. +// The twist: only one FE_Nothing element dominates. +// +// Note: output corresponds to current results, which are not correct. + + +#include +#include +#include + +#include +#include +#include + +#include +#include + +#include "../tests.h" + +#include "fe_nothing_dominates.h" + + +template +void +test() +{ + Functions::ConstantFunction function(1., 2); + hp::QCollection q_collection(QGauss(2), QGauss(2)); + + { + deallog << "(FE_Nothing(true)xFE_Q, FE_QxFE_Nothing(false))" << std::endl; + hp::FECollection fe_collection( + FESystem(FE_Nothing(1, true), FE_Q(1)), + FESystem(FE_Q(1), FE_Nothing(1, false))); + project(fe_collection, q_collection, function); + } + { + deallog << "(FE_QxFE_Nothing(false), FE_Nothing(true)xFE_Q)" << std::endl; + hp::FECollection fe_collection( + FESystem(FE_Q(1), FE_Nothing(1, false)), + FESystem(FE_Nothing(1, true), FE_Q(1))); + project(fe_collection, q_collection, function); + } +} + + +int +main() +{ + initlog(); + + deallog.push("1d"); + test<1>(); + deallog.pop(); + deallog.push("2d"); + test<2>(); + deallog.pop(); + deallog.push("3d"); + test<3>(); + deallog.pop(); +} diff --git a/tests/hp/fe_nothing_dominates_03.output b/tests/hp/fe_nothing_dominates_03.output new file mode 100644 index 0000000000..545c03e92c --- /dev/null +++ b/tests/hp/fe_nothing_dominates_03.output @@ -0,0 +1,49 @@ + +DEAL:1d::(FE_Nothing(true)xFE_Q, FE_QxFE_Nothing(false)) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: 1 1 +DEAL:1d:: cell 1: 1 1 +DEAL:1d::OK +DEAL:1d::(FE_QxFE_Nothing(false), FE_Nothing(true)xFE_Q) +DEAL:1d::constraints: +DEAL:1d::dof values: +DEAL:1d:: cell 0: 1 1 +DEAL:1d:: cell 1: 1 1 +DEAL:1d::OK +DEAL:2d::(FE_Nothing(true)xFE_Q, FE_QxFE_Nothing(false)) +DEAL:2d::constraints: + 4 = 0 + 6 = 0 +DEAL:2d::dof values: +DEAL:2d:: cell 0: 1 1 1 1 +DEAL:2d:: cell 1: 0 1.5 0 1.5 +DEAL:2d::OK +DEAL:2d::(FE_QxFE_Nothing(false), FE_Nothing(true)xFE_Q) +DEAL:2d::constraints: + 1 = 0 + 3 = 0 +DEAL:2d::dof values: +DEAL:2d:: cell 0: 1.5 0 1.5 0 +DEAL:2d:: cell 1: 1 1 1 1 +DEAL:2d::OK +DEAL:3d::(FE_Nothing(true)xFE_Q, FE_QxFE_Nothing(false)) +DEAL:3d::constraints: + 8 = 0 + 10 = 0 + 12 = 0 + 14 = 0 +DEAL:3d::dof values: +DEAL:3d:: cell 0: 1 1 1 1 1 1 1 1 +DEAL:3d:: cell 1: 0 1.5 0 1.5 0 1.5 0 1.5 +DEAL:3d::OK +DEAL:3d::(FE_QxFE_Nothing(false), FE_Nothing(true)xFE_Q) +DEAL:3d::constraints: + 1 = 0 + 3 = 0 + 5 = 0 + 7 = 0 +DEAL:3d::dof values: +DEAL:3d:: cell 0: 1.5 0 1.5 0 1.5 0 1.5 0 +DEAL:3d:: cell 1: 1 1 1 1 1 1 1 1 +DEAL:3d::OK -- 2.39.5