From 52b214cf9c5d99edf1ca58d54eeb63366cc5a318 Mon Sep 17 00:00:00 2001 From: Magdalena Schreter Date: Mon, 19 Feb 2024 19:14:25 +0100 Subject: [PATCH] fix make_flux_sparsity_pattern for FENothing --- .../changes/minor/20240219SchreterRitthaler | 3 + .../lac/affine_constraints.templates.h | 7 +- tests/hp/flux_sparsity_03.cc | 144 ++++++++++++++++++ tests/hp/flux_sparsity_03.output | 15 ++ 4 files changed, 168 insertions(+), 1 deletion(-) create mode 100644 doc/news/changes/minor/20240219SchreterRitthaler create mode 100644 tests/hp/flux_sparsity_03.cc create mode 100644 tests/hp/flux_sparsity_03.output diff --git a/doc/news/changes/minor/20240219SchreterRitthaler b/doc/news/changes/minor/20240219SchreterRitthaler new file mode 100644 index 0000000000..bf1f8986b0 --- /dev/null +++ b/doc/news/changes/minor/20240219SchreterRitthaler @@ -0,0 +1,3 @@ +Fixed: Fix make_flux_sparsity_pattern for a FECollection containing FENothing elements. +
+(Magdalena Schreter-Fleischhacker, Andreas Ritthaler, 2024/02/19) diff --git a/include/deal.II/lac/affine_constraints.templates.h b/include/deal.II/lac/affine_constraints.templates.h index c6d24ec496..359196623f 100644 --- a/include/deal.II/lac/affine_constraints.templates.h +++ b/include/deal.II/lac/affine_constraints.templates.h @@ -4714,12 +4714,17 @@ AffineConstraints::add_entries_local_to_global( const size_type n_local_rows = row_indices.size(); const size_type n_local_cols = col_indices.size(); + // Early return if the length of row and column indices is zero, relevant + // for the usage with FENothing. + if (n_local_cols == 0 && n_local_rows == 0) + return; + typename internal::AffineConstraints::ScratchDataAccessor scratch_data(this->scratch_data); std::vector> &cell_entries = scratch_data->new_entries; cell_entries.resize(0); - cell_entries.reserve(row_indices.size() * col_indices.size()); + cell_entries.reserve(n_local_rows * n_local_cols); // if constrained entries should be kept, need to add rows and columns of // those to the sparsity pattern diff --git a/tests/hp/flux_sparsity_03.cc b/tests/hp/flux_sparsity_03.cc new file mode 100644 index 0000000000..73e3028448 --- /dev/null +++ b/tests/hp/flux_sparsity_03.cc @@ -0,0 +1,144 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 - 2020 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. +// +// --------------------------------------------------------------------- + +// A test that checks make_flux_sparsity_pattern +// with a FECollection containing FE_Q and FE_Nothing elements. + +// x__________x__________ _________ +// | | | | +// | FE_Q | FE_N | FE_N | x: constrained DoFs +// | | | | +// |__________|__________|_________| +// x x + +#include + +#include + +#include +#include + +#include +#include +#include +#include + +#include + +#include +#include + +#include +#include + +#include "../tests.h" + +template +void +check() +{ + // create triangulation + Triangulation triangulation; + std::vector subdivisions(dim, 1U); + subdivisions[0] = 3; + subdivisions[1] = 1; + if (dim == 3) + subdivisions[2] = 1; + Point p1, p2; + switch (dim) + { + case 2: + p1[0] = p1[1] = 0.0; + p2[0] = 3.0; + p2[1] = 1.0; + break; + case 3: + p1[0] = p1[1] = p1[2] = 0.0; + p2[0] = 3.0; + p2[1] = p2[2] = 1.0; + break; + default: + DEAL_II_NOT_IMPLEMENTED(); + } + + GridGenerator::subdivided_hyper_rectangle(triangulation, + subdivisions, + p1, + p2); + + // create FE Collection and insert two FE objects + // Q1 and one FE object Nothing + hp::FECollection fe_collection; + + fe_collection.push_back(FE_Q(1)); + fe_collection.push_back(FE_Nothing()); + + // set-up DoFHandler + DoFHandler dof_handler(triangulation); + + // set active fe index + unsigned int i = 0; + + for (const auto &cell : dof_handler.active_cell_iterators()) + { + if (i == 2) + { + cell->set_active_fe_index(0); + ++i; + } + else + { + cell->set_active_fe_index(1); + ++i; + } + } + + dof_handler.distribute_dofs(fe_collection); + + { + // set constraints + AffineConstraints constraints; + + for (unsigned i = 0; i < Utilities::pow(2, dim); ++i) + constraints.add_constraint(i, {}, i); + + constraints.close(); + + // setup sparsity pattern + DynamicSparsityPattern dsp(dof_handler.n_dofs(), dof_handler.n_dofs()); + DoFTools::make_flux_sparsity_pattern(dof_handler, dsp, constraints); + dsp.compress(); + + // print sparsity pattern + deallog << "nonzero matrix elements: " << dsp.n_nonzero_elements() + << std::endl; + dsp.print(deallog.get_file_stream()); + } +} + + + +int +main() +{ + initlog(); + + deallog.push("2d"); + check<2>(); + deallog.pop(); + deallog.push("3d"); + check<3>(); + deallog.pop(); +} diff --git a/tests/hp/flux_sparsity_03.output b/tests/hp/flux_sparsity_03.output new file mode 100644 index 0000000000..adbee02e6b --- /dev/null +++ b/tests/hp/flux_sparsity_03.output @@ -0,0 +1,15 @@ + +DEAL:2d::nonzero matrix elements: 16 +[0,0,1,2,3] +[1,0,1,2,3] +[2,0,1,2,3] +[3,0,1,2,3] +DEAL:3d::nonzero matrix elements: 64 +[0,0,1,2,3,4,5,6,7] +[1,0,1,2,3,4,5,6,7] +[2,0,1,2,3,4,5,6,7] +[3,0,1,2,3,4,5,6,7] +[4,0,1,2,3,4,5,6,7] +[5,0,1,2,3,4,5,6,7] +[6,0,1,2,3,4,5,6,7] +[7,0,1,2,3,4,5,6,7] -- 2.39.5