From 1ddaa3c71724dd5f4620d9174fde7ef7487a3ed7 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 24 Apr 2001 12:33:45 +0000 Subject: [PATCH] Implement fully DoFTools::map_dofs_to_support_points. git-svn-id: https://svn.dealii.org/trunk@4475 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/dofs/dof_tools.cc | 23 ++++++++++++++--------- 1 file changed, 14 insertions(+), 9 deletions(-) diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index cdb2b67540..d92f5bc05a 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -1880,10 +1880,17 @@ DoFTools::map_dofs_to_support_points (const Mapping &mapping, // now loop over all cells and // enquire the support points on - // each of these - QMidpoint q_dummy; + // 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 + Quadrature q_dummy(dof_handler.get_fe().get_unit_support_points(), + std::vector (dofs_per_cell, + 1./dofs_per_cell)); FEValues fe_values (mapping, dof_handler.get_fe(), - q_dummy, update_support_points); + q_dummy, update_q_points); typename DoFHandler::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); @@ -1893,12 +1900,10 @@ DoFTools::map_dofs_to_support_points (const Mapping &mapping, { fe_values.reinit (cell); cell->get_dof_indices (local_dof_indices); -//TODO:[WB] Uncomment code once FEValues::get_support_points is available. Add a ChangeLog entry then. - Assert (false, ExcNotImplemented()); -// const std::vector > & support_points -// = fe_values.get_support_points (); -// for (unsigned int i=0; i > & points + = fe_values.get_quadrature_points (); + for (unsigned int i=0; i