From 23ccb5921ce1ddb565d2f035f14f9198542288fb Mon Sep 17 00:00:00 2001 From: bangerth Date: Mon, 14 Aug 2006 23:25:33 +0000 Subject: [PATCH] Change a few function parameters from pointer to FullMatrix to reference to array of FullMatrix git-svn-id: https://svn.dealii.org/trunk@13705 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_tools.h | 44 ++++++++++--------- deal.II/deal.II/source/fe/fe_abf.cc | 4 +- deal.II/deal.II/source/fe/fe_dgp.cc | 4 +- deal.II/deal.II/source/fe/fe_dgp_monomial.cc | 4 +- deal.II/deal.II/source/fe/fe_dgq.cc | 8 ++-- .../deal.II/source/fe/fe_raviart_thomas.cc | 12 ++--- .../source/fe/fe_raviart_thomas_nodal.cc | 12 ++--- deal.II/deal.II/source/fe/fe_tools.cc | 20 ++++----- deal.II/doc/news/changes.html | 17 +++++++ 9 files changed, 73 insertions(+), 52 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h index a97f56d041..faca7d90de 100644 --- a/deal.II/deal.II/include/fe/fe_tools.h +++ b/deal.II/deal.II/include/fe/fe_tools.h @@ -264,15 +264,16 @@ class FETools * class for which we compute the * embedding matrices. * @param matrices A pointer to - * 2dim FullMatrix + * GeometryInfo::children_per_cell=2dim FullMatrix * objects. This is the format * used in FiniteElement, * where we want to use ths * function mostly. */ template - static void compute_embedding_matrices(const FiniteElement &fe, - FullMatrix* matrices); + static void + compute_embedding_matrices(const FiniteElement &fe, + FullMatrix (&matrices)[GeometryInfo::children_per_cell]); /** * Compute the embedding matrices @@ -281,17 +282,18 @@ class FETools * * @param fe The finite element * for which to compute these - * matrices. - * @param matrices An array of - * 2dim-1 FullMatrix - * objects,holding the embedding - * matrix for each subface. - * @param face_coarse The number - * of the face on the coarse side - * of the face for which this is - * computed. - * @param face_fine The number - * of the face on the refined side + * matrices. @param matrices An + * array of + * GeometryInfo::subfaces_per_face + * = 2dim-1 + * FullMatrix objects,holding the + * embedding matrix for each + * subface. @param face_coarse + * The number of the face on the + * coarse side of the face for + * which this is computed. + * @param face_fine The number of + * the face on the refined side * of the face for which this is * computed. * @@ -301,10 +303,11 @@ class FETools * sufficiently tested yet. */ template - static void compute_face_embedding_matrices(const FiniteElement& fe, - FullMatrix* matrices, - unsigned int face_coarse, - unsigned int face_fine); + static void + compute_face_embedding_matrices(const FiniteElement& fe, + FullMatrix (&matrices)[GeometryInfo::subfaces_per_face], + const unsigned int face_coarse, + const unsigned int face_fine); /** * Compute the @@ -321,8 +324,9 @@ class FETools * want to use this function mostly. */ template - static void compute_projection_matrices(const FiniteElement &fe, - FullMatrix* matrices); + static void + compute_projection_matrices(const FiniteElement &fe, + FullMatrix (&matrices)[GeometryInfo::children_per_cell]); //TODO:[WB] Replace this documentation by something comprehensible diff --git a/deal.II/deal.II/source/fe/fe_abf.cc b/deal.II/deal.II/source/fe/fe_abf.cc index c8eb6b843c..c16dc5a6fa 100644 --- a/deal.II/deal.II/source/fe/fe_abf.cc +++ b/deal.II/deal.II/source/fe/fe_abf.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2003, 2004, 2005 by the deal.II authors +// Copyright (C) 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 @@ -79,7 +79,7 @@ FE_ABF::FE_ABF (const unsigned int deg) this->restriction[i].reinit (n_dofs, n_dofs); } - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); initialize_restriction (); std::vector > diff --git a/deal.II/deal.II/source/fe/fe_dgp.cc b/deal.II/deal.II/source/fe/fe_dgp.cc index e2a987bf95..6736e464e0 100644 --- a/deal.II/deal.II/source/fe/fe_dgp.cc +++ b/deal.II/deal.II/source/fe/fe_dgp.cc @@ -30,12 +30,12 @@ FE_DGP::FE_DGP (const unsigned int degree) for (unsigned int i=0; i::children_per_cell; ++i) this->prolongation[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); // Fill restriction matrices with L2-projection for (unsigned int i=0; i::children_per_cell; ++i) this->restriction[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_projection_matrices (*this, &this->restriction[0]); + FETools::compute_projection_matrices (*this, this->restriction); } diff --git a/deal.II/deal.II/source/fe/fe_dgp_monomial.cc b/deal.II/deal.II/source/fe/fe_dgp_monomial.cc index dc192ea696..55c267f935 100644 --- a/deal.II/deal.II/source/fe/fe_dgp_monomial.cc +++ b/deal.II/deal.II/source/fe/fe_dgp_monomial.cc @@ -133,12 +133,12 @@ FE_DGPMonomial::FE_DGPMonomial (const unsigned int degree) for (unsigned int i=0; i::children_per_cell; ++i) this->prolongation[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); // Fill restriction matrices with L2-projection for (unsigned int i=0; i::children_per_cell; ++i) this->restriction[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_projection_matrices (*this, &this->restriction[0]); + FETools::compute_projection_matrices (*this, this->restriction); } diff --git a/deal.II/deal.II/source/fe/fe_dgq.cc b/deal.II/deal.II/source/fe/fe_dgq.cc index e14c20e524..05a9d8929b 100644 --- a/deal.II/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/deal.II/source/fe/fe_dgq.cc @@ -140,12 +140,12 @@ FE_DGQ::FE_DGQ (const unsigned int degree) for (unsigned int i=0; i::children_per_cell; ++i) this->prolongation[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); // Fill restriction matrices with L2-projection for (unsigned int i=0; i::children_per_cell; ++i) this->restriction[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_projection_matrices (*this, &this->restriction[0]); + FETools::compute_projection_matrices (*this, this->restriction); // finally fill in support points @@ -203,12 +203,12 @@ FE_DGQ::FE_DGQ (const Quadrature<1>& points) for (unsigned int i=0; i::children_per_cell; ++i) this->prolongation[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); // Fill restriction matrices with L2-projection for (unsigned int i=0; i::children_per_cell; ++i) this->restriction[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_projection_matrices (*this, &this->restriction[0]); + FETools::compute_projection_matrices (*this, this->restriction); // Compute support points, whivh // are the tensor product of the diff --git a/deal.II/deal.II/source/fe/fe_raviart_thomas.cc b/deal.II/deal.II/source/fe/fe_raviart_thomas.cc index 0190e4849e..788b303a18 100644 --- a/deal.II/deal.II/source/fe/fe_raviart_thomas.cc +++ b/deal.II/deal.II/source/fe/fe_raviart_thomas.cc @@ -69,17 +69,17 @@ FE_RaviartThomas::FE_RaviartThomas (const unsigned int deg) this->restriction[i].reinit (n_dofs, n_dofs); } - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); initialize_restriction(); - std::vector > - face_embeddings(1<<(dim-1), FullMatrix(this->dofs_per_face, - this->dofs_per_face)); - FETools::compute_face_embedding_matrices(*this, &face_embeddings[0], 0, 0); + FullMatrix face_embeddings[GeometryInfo::subfaces_per_face]; + for (unsigned int i=0; i::subfaces_per_face; ++i) + face_embeddings[i].reinit (this->dofs_per_face, this->dofs_per_face); + FETools::compute_face_embedding_matrices(*this, face_embeddings, 0, 0); this->interface_constraints.reinit((1<<(dim-1)) * this->dofs_per_face, this->dofs_per_face); unsigned int target_row=0; - for (unsigned int d=0;d::subfaces_per_face;++d) for (unsigned int i=0;i::FE_RaviartThomasNodal (const unsigned int deg) for (unsigned int i=0; i::children_per_cell; ++i) this->prolongation[i].reinit (this->dofs_per_cell, this->dofs_per_cell); - FETools::compute_embedding_matrices (*this, &this->prolongation[0]); + FETools::compute_embedding_matrices (*this, this->prolongation); - std::vector > - face_embeddings(1<<(dim-1), FullMatrix(this->dofs_per_face, - this->dofs_per_face)); - FETools::compute_face_embedding_matrices(*this, &face_embeddings[0], 0, 0); + FullMatrix face_embeddings[GeometryInfo::subfaces_per_face]; + for (unsigned int i=0; i::subfaces_per_face; ++i) + face_embeddings[i].reinit (this->dofs_per_face, this->dofs_per_face); + FETools::compute_face_embedding_matrices(*this, face_embeddings, 0, 0); this->interface_constraints.reinit((1<<(dim-1)) * this->dofs_per_face, this->dofs_per_face); unsigned int target_row=0; - for (unsigned int d=0;d::subfaces_per_face;++d) for (unsigned int i=0;i void FETools::compute_embedding_matrices(const FiniteElement& fe, - FullMatrix* matrices) + FullMatrix (&matrices)[GeometryInfo::children_per_cell]) { const unsigned int nc = GeometryInfo::children_per_cell; const unsigned int n = fe.dofs_per_cell; @@ -592,9 +592,9 @@ FETools::compute_embedding_matrices(const FiniteElement& fe, template void FETools::compute_face_embedding_matrices(const FiniteElement& fe, - FullMatrix* matrices, - unsigned int face_coarse, - unsigned int face_fine) + FullMatrix (&matrices)[GeometryInfo::subfaces_per_face], + const unsigned int face_coarse, + const unsigned int face_fine) { const unsigned int nc = GeometryInfo::subfaces_per_face; const unsigned int n = fe.dofs_per_face; @@ -607,7 +607,7 @@ FETools::compute_face_embedding_matrices(const FiniteElement& fe, Assert(matrices[i].m() == n, ExcDimensionMismatch(matrices[i].m(),n)); } - // Set up a meshes, one with a single + // Set up meshes, one with a single // reference cell and refine it once Triangulation tria; GridGenerator::hyper_cube (tria, 0, 1); @@ -711,7 +711,7 @@ FETools::compute_face_embedding_matrices(const FiniteElement& fe, typename Triangulation::active_cell_iterator fine_cell = tria.begin(0)->child(GeometryInfo::child_cell_on_face(face_coarse, - cell_number)); + cell_number)); fine.reinit(fine_cell); coarse.reinit(tria.begin(0)); @@ -758,7 +758,7 @@ FETools::compute_face_embedding_matrices(const FiniteElement& fe, template void FETools::compute_projection_matrices(const FiniteElement& fe, - FullMatrix* matrices) + FullMatrix (&matrices)[GeometryInfo::children_per_cell]) { const unsigned int nc = GeometryInfo::children_per_cell; const unsigned int n = fe.dofs_per_cell; @@ -1881,16 +1881,16 @@ void FETools::get_projection_matrix template void FETools::compute_embedding_matrices -(const FiniteElement &, FullMatrix*); +(const FiniteElement &, FullMatrix (&)[GeometryInfo::children_per_cell]); template void FETools::compute_face_embedding_matrices -(const FiniteElement &, FullMatrix*, +(const FiniteElement &, FullMatrix (&matrices)[GeometryInfo::subfaces_per_face], unsigned int, unsigned int); template void FETools::compute_projection_matrices -(const FiniteElement &, FullMatrix*); +(const FiniteElement &, FullMatrix (&)[GeometryInfo::children_per_cell]); template void FETools::interpolate diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 1ece207c65..ea11b73ace 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -794,6 +794,23 @@ inconvenience this causes.

deal.II

    +
  1. + Changed: The functions FETools::compute_embedding_matrices, + FETools::compute_face_embedding_matrices, and + FETools::compute_projection_matrices + (mostly used in internal computations in setting up finite + element objects) previously took pointers to the first element + of an array of matrices as arguments. This isn't type-safe, and + in particular did not allow to check for the number of matrices + in the array. The functions now take a reference to an array of + the correct length. +
    + (WB 2006/08/14) +

    +
  2. Extended: The VectorTools::project functions are now also implemented for 1d. -- 2.39.5