From 863c71bd8944f52bf1159d6f45c674ae42059673 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 18 Apr 2006 21:52:12 +0000 Subject: [PATCH] Make VectorTools::project fit for other types of vectors as well. git-svn-id: https://svn.dealii.org/trunk@12859 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/numerics/vectors.h | 17 ++++---- .../include/numerics/vectors.templates.h | 39 +++++++++++++------ deal.II/deal.II/source/numerics/vectors.cc | 23 +---------- .../source/numerics/vectors.instance.h | 23 +++++++++++ deal.II/doc/news/changes.html | 8 ++++ tests/deal.II/Makefile | 3 +- 6 files changed, 70 insertions(+), 43 deletions(-) diff --git a/deal.II/deal.II/include/numerics/vectors.h b/deal.II/deal.II/include/numerics/vectors.h index 7d17b68c68..66a4e4190b 100644 --- a/deal.II/deal.II/include/numerics/vectors.h +++ b/deal.II/deal.II/include/numerics/vectors.h @@ -25,7 +25,7 @@ template class Point; template class Function; template class FunctionMap; template class Quadrature; -template class QGauss2; +template class QGauss; template class Vector; template class FullMatrix; @@ -491,15 +491,15 @@ class VectorTools * See the general documentation of this * class for further information. */ - template + template static void project (const Mapping &mapping, const DoFHandler &dof, const ConstraintMatrix &constraints, const Quadrature &quadrature, const Function &function, - Vector &vec, + VECTOR &vec, const bool enforce_zero_boundary = false, - const Quadrature &q_boundary = QGauss2(), + const Quadrature &q_boundary = QGauss(2), const bool project_to_boundary_first = false); /** @@ -512,12 +512,13 @@ class VectorTools * quadrature formula is an invalid * object since it makes no sense in 1d. */ + template static void project (const Mapping<1> &mapping, const DoFHandler<1> &dof, const ConstraintMatrix &constraints, const Quadrature<1> &quadrature, const Function<1> &function, - Vector &vec, + VECTOR &vec, const bool enforce_zero_boundary = false, const Quadrature<0> &q_boundary = *invalid_face_quadrature, const bool project_to_boundary_first = false); @@ -527,14 +528,14 @@ class VectorTools * function, see above, with * mapping=MappingQ1@(). */ - template + template static void project (const DoFHandler &dof, const ConstraintMatrix &constraints, const Quadrature &quadrature, const Function &function, - Vector &vec, + VECTOR &vec, const bool enforce_zero_boundary = false, - const Quadrature &q_boundary = QGauss2(), + const Quadrature &q_boundary = QGauss(2), const bool project_to_boundary_first = false); /** diff --git a/deal.II/deal.II/include/numerics/vectors.templates.h b/deal.II/deal.II/include/numerics/vectors.templates.h index c1e57726a5..116bfe9bc0 100644 --- a/deal.II/deal.II/include/numerics/vectors.templates.h +++ b/deal.II/deal.II/include/numerics/vectors.templates.h @@ -273,22 +273,25 @@ VectorTools::interpolate (const DoFHandler &dof_1, #if deal_II_dimension == 1 +template void VectorTools::project (const Mapping<1> &, const DoFHandler<1> &, const ConstraintMatrix &, const Quadrature<1> &, const Function<1> &, - Vector &, + VECTOR &, const bool , const Quadrature<0> &, const bool ) { - // this function should easily be implemented - // using the template below. However some - // changes have to be made since faces don't - // exist in 1D. Maybe integrate the creation of - // zero boundary values into the - // project_boundary_values function? + // this function should easily be + // implemented using the template + // below. However some changes have + // to be made since faces don't + // exist in 1D. Maybe integrate the + // creation of zero boundary values + // into the project_boundary_values + // function? Assert (false, ExcNotImplemented()); } @@ -296,19 +299,22 @@ void VectorTools::project (const Mapping<1> &, #endif -template +template void VectorTools::project (const Mapping &mapping, const DoFHandler &dof, const ConstraintMatrix &constraints, const Quadrature &quadrature, const Function &function, - Vector &vec, + VECTOR &vec_result, const bool enforce_zero_boundary, const Quadrature &q_boundary, const bool project_to_boundary_first) { Assert (dof.get_fe().n_components() == function.n_components, ExcInvalidFE()); + + Assert (vec_result.size() == dof.n_dofs(), + ExcDimensionMismatch (vec_result.size(), dof.n_dofs())); const FiniteElement &fe = dof.get_fe(); @@ -377,7 +383,7 @@ void VectorTools::project (const Mapping &mapping, // set up mass matrix and right hand side - vec.reinit (dof.n_dofs()); + Vector vec (dof.n_dofs()); SparsityPattern sparsity(dof.n_dofs(), dof.n_dofs(), dof.max_couplings_between_dofs()); @@ -409,15 +415,24 @@ void VectorTools::project (const Mapping &mapping, // distribute solution constraints.distribute (vec); + + // copy vec into vec_result. we + // can't use ve_result itself + // above, since it may be of + // another type than Vector + // and that wouldn't necessarily go + // together with the matrix and + // other functions + std::copy (vec.begin(), vec.end(), vec_result.begin()); } -template +template void VectorTools::project (const DoFHandler &dof, const ConstraintMatrix &constraints, const Quadrature &quadrature, const Function &function, - Vector &vec, + VECTOR &vec, const bool enforce_zero_boundary, const Quadrature &q_boundary, const bool project_to_boundary_first) diff --git a/deal.II/deal.II/source/numerics/vectors.cc b/deal.II/deal.II/source/numerics/vectors.cc index bf045fa4e7..a69fb010fb 100644 --- a/deal.II/deal.II/source/numerics/vectors.cc +++ b/deal.II/deal.II/source/numerics/vectors.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -88,16 +88,6 @@ namespace #include "vectors.instance.h" #undef VEC -template -void VectorTools::project -(const DoFHandler &, - const ConstraintMatrix &, - const Quadrature &, - const Function &, - Vector &, - const bool, - const Quadrature &, - const bool); template void VectorTools::create_right_hand_side @@ -193,17 +183,6 @@ void VectorTools::interpolate_boundary_values std::map &, const std::vector &); -template -void VectorTools::project -(const Mapping &, - const DoFHandler &, - const ConstraintMatrix &, - const Quadrature &, - const Function &, - Vector &, - const bool, - const Quadrature &, - const bool); #endif diff --git a/deal.II/deal.II/source/numerics/vectors.instance.h b/deal.II/deal.II/source/numerics/vectors.instance.h index 4faf6c45a1..7178be51a1 100644 --- a/deal.II/deal.II/source/numerics/vectors.instance.h +++ b/deal.II/deal.II/source/numerics/vectors.instance.h @@ -114,3 +114,26 @@ double VectorTools::compute_mean_value const Quadrature&, const VEC&, const unsigned int); + +template +void VectorTools::project +(const Mapping &, + const DoFHandler &, + const ConstraintMatrix &, + const Quadrature &, + const Function &, + VEC &, + const bool, + const Quadrature &, + const bool); + +template +void VectorTools::project +(const DoFHandler &, + const ConstraintMatrix &, + const Quadrature &, + const Function &, + VEC &, + const bool, + const Quadrature &, + const bool); diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index d6750f954d..021e7db38d 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -583,6 +583,14 @@ inconvenience this causes.

deal.II

    +
  1. + New: The function VectorTools::project functions can now + also be used for vector arguments of type other than + Vector<double>. +
    + (WB 2006/04/17) +

    +
  2. New: The function GridTools::get_finest_common_cells can be used to find those cells that are shared by two triangulations that are diff --git a/tests/deal.II/Makefile b/tests/deal.II/Makefile index 31c3447706..9d837e4ff2 100644 --- a/tests/deal.II/Makefile +++ b/tests/deal.II/Makefile @@ -53,7 +53,8 @@ tests_x = block_matrices \ subcelldata \ interpolate_boundary_values_* \ have_same_coarse_mesh_* \ - get_finest_common_cells_* + get_finest_common_cells_* \ + project_* # from above list of regular expressions, generate the real set of # tests -- 2.39.5