From db6008691d06a82e36fc35f3a49b6df5d36321fd Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 8 Jan 2009 16:01:33 +0000 Subject: [PATCH] Fix a snafu from last night: the function I moved into an internal namespace is still used from another file. Unfortunately, I forgot to also remove the declaration of the function, so the code continued to compile, but wouldn't link. Fix this by re-surrecting the function. By tweaking the interface a bit it actually made for a nice function that could as well be public. git-svn-id: https://svn.dealii.org/trunk@18136 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 28 +++--- deal.II/deal.II/source/dofs/dof_tools.cc | 94 +++++++++---------- .../deal.II/source/multigrid/mg_dof_tools.cc | 17 ++-- 3 files changed, 63 insertions(+), 76 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index b76e1b40c2..1a928aad0d 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -1768,6 +1768,20 @@ class DoFTools convert_couplings_to_blocks (const DoFHandler& dof_handler, const Table<2, Coupling>& table_by_component, std::vector >& tables_by_block); + + /** + * Given a finite element and a table how + * the vector components of it couple + * with each other, compute and return a + * table that describes how the + * individual shape functions couple with + * each other. + */ + template + static + Table<2,Coupling> + dof_couplings_from_component_couplings (const FiniteElement &fe, + const Table<2,Coupling> &component_couplings); /** * Exception @@ -1813,20 +1827,6 @@ class DoFTools * Exception */ DeclException0 (ExcInvalidBoundaryIndicator); - - private: - - /** - * Compute coupling of dofs from - * coupling of components. - */ - template - static void compute_dof_couplings ( - Table<2,Coupling>& dof_couplings, - const Table<2,Coupling>& component_couplings, - const FiniteElement& fe); - - friend class MGTools; }; diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 58e319be5e..1c38e13248 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -681,65 +681,57 @@ DoFTools::make_flux_sparsity_pattern ( #endif -namespace internal +template +Table<2,DoFTools::Coupling> +DoFTools::dof_couplings_from_component_couplings +(const FiniteElement &fe, + const Table<2,Coupling> &component_couplings) { - namespace - { - template - void - compute_dof_couplings(Table<2,DoFTools::Coupling>& dof_couplings, - const Table<2,DoFTools::Coupling>& component_couplings, - const FiniteElement& fe) - { - Assert(component_couplings.n_rows() == fe.n_components(), - ExcDimensionMismatch(component_couplings.n_rows(), - fe.n_components())); - Assert(component_couplings.n_cols() == fe.n_components(), - ExcDimensionMismatch(component_couplings.n_cols(), - fe.n_components())); + Assert(component_couplings.n_rows() == fe.n_components(), + ExcDimensionMismatch(component_couplings.n_rows(), + fe.n_components())); + Assert(component_couplings.n_cols() == fe.n_components(), + ExcDimensionMismatch(component_couplings.n_cols(), + fe.n_components())); - const unsigned int n_dofs = fe.dofs_per_cell; + const unsigned int n_dofs = fe.dofs_per_cell; - Assert(dof_couplings.n_rows() == n_dofs, - ExcDimensionMismatch(dof_couplings.n_rows(), n_dofs)); - Assert(dof_couplings.n_cols() == n_dofs, - ExcDimensionMismatch(dof_couplings.n_cols(), n_dofs)); + Table<2,DoFTools::Coupling> dof_couplings (n_dofs, n_dofs); - for (unsigned int i=0; i int_dof_mask(total_dofs, total_dofs); - Table<2,Coupling> flux_dof_mask(total_dofs, total_dofs); - - internal::compute_dof_couplings(int_dof_mask, int_mask, fe); - internal::compute_dof_couplings(flux_dof_mask, flux_mask, fe); + const Table<2,Coupling> + int_dof_mask = dof_couplings_from_component_couplings(fe, int_mask), + flux_dof_mask = dof_couplings_from_component_couplings(fe, flux_mask); for (unsigned int i=0; i::faces_per_cell;++f) diff --git a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc index 8d04955e8c..742e9d4d6d 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -685,11 +685,9 @@ MGTools::make_flux_sparsity_pattern ( typename MGDoFHandler::cell_iterator cell = dof.begin(level), endc = dof.end(level); - Table<2,DoFTools::Coupling> int_dof_mask(total_dofs, total_dofs); - Table<2,DoFTools::Coupling> flux_dof_mask(total_dofs, total_dofs); - - DoFTools::compute_dof_couplings(int_dof_mask, int_mask, fe); - DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe); + const Table<2,DoFTools::Coupling> + int_dof_mask = DoFTools::dof_couplings_from_component_couplings(fe, int_mask), + flux_dof_mask = DoFTools::dof_couplings_from_component_couplings(fe, flux_mask); for (unsigned int i=0; i::faces_per_cell;++f) @@ -864,11 +862,10 @@ MGTools::make_flux_sparsity_pattern_edge ( Table<2,bool> support_on_face(dofs_per_cell, GeometryInfo::faces_per_cell); typename MGDoFHandler::cell_iterator cell = dof.begin(level), - endc = dof.end(level); - - Table<2,DoFTools::Coupling> flux_dof_mask(dofs_per_cell, dofs_per_cell); + endc = dof.end(level); - DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe); + const Table<2,DoFTools::Coupling> flux_dof_mask + = DoFTools::dof_couplings_from_component_couplings(fe, flux_mask); for (unsigned int i=0; i::faces_per_cell;++f) -- 2.39.5