From cbedb2ed2a8cce0a002c2ec22f42cde581db9e0f Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 20 Jun 2005 20:18:36 +0000 Subject: [PATCH] Several simplifications in FETools::get_projection_matrix. git-svn-id: https://svn.dealii.org/trunk@10893 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/source/fe/fe_tools.cc | 76 ++++++++++++--------------- 1 file changed, 34 insertions(+), 42 deletions(-) diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 66e21d7976..dc14f1ff18 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -457,15 +457,14 @@ void FETools::get_projection_matrix (const FiniteElement &fe1, { b = 0.; for (unsigned int i=0;i& fe, Assert(matrices[i].m() == n, ExcDimensionMismatch(matrices[i].m(),n)); } - /* - * Set up two meshes, one with a - * single reference cell and the - * other refined, together with - * DoFHandler and finite elements. - */ - Triangulation tr_coarse; - Triangulation tr_fine; - GridGenerator::hyper_cube (tr_coarse, 0, 1); - GridGenerator::hyper_cube (tr_fine, 0, 1); - tr_fine.refine_global(1); - DoFHandler dof_coarse(tr_coarse); - dof_coarse.distribute_dofs(fe); - DoFHandler dof_fine(tr_fine); - dof_fine.distribute_dofs(fe); + // Set up a meshes, one with a single + // reference cell and refine it once + Triangulation tria; + GridGenerator::hyper_cube (tria, 0, 1); + tria.refine_global(1); MappingCartesian mapping; QGauss q_fine(degree+1); const unsigned int nq = q_fine.n_quadrature_points; FEValues fine (mapping, fe, q_fine, - update_q_points | update_JxW_values | update_values); - - typename DoFHandler::active_cell_iterator coarse_cell - = dof_coarse.begin_active(); - typename DoFHandler::active_cell_iterator fine_cell; + update_q_points | update_JxW_values | update_values); /** * We search for the polynomial on @@ -532,7 +517,7 @@ FETools::compute_embedding_matrices(const FiniteElement& fe, * This matrix is the same for all * children. */ - fine.reinit(dof_fine.begin_active()); + fine.reinit(tria.begin_active()); FullMatrix A(nq*nd, n); for (unsigned int d=0;d& fe, Vector v_fine(n); unsigned int cell_number = 0; - for (fine_cell = dof_fine.begin_active(); - fine_cell != dof_fine.end(); - ++fine_cell, ++cell_number) + for (typename Triangulation::active_cell_iterator fine_cell + = tria.begin_active(); + fine_cell != tria.end(); ++fine_cell, ++cell_number) { - FullMatrix& matrix = matrices[cell_number]; fine.reinit(fine_cell); - Quadrature q_coarse (fine.get_quadrature_points(), - fine.get_JxW_values()); + + // evaluate on the coarse cell (which + // is the first -- inactive -- cell on + // the lowest level of the + // triangulation we have created) + const Quadrature q_coarse (fine.get_quadrature_points(), + fine.get_JxW_values()); FEValues coarse (mapping, fe, q_coarse, update_values); - coarse.reinit(coarse_cell); + coarse.reinit(tria.begin(0)); + FullMatrix &this_matrix = matrices[cell_number]; + // Compute this once for each // coarse grid basis function for (unsigned int i=0;i& fe, // function, the columns // are fine grid. for (unsigned int j=0;j void FETools::compute_projection_matrices(const FiniteElement& fe, -- 2.39.5