* 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
*
* @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.
*
* 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
* 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
// $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
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> >
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);
}
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);
}
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
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
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)
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)
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;
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;
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<dim> tria;
GridGenerator::hyper_cube (tria, 0, 1);
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));
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;
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>
<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.