From: Wolfgang Bangerth Date: Wed, 9 Feb 2011 15:30:57 +0000 (+0000) Subject: Compute restriction and prolongation matrices for DGQ elements also if dimSpecific improvements
    +
  1. Fixed: Prolongation and restriction matrices were not computed at all +for elements of type FE_DGQ if dim@. Consequently, +consumers of this information, such as the SolutionTransfer class or +the DoFCellAccess::set_dof_values_by_interpolation function did not +work either and simply returned zero results. This is now fixed. +
    +(M. Sebastian Pauletti, Wolfgang Bangerth, 2011/02/09) +
  2. +
  3. Fixed: When refining a mesh with codimension one, edges were refined using the same manifold description as adjacent cells, but this ignored that a boundary indicator might have been purposefully set for edges that are truly at diff --git a/deal.II/source/fe/fe_dgq.cc b/deal.II/source/fe/fe_dgq.cc index 3577c85c56..11d3216cca 100644 --- a/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/source/fe/fe_dgq.cc @@ -142,12 +142,32 @@ FE_DGQ::FE_DGQ (const unsigned int degree) // restriction and prolongation // matrices to the right sizes this->reinit_restriction_and_prolongation_matrices(); - // Fill prolongation matrices with embedding operators - if (dim ==spacedim) - FETools::compute_embedding_matrices (*this, this->prolongation); - // Fill restriction matrices with L2-projection - if (dim ==spacedim) - FETools::compute_projection_matrices (*this, this->restriction); + + // Fill prolongation matrices with + // embedding operators and + // restriction with L2 projection + // + // we can call the respective + // functions in the case + // dim==spacedim, but they are not + // currently implemented for the + // codim-1 case. there's no harm + // here, though, since these + // matrices are independent of the + // space dimension and so we can + // copy them from the respective + // elements (not cheap, but works) + if (dim == spacedim) + { + FETools::compute_embedding_matrices (*this, this->prolongation); + FETools::compute_projection_matrices (*this, this->restriction); + } + else + { + FE_DGQ tmp (degree); + this->prolongation = tmp.prolongation; + this->restriction = tmp.restriction; + } // finally fill in support points diff --git a/tests/codim_one/fe_dgq_prolongation_01.cc b/tests/codim_one/fe_dgq_prolongation_01.cc new file mode 100644 index 0000000000..567cdf9a6a --- /dev/null +++ b/tests/codim_one/fe_dgq_prolongation_01.cc @@ -0,0 +1,45 @@ +//---------------------------- fe_dgq_prolongation_01.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2010, 2011 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- fe_dgq_prolongation_01.cc --------------------------- + +// this is the root cause for solution_tranfer_01: the prolongation +// matrices for FE_DGQ were not computed at all + +#include "../tests.h" +#include + + +#include + +using namespace dealii; + + +int main () +{ + std::ofstream logfile("fe_dgq_prolongation_01/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + const unsigned int spacedim = 2; + const unsigned int dim = spacedim-1; + + for (unsigned int degree=0; degree<3; ++degree) + { + deallog << "Degree=" << degree << std::endl; + FE_DGQ fe(degree); + for (unsigned int i=0; i +#include +#include + +#include +#include +#include +#include + +#include + +#include +#include + +#include + +using namespace dealii; + + +int main () +{ + std::ofstream logfile("solution_transfer_01/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + const unsigned int spacedim = 2; + const unsigned int dim = spacedim-1; + + Triangulation boundary_mesh; + + // create a line in 2d + GridGenerator::hyper_cube (boundary_mesh); + + // attach a piecewise constant + // element to this one cell + FE_DGQ fe(0); + DoFHandler dh (boundary_mesh); + dh.distribute_dofs(fe); + + Vector solution(dh.n_dofs()); + solution = 1.0; + + deallog << "Old values:" << std::endl; + for (unsigned int i=0; iset_refine_flag (); + + SolutionTransfer, DoFHandler > + soltrans(dh); + + boundary_mesh.prepare_coarsening_and_refinement(); + + soltrans.prepare_for_coarsening_and_refinement(solution); + boundary_mesh.execute_coarsening_and_refinement (); + + dh.distribute_dofs(fe); + + // get the interpolated solution + // back + Vector tmp(dh.n_dofs()); + tmp = 2; + soltrans.interpolate(solution, tmp); + + deallog << "New values:" << std::endl; + for (unsigned int i=0; i