From 194638b50bc3288e5d2d324f6f583106d8187945 Mon Sep 17 00:00:00 2001 From: bangerth Date: Wed, 8 Jan 2014 22:54:05 +0000 Subject: [PATCH] Fix a testcase by Minh Do-Quang: FESystem got a bit confused with FE_Nothing elements. git-svn-id: https://svn.dealii.org/trunk@32178 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 9 ++- deal.II/source/fe/fe_system.cc | 14 +++- tests/hp/solution_transfer_04.cc | 104 +++++++++++++++++++++++++++ tests/hp/solution_transfer_04.output | 2 + 4 files changed, 125 insertions(+), 4 deletions(-) create mode 100644 tests/hp/solution_transfer_04.cc create mode 100644 tests/hp/solution_transfer_04.output diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 8b705acd64..444e2a9d14 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2013 by the deal.II authors +// Copyright (C) 2013, 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -85,6 +85,13 @@ inconvenience this causes.
    +
  1. Fixed: FESystem::get_interpolation_matrix, a function that is among + other places used by SolutionTransfer, had a bug that prevented it from + running correctly in some situations where one uses FE_Nothing. +
    + (Minh Do-Quang, Wolfgang Bangerth, 2014/01/08) +
  2. +
  3. Improved: When you call WorkStream::run with an empty function object for the copier, operations on individual cells are essentially all independent. In other words, you have a massively parallel collection of jobs. In this diff --git a/deal.II/source/fe/fe_system.cc b/deal.II/source/fe/fe_system.cc index bf70ec76a6..8683136e66 100644 --- a/deal.II/source/fe/fe_system.cc +++ b/deal.II/source/fe/fe_system.cc @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 1999 - 2013 by the deal.II authors +// Copyright (C) 1999 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -528,10 +528,18 @@ FESystem::get_interpolation_matrix ( const FiniteElement &x_source_fe, FullMatrix &interpolation_matrix) const { - Assert (interpolation_matrix.m() == this->dofs_per_cell, + // check that the size of the matrices is correct. for historical + // reasons, if you call matrix.reinit(8,0), it sets the sizes + // to m==n==0 internally. this may happen when we use a FE_Nothing, + // so write the test in a more lenient way + Assert ((interpolation_matrix.m() == this->dofs_per_cell) + || + (x_source_fe.dofs_per_cell == 0), ExcDimensionMismatch (interpolation_matrix.m(), this->dofs_per_cell)); - Assert (interpolation_matrix.n() == x_source_fe.dofs_per_cell, + Assert ((interpolation_matrix.n() == x_source_fe.dofs_per_cell) + || + (this->dofs_per_cell == 0), ExcDimensionMismatch (interpolation_matrix.m(), x_source_fe.dofs_per_cell)); diff --git a/tests/hp/solution_transfer_04.cc b/tests/hp/solution_transfer_04.cc new file mode 100644 index 0000000000..985f6c6829 --- /dev/null +++ b/tests/hp/solution_transfer_04.cc @@ -0,0 +1,104 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2011 - 2014 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. +// +// --------------------------------------------------------------------- + + + +// testcase by Minh Do-Quang: a case where SolutionTransfer got into trouble +// in a couple of places when using FE_Nothing and FESystem. + +#include "../tests.h" +#include +#include +#include + +#include + +#include + +#include +#include +#include +#include + +#include +#include +#include +#include + +#include + +#include + +#include + +using namespace dealii; + +int main() +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + Triangulation<2> triangulation; + GridGenerator::hyper_cube (triangulation); + triangulation.refine_global (1); + + + hp::FECollection<2> fe_collection; + + fe_collection.push_back(FESystem<2>(FE_Q<2>(1), 1, FE_Q<2>(1), 1)); + fe_collection.push_back(FESystem<2>(FE_Nothing<2>(), 1, FE_Nothing<2>(), 1)); + + hp::DoFHandler<2> dof_handler(triangulation); + + // Assign FE to cells + hp::DoFHandler<2>::active_cell_iterator cell; + hp::DoFHandler<2>::active_cell_iterator endc = dof_handler.end(); + + + cell = dof_handler.begin_active(); + cell->set_active_fe_index(1); + cell++; + cell->set_active_fe_index(1); + cell++; + cell->set_active_fe_index(0); + cell++; + cell->set_active_fe_index(0); + + + dof_handler.distribute_dofs (fe_collection); + + // Init solution + Vector solution(dof_handler.n_dofs()); + solution = 1.0; + + + SolutionTransfer<2, Vector, hp::DoFHandler<2> > solultion_trans(dof_handler); + solultion_trans.prepare_for_coarsening_and_refinement(solution); + + triangulation.execute_coarsening_and_refinement (); + dof_handler.distribute_dofs (fe_collection); + + Vector new_solution(dof_handler.n_dofs()); + solultion_trans.interpolate(solution, new_solution); + + // we are good if we made it to here + deallog << "OK" << std::endl; +} + + + diff --git a/tests/hp/solution_transfer_04.output b/tests/hp/solution_transfer_04.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/hp/solution_transfer_04.output @@ -0,0 +1,2 @@ + +DEAL::OK -- 2.39.5