From 4ac80fa51701fba16a275d894e944f6abc2f0ba4 Mon Sep 17 00:00:00 2001 From: goll Date: Thu, 2 Feb 2012 09:25:53 +0000 Subject: [PATCH] DoFTools::map_dofs_to_support_points now works within the hp framekwork. git-svn-id: https://svn.dealii.org/trunk@24975 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 4 + deal.II/examples/step-27/step-27.cc | 9 ++ deal.II/include/deal.II/dofs/dof_tools.h | 12 +++ deal.II/source/dofs/dof_tools.cc | 126 ++++++++++++++++------- deal.II/source/dofs/dof_tools.inst.in | 15 +++ 5 files changed, 130 insertions(+), 36 deletions(-) diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index d5d7e084b8..b63d73834a 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -81,6 +81,10 @@ enabled due to a missing include file in file

Specific improvements

    +
  1. Improved: DoFTools::map_dofs_to_support_points() now also works within the hp framework. +
    +(Christian Goll, 2012/02/02) +
  2. Improved: There is now a constructor for FESystem that allows to create collections of finite elements of arbitrary length.
    diff --git a/deal.II/examples/step-27/step-27.cc b/deal.II/examples/step-27/step-27.cc index a3c091558d..8c8d5ec7b3 100644 --- a/deal.II/examples/step-27/step-27.cc +++ b/deal.II/examples/step-27/step-27.cc @@ -454,6 +454,15 @@ namespace Step27 solution, estimated_error_per_cell); + std::vector > support_points(dof_handler.n_dofs()); + + DoFTools::map_dofs_to_support_points(hp::StaticMappingQ1::mapping_collection, + dof_handler, support_points); + + + dealii::Point p(0.5,0.75); + std::cout << VectorTools::point_value(dof_handler, solution, p) << std::endl; + Vector smoothness_indicators (triangulation.n_active_cells()); estimate_smoothness (smoothness_indicators); diff --git a/deal.II/include/deal.II/dofs/dof_tools.h b/deal.II/include/deal.II/dofs/dof_tools.h index e5202f3d36..d8d3d4e6a1 100644 --- a/deal.II/include/deal.II/dofs/dof_tools.h +++ b/deal.II/include/deal.II/dofs/dof_tools.h @@ -23,6 +23,7 @@ #include #include #include +#include #include #include @@ -39,6 +40,7 @@ template class DoFHandler; namespace hp { template class DoFHandler; + template class MappingCollection; } template class MGDoFHandler; class ConstraintMatrix; @@ -2192,6 +2194,16 @@ namespace DoFTools const DoFHandler &dof_handler, std::vector > &support_points); + /** + * Same as above for the hp case. + */ + + template + void + map_dofs_to_support_points (const dealii::hp::MappingCollection &mapping, + const hp::DoFHandler &dof_handler, + std::vector > &support_points); + /** * This is the opposite function * to the one above. It generates diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index a071a393d9..5adc456716 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -36,6 +36,8 @@ #include #include #include +#include +#include #include #include @@ -5688,51 +5690,103 @@ namespace DoFTools dof_handler.n_boundary_dofs (boundary_indicators)); } + namespace internal + { + namespace + { + template + void map_dofs_to_support_points_internal( + const hp::MappingCollection & mapping, + const DH &dof_handler, + std::vector > &support_points) + { + AssertDimension(support_points.size(), dof_handler.n_dofs()); + + const unsigned int dim = DH::dimension; + const unsigned int spacedim = DH::space_dimension; + + hp::FECollection fe_collection(dof_handler.get_fe()); + hp::QCollection q_coll_dummy; + + for (unsigned int fe_index = 0; fe_index < fe_collection.size(); ++fe_index) + { + // check whether every fe in the collection + // has support points + Assert(fe_collection[fe_index].has_support_points(), + typename FiniteElement::ExcFEHasNoSupportPoints()); + q_coll_dummy.push_back( + Quadrature ( + fe_collection[fe_index].get_unit_support_points())); + } + // now loop over all cells and + // enquire the support points on + // each of these. we use dummy + // quadrature formulas where the + // quadrature points are located at + // the unit support points to + // enquire the location of the + // support points in real space + // + // the weights of the quadrature + // rule have been set to invalid values + // by the used constructor. + hp::FEValues hp_fe_values(mapping, fe_collection, + q_coll_dummy, update_quadrature_points); + typename DH::active_cell_iterator cell = + dof_handler.begin_active(), endc = dof_handler.end(); + + std::vector local_dof_indices; + for (; cell != endc; ++cell) + { + hp_fe_values.reinit(cell); + const FEValues &fe_values = hp_fe_values.get_present_fe_values(); + + local_dof_indices.resize(cell->get_fe().dofs_per_cell); + cell->get_dof_indices(local_dof_indices); + + const std::vector > & points = + fe_values.get_quadrature_points(); + for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i) + support_points[local_dof_indices[i]] = points[i]; + }; + } + } + } template void map_dofs_to_support_points (const Mapping &mapping, - const DoFHandler &dof_handler, - std::vector > &support_points) + const DoFHandler &dof_handler, + std::vector > &support_points) { - const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell; + AssertDimension(support_points.size(), dof_handler.n_dofs()); - // check whether fe has support - // points - Assert (dof_handler.get_fe().has_support_points(), - typename FiniteElement::ExcFEHasNoSupportPoints()); - AssertDimension (support_points.size(), dof_handler.n_dofs()); + //Let the internal function do all the work, + //just make sure that it gets a MappingCollection + const hp::MappingCollection mapping_collection(mapping); - // now loop over all cells and - // enquire the support points on - // each of these. use a dummy - // quadrature formula where the - // quadrature points are located at - // the unit support points to - // enquire the location of the - // support points in real space - // - // the weights of the quadrature - // rule are set to invalid values - // by the used constructor. - Quadrature q_dummy(dof_handler.get_fe().get_unit_support_points()); - FEValues fe_values (mapping, dof_handler.get_fe(), - q_dummy, update_quadrature_points); - typename DoFHandler::active_cell_iterator - cell = dof_handler.begin_active(), - endc = dof_handler.end(); + internal::map_dofs_to_support_points_internal( + mapping_collection, + dof_handler, support_points); + } - std::vector local_dof_indices (dofs_per_cell); - for (; cell!=endc; ++cell) - { - fe_values.reinit (cell); - cell->get_dof_indices (local_dof_indices); - const std::vector > & points - = fe_values.get_quadrature_points (); - for (unsigned int i=0; i + void + map_dofs_to_support_points( + const hp::MappingCollection &mapping, + const hp::DoFHandler &dof_handler, + std::vector > &support_points) + { + AssertDimension(support_points.size(), dof_handler.n_dofs()); + + //Let the internal function do all the work, + //just make sure that it gets a MappingCollection + internal::map_dofs_to_support_points_internal( + mapping, + dof_handler, + support_points); } diff --git a/deal.II/source/dofs/dof_tools.inst.in b/deal.II/source/dofs/dof_tools.inst.in index 996ba002a9..af7cdb324f 100644 --- a/deal.II/source/dofs/dof_tools.inst.in +++ b/deal.II/source/dofs/dof_tools.inst.in @@ -618,6 +618,14 @@ DoFTools::map_dofs_to_support_points (const Mapping&, const DoFHandler&, std::vector >&); + + +template +void +DoFTools::map_dofs_to_support_points +(const hp::MappingCollection&, + const hp::DoFHandler&, + std::vector >&); #if deal_II_dimension < 3 @@ -627,6 +635,13 @@ DoFTools::map_dofs_to_support_points (const Mapping&, const DoFHandler&, std::vector >&); + +template +void +DoFTools::map_dofs_to_support_points +(const hp::MappingCollection&, + const hp::DoFHandler&, + std::vector >&); template void -- 2.39.5