From: hartmann Date: Wed, 22 Oct 2003 15:52:34 +0000 (+0000) Subject: Templatize the create_mass_matrix functions also on the template argument of the... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=84577e41137034d05111727adae53f854b3115d3;p=dealii-svn.git Templatize the create_mass_matrix functions also on the template argument of the SparseMatrix class. git-svn-id: https://svn.dealii.org/trunk@8125 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/numerics/matrices.h b/deal.II/deal.II/include/numerics/matrices.h index 2dd2d7cd06..582a195655 100644 --- a/deal.II/deal.II/include/numerics/matrices.h +++ b/deal.II/deal.II/include/numerics/matrices.h @@ -203,11 +203,11 @@ class MatrixCreator * See the general doc of this class * for more information. */ - template + template static void create_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const a = 0); /** @@ -215,10 +215,10 @@ class MatrixCreator * function, see above, with * @p{mapping=MappingQ1()}. */ - template + template static void create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const a = 0); /** @@ -234,11 +234,11 @@ class MatrixCreator * See the general doc of this * class for more information. */ - template + template static void create_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const a = 0); @@ -248,10 +248,10 @@ class MatrixCreator * function, see above, with * @p{mapping=MappingQ1()}. */ - template + template static void create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const a = 0); @@ -445,12 +445,12 @@ class MatrixCreator * used to synchronise access to * the matrix. */ - template + template static void create_mass_matrix_1 (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const a, const IteratorRange range, Threads::ThreadMutex &mutex); @@ -465,12 +465,12 @@ class MatrixCreator * used to synchronise access to * the matrix. */ - template + template static void create_mass_matrix_2 (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const a, diff --git a/deal.II/deal.II/source/numerics/matrices.cc b/deal.II/deal.II/source/numerics/matrices.cc index faa104f129..3e146b28ff 100644 --- a/deal.II/deal.II/source/numerics/matrices.cc +++ b/deal.II/deal.II/source/numerics/matrices.cc @@ -58,11 +58,11 @@ MatrixCreator::IteratorRange::IteratorRange (const iterator_pair &ip) -template +template void MatrixCreator::create_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const coefficient) { Assert (matrix.m() == dof.n_dofs(), @@ -88,11 +88,11 @@ void MatrixCreator::create_mass_matrix (const Mapping &mapping, typedef void (*create_mass_matrix_1_t) (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const coefficient, const IteratorRange range, Threads::ThreadMutex &mutex); - create_mass_matrix_1_t p = &MatrixCreator::template create_mass_matrix_1; + create_mass_matrix_1_t p = &MatrixCreator::template create_mass_matrix_1; for (unsigned int thread=0; thread &mapping, -template +template void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const coefficient, const IteratorRange range, Threads::ThreadMutex &mutex) @@ -125,7 +125,7 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, coefficient->n_components==1 || coefficient->n_components==n_components, ExcComponentMismatch()); - FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); std::vector coefficient_values (n_q_points); std::vector > coefficient_vector_values (n_q_points, Vector (n_components)); @@ -222,10 +222,10 @@ void MatrixCreator::create_mass_matrix_1 (const Mapping &mapping, } -template +template void MatrixCreator::create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function * const coefficient) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); @@ -234,11 +234,11 @@ void MatrixCreator::create_mass_matrix (const DoFHandler &dof, } -template +template void MatrixCreator::create_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const coefficient) @@ -266,13 +266,13 @@ void MatrixCreator::create_mass_matrix (const Mapping &mapping, typedef void (*create_mass_matrix_2_t) (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const coefficient, const IteratorRange range, Threads::ThreadMutex &mutex); - create_mass_matrix_2_t p = &MatrixCreator::template create_mass_matrix_2; + create_mass_matrix_2_t p = &MatrixCreator::template create_mass_matrix_2; for (unsigned int thread=0; thread &mapping, -template +template void MatrixCreator::create_mass_matrix_2 (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const coefficient, @@ -311,7 +311,7 @@ MatrixCreator::create_mass_matrix_2 (const Mapping &mapping, coefficient->n_components==1 || coefficient->n_components==n_components, ExcComponentMismatch()); - FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); Vector local_rhs (dofs_per_cell); std::vector rhs_values (fe_values.n_quadrature_points); std::vector coefficient_values (n_q_points); @@ -418,10 +418,10 @@ MatrixCreator::create_mass_matrix_2 (const Mapping &mapping, } -template +template void MatrixCreator::create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix &matrix, const Function &rhs, Vector &rhs_vector, const Function * const coefficient) @@ -1297,6 +1297,39 @@ void MatrixCreator::create_mass_matrix Vector &rhs_vector, const Function * const coefficient); + +template +void MatrixCreator::create_mass_matrix +(const Mapping &mapping, + const DoFHandler &dof, + const Quadrature &q, + SparseMatrix &matrix, + const Function * const coefficient); +template +void MatrixCreator::create_mass_matrix +(const DoFHandler &dof, + const Quadrature &q, + SparseMatrix &matrix, + const Function * const coefficient); +template +void MatrixCreator::create_mass_matrix +(const Mapping &mapping, + const DoFHandler &dof, + const Quadrature &q, + SparseMatrix &matrix, + const Function &rhs, + Vector &rhs_vector, + const Function * const coefficient); +template +void MatrixCreator::create_mass_matrix +(const DoFHandler &dof, + const Quadrature &q, + SparseMatrix &matrix, + const Function &rhs, + Vector &rhs_vector, + const Function * const coefficient); + + #if deal_II_dimension != 1 template void diff --git a/deal.II/doc/news/c-4-0.html b/deal.II/doc/news/c-4-0.html index 982686e07b..dba522d046 100644 --- a/deal.II/doc/news/c-4-0.html +++ b/deal.II/doc/news/c-4-0.html @@ -256,6 +256,18 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK

deal.II

    +
  1. + Improved: The MatrixCreator::create_mass_matrix + functions are now templatized also on the template argument of + the SparseMatrix class. This allows + invoking this function for SparseMatrix<double> and SparseMatrix<float> objects. +
    + (RH 2003/10/22) +

    +
  2. New: There is now also a function MGDoFCellAccessor::neighbor_child_on_subface