]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Change a few function parameters from pointer to FullMatrix to reference to array...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 14 Aug 2006 23:25:33 +0000 (23:25 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 14 Aug 2006 23:25:33 +0000 (23:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@13705 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/source/fe/fe_abf.cc
deal.II/deal.II/source/fe/fe_dgp.cc
deal.II/deal.II/source/fe/fe_dgp_monomial.cc
deal.II/deal.II/source/fe/fe_dgq.cc
deal.II/deal.II/source/fe/fe_raviart_thomas.cc
deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc
deal.II/deal.II/source/fe/fe_tools.cc
deal.II/doc/news/changes.html

index a97f56d041783b4f34c4cca0bab4ddf586bceb6f..faca7d90de09d2f8c8fdb9d8d0f9b62b85e35eee 100644 (file)
@@ -264,15 +264,16 @@ class FETools
                                      * class for which we compute the
                                      * embedding matrices.
                                      * @param matrices A pointer to
-                                     * <i>2<sup>dim</sup></i> FullMatrix
+                                     * <i>GeometryInfo::children_per_cell=2<sup>dim</sup></i> FullMatrix
                                      * objects. This is the format
                                      * used in FiniteElement,
                                      * where we want to use ths
                                      * function mostly.
                                      */
     template <int dim, typename number>
-    static void compute_embedding_matrices(const FiniteElement<dim> &fe,
-                                          FullMatrix<number>* matrices);
+    static void
+    compute_embedding_matrices(const FiniteElement<dim> &fe,
+                              FullMatrix<number> (&matrices)[GeometryInfo<dim>::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
-                                     * <i>2<sup>dim-1</sup></i> 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
+                                     * <i>GeometryInfo<dim>::subfaces_per_face
+                                     * = 2<sup>dim-1</sup></i>
+                                     * 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<int dim, typename number>
-    static void compute_face_embedding_matrices(const FiniteElement<dim>& fe,
-                                               FullMatrix<number>* matrices,
-                                               unsigned int face_coarse,
-                                               unsigned int face_fine);
+    static void
+    compute_face_embedding_matrices(const FiniteElement<dim>& fe,
+                                   FullMatrix<number> (&matrices)[GeometryInfo<dim>::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 <int dim, typename number>
-    static void compute_projection_matrices(const FiniteElement<dim> &fe,
-                                           FullMatrix<number>* matrices);
+    static void
+    compute_projection_matrices(const FiniteElement<dim> &fe,
+                               FullMatrix<number> (&matrices)[GeometryInfo<dim>::children_per_cell]);
 
 //TODO:[WB] Replace this documentation by something comprehensible
     
index c8eb6b843cbef17d7072fa9ba78603ef0b58309c..c16dc5a6fa1c4be85d1d3ff109b5c5d46fb9750c 100644 (file)
@@ -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<dim>::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<FullMatrix<double> >
index e2a987bf9596ff1efe893255b11106f24136627d..6736e464e0a0bb5c48bb39e5caee890da98b7c6e 100644 (file)
@@ -30,12 +30,12 @@ FE_DGP<dim>::FE_DGP (const unsigned int degree)
   for (unsigned int i=0; i<GeometryInfo<dim>::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<GeometryInfo<dim>::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);
 }
 
 
index dc192ea69609121fe2ea507895bb6f7ce0c77b0a..55c267f9352f1c8525fe4e80bb3ad2f31c1ff88d 100644 (file)
@@ -133,12 +133,12 @@ FE_DGPMonomial<dim>::FE_DGPMonomial (const unsigned int degree)
   for (unsigned int i=0; i<GeometryInfo<dim>::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<GeometryInfo<dim>::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);
 }
 
 
index e14c20e5240b2110c277390e424b2666f5e37735..05a9d8929ba7e415024043d501ee063bbc6a2cf5 100644 (file)
@@ -140,12 +140,12 @@ FE_DGQ<dim>::FE_DGQ (const unsigned int degree)
   for (unsigned int i=0; i<GeometryInfo<dim>::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<GeometryInfo<dim>::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<dim>::FE_DGQ (const Quadrature<1>& points)
   for (unsigned int i=0; i<GeometryInfo<dim>::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<GeometryInfo<dim>::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
index 0190e4849e401eaebe2dab7fecf53967220bc580..788b303a18b23e618234a9aff9e9e5214e38ca17 100644 (file)
@@ -69,17 +69,17 @@ FE_RaviartThomas<dim>::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<FullMatrix<double> >
-    face_embeddings(1<<(dim-1), FullMatrix<double>(this->dofs_per_face,
-                                                  this->dofs_per_face));
-  FETools::compute_face_embedding_matrices(*this, &face_embeddings[0], 0, 0);
+  FullMatrix<double> face_embeddings[GeometryInfo<dim>::subfaces_per_face];
+  for (unsigned int i=0; i<GeometryInfo<dim>::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<face_embeddings.size();++d)
+  for (unsigned int d=0;d<GeometryInfo<dim>::subfaces_per_face;++d)
     for (unsigned int i=0;i<face_embeddings[d].m();++i)
       {
        for (unsigned int j=0;j<face_embeddings[d].n();++j)
index acf927e5d713045cbe434ffecff6965b9d784c10..c27861e07478622b9bbe65ca84ff5da23b88fb5d 100644 (file)
@@ -68,16 +68,16 @@ FE_RaviartThomasNodal<dim>::FE_RaviartThomasNodal (const unsigned int deg)
   for (unsigned int i=0; i<GeometryInfo<dim>::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<FullMatrix<double> >
-    face_embeddings(1<<(dim-1), FullMatrix<double>(this->dofs_per_face,
-                                                  this->dofs_per_face));
-  FETools::compute_face_embedding_matrices(*this, &face_embeddings[0], 0, 0);
+  FullMatrix<double> face_embeddings[GeometryInfo<dim>::subfaces_per_face];
+  for (unsigned int i=0; i<GeometryInfo<dim>::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<face_embeddings.size();++d)
+  for (unsigned int d=0;d<GeometryInfo<dim>::subfaces_per_face;++d)
     for (unsigned int i=0;i<face_embeddings[d].m();++i)
       {
        for (unsigned int j=0;j<face_embeddings[d].n();++j)
index 5091d1aa75381cfb9c1df3102d4546ef01904e33..2baa4c5ae621b5a34afd4fefa79b183de1e14306 100644 (file)
@@ -477,7 +477,7 @@ FETools::compute_node_matrix(
 template<int dim, typename number>
 void
 FETools::compute_embedding_matrices(const FiniteElement<dim>& fe,
-                                   FullMatrix<number>* matrices)
+                                   FullMatrix<number> (&matrices)[GeometryInfo<dim>::children_per_cell])
 {
   const unsigned int nc = GeometryInfo<dim>::children_per_cell;
   const unsigned int n  = fe.dofs_per_cell;
@@ -592,9 +592,9 @@ FETools::compute_embedding_matrices(const FiniteElement<dim>& fe,
 template<int dim, typename number>
 void
 FETools::compute_face_embedding_matrices(const FiniteElement<dim>& fe,
-                                        FullMatrix<number>* matrices,
-                                        unsigned int face_coarse,
-                                        unsigned int face_fine)
+                                        FullMatrix<number> (&matrices)[GeometryInfo<dim>::subfaces_per_face],
+                                        const unsigned int face_coarse,
+                                        const unsigned int face_fine)
 {
   const unsigned int nc = GeometryInfo<dim>::subfaces_per_face;
   const unsigned int n  = fe.dofs_per_face;
@@ -607,7 +607,7 @@ FETools::compute_face_embedding_matrices(const FiniteElement<dim>& fe,
       Assert(matrices[i].m() == n, ExcDimensionMismatch(matrices[i].m(),n));
     }
   
-                                   // Set up meshes, one with a single
+                                   // Set up meshes, one with a single
                                    // reference cell and refine it once
   Triangulation<dim> tria;
   GridGenerator::hyper_cube (tria, 0, 1);
@@ -711,7 +711,7 @@ FETools::compute_face_embedding_matrices(const FiniteElement<dim>& fe,
       
       typename Triangulation<dim>::active_cell_iterator fine_cell
        = tria.begin(0)->child(GeometryInfo<dim>::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<dim>& fe,
 template<int dim, typename number>
 void
 FETools::compute_projection_matrices(const FiniteElement<dim>& fe,
-                                    FullMatrix<number>* matrices)
+                                    FullMatrix<number> (&matrices)[GeometryInfo<dim>::children_per_cell])
 {
   const unsigned int nc = GeometryInfo<dim>::children_per_cell;
   const unsigned int n  = fe.dofs_per_cell;
@@ -1881,16 +1881,16 @@ void FETools::get_projection_matrix<deal_II_dimension>
 
 template
 void FETools::compute_embedding_matrices<deal_II_dimension>
-(const FiniteElement<deal_II_dimension> &, FullMatrix<double>*);
+(const FiniteElement<deal_II_dimension> &, FullMatrix<double> (&)[GeometryInfo<deal_II_dimension>::children_per_cell]);
 
 template
 void FETools::compute_face_embedding_matrices<deal_II_dimension>
-(const FiniteElement<deal_II_dimension> &, FullMatrix<double>*,
+(const FiniteElement<deal_II_dimension> &, FullMatrix<double> (&matrices)[GeometryInfo<deal_II_dimension>::subfaces_per_face],
  unsigned int, unsigned int);
 
 template
 void FETools::compute_projection_matrices<deal_II_dimension>
-(const FiniteElement<deal_II_dimension> &, FullMatrix<double>*);
+(const FiniteElement<deal_II_dimension> &, FullMatrix<double> (&)[GeometryInfo<deal_II_dimension>::children_per_cell]);
 
 template
 void FETools::interpolate<deal_II_dimension>
index 1ece207c65f7cda0164cf9db642d543b24bbff19..ea11b73ace7d875b6000c449aeba1aef12267c2d 100644 (file)
@@ -794,6 +794,23 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       Changed: The functions <code
+       class="member">FETools::compute_embedding_matrices</code>,
+       <code
+       class="member">FETools::compute_face_embedding_matrices</code>, and
+       <code
+       class="member">FETools::compute_projection_matrices</code>
+       (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.
+       <br> 
+       (WB 2006/08/14)
+       </p>
+
   <li> <p>
        Extended: The <code class="member">VectorTools::project</code> functions
        are now also implemented for 1d.

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.