From db3fcb136183813147197287e70cdcd81dac7f22 Mon Sep 17 00:00:00 2001 From: bangerth Date: Sun, 26 Jan 2014 18:47:48 +0000 Subject: [PATCH] Fix a forgotten case with FENothing (thanks to a testcase by Krzysztof Bzowski). git-svn-id: https://svn.dealii.org/trunk@32316 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/source/dofs/dof_accessor_set.cc | 4 +- tests/hp/solution_transfer_09.cc | 105 +++++++++++++++++++++ tests/hp/solution_transfer_09.output | 4 + tests/hp/solution_transfer_10.cc | 119 ++++++++++++++++++++++++ tests/hp/solution_transfer_10.output | 2 + 5 files changed, 233 insertions(+), 1 deletion(-) create mode 100644 tests/hp/solution_transfer_09.cc create mode 100644 tests/hp/solution_transfer_09.output create mode 100644 tests/hp/solution_transfer_10.cc create mode 100644 tests/hp/solution_transfer_10.output diff --git a/deal.II/source/dofs/dof_accessor_set.cc b/deal.II/source/dofs/dof_accessor_set.cc index b4fef1b51f..5a63df6a36 100644 --- a/deal.II/source/dofs/dof_accessor_set.cc +++ b/deal.II/source/dofs/dof_accessor_set.cc @@ -112,7 +112,9 @@ set_dof_values_by_interpolation (const Vector &local_values, for (unsigned int child=0; childn_children(); ++child) { - fe.get_prolongation_matrix(child, this->refinement_case()).vmult (tmp, local_values); + if (tmp.size() > 0) + fe.get_prolongation_matrix(child, this->refinement_case()) + .vmult (tmp, local_values); this->child(child)->set_dof_values_by_interpolation (tmp, values, fe_index); } } diff --git a/tests/hp/solution_transfer_09.cc b/tests/hp/solution_transfer_09.cc new file mode 100644 index 0000000000..1527b87d5a --- /dev/null +++ b/tests/hp/solution_transfer_09.cc @@ -0,0 +1,105 @@ +// --------------------------------------------------------------------- +// $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. +// +// --------------------------------------------------------------------- + + + +// Like _07 but with all same fe indices. This triggered yet another place +// where we had the same kind of error. + +#include "../tests.h" +#include +#include +#include + +#include + +#include + +#include +#include +#include +#include + +#include +#include +#include +#include + +#include + +#include + +#include + +using namespace dealii; + + +template +void test () +{ + Triangulation triangulation; + GridGenerator::hyper_cube (triangulation); + + hp::FECollection fe_collection; + fe_collection.push_back(FE_Q(1)); + fe_collection.push_back(FE_Q(2)); + + hp::DoFHandler dof_handler(triangulation); + dof_handler.begin_active()->set_active_fe_index(0); + dof_handler.distribute_dofs (fe_collection); + + // Init solution + Vector solution(dof_handler.n_dofs()); + solution = 1.0; + + + // set refine flag for the only cell we have, then do the refinement + SolutionTransfer, hp::DoFHandler > + solution_trans(dof_handler); + dof_handler.begin_active()->set_refine_flag (); + solution_trans.prepare_for_pure_refinement(); + triangulation.execute_coarsening_and_refinement (); + + // now set the active_fe_index flags on the new set of fine level cells + for (unsigned int c=0; cn_children(); ++c) + dof_handler.begin(0)->child(c)->set_active_fe_index(0); + + // distribute dofs and transfer solution there + dof_handler.distribute_dofs (fe_collection); + + Vector new_solution(dof_handler.n_dofs()); + solution_trans.refine_interpolate(solution, new_solution); + + // we should now have only 1s in the new_solution vector + for (unsigned int i=0; i (); + test<2> (); + test<3> (); +} diff --git a/tests/hp/solution_transfer_09.output b/tests/hp/solution_transfer_09.output new file mode 100644 index 0000000000..fb71de2867 --- /dev/null +++ b/tests/hp/solution_transfer_09.output @@ -0,0 +1,4 @@ + +DEAL::OK +DEAL::OK +DEAL::OK diff --git a/tests/hp/solution_transfer_10.cc b/tests/hp/solution_transfer_10.cc new file mode 100644 index 0000000000..75ac94f5fb --- /dev/null +++ b/tests/hp/solution_transfer_10.cc @@ -0,0 +1,119 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 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. +// +// --------------------------------------------------------------------- + +// A test by Krzysztof Bzowski that verifies something in SolutionTransfer +// that didn't work for a few days + +#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(Triangulation<2>::none); + GridGenerator::hyper_cube (triangulation); + triangulation.refine_global(1); + + hp::FECollection<2> fe_collection; + fe_collection.push_back(FE_Q<2>(1)); + fe_collection.push_back(FE_Nothing<2>()); + + hp::DoFHandler<2> dof_handler(triangulation); + + // Assign FEQ to all cells + hp::DoFHandler<2>::active_cell_iterator cell = dof_handler.begin_active(); + hp::DoFHandler<2>::active_cell_iterator endc = dof_handler.end(); + + + + /* + * ----------- + * | 0 | 0 | + * ----------- + * | 1 | 1 | 0 - FEQ, 1 - FE_Nothing + * ----------- + */ + + 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; + + + /* Set refine flags: + * ----------- + * | | R | + * ----------- + * | R | | + * ----------- + */ + + + cell = dof_handler.begin_active(); + cell->set_refine_flag(); + cell++; + cell++; + cell++; + cell->set_refine_flag(); + + triangulation.prepare_coarsening_and_refinement(); + + // Interpolate solution + 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); + + deallog << "OK" << std::endl; +} + diff --git a/tests/hp/solution_transfer_10.output b/tests/hp/solution_transfer_10.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/hp/solution_transfer_10.output @@ -0,0 +1,2 @@ + +DEAL::OK -- 2.39.5