Functions

MatrixCreator Namespace Reference

Functions

template<int dim, typename number , int spacedim>
void create_mass_matrix (const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const hp::MappingCollection< dim, spacedim > &mapping, const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const hp::MappingCollection< dim, spacedim > &mapping, const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, typename number , int spacedim>
void create_mass_matrix (const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_boundary_mass_matrix (const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim-1 > &q, SparseMatrix< double > &matrix, const typename FunctionMap< spacedim >::type &boundary_functions, Vector< double > &rhs_vector, std::vector< unsigned int > &dof_to_boundary_mapping, const Function< spacedim > *const weight=0, std::vector< unsigned int > component_mapping=std::vector< unsigned int >())
template<int dim, int spacedim>
void create_boundary_mass_matrix (const DoFHandler< dim, spacedim > &dof, const Quadrature< dim-1 > &q, SparseMatrix< double > &matrix, const typename FunctionMap< spacedim >::type &boundary_functions, Vector< double > &rhs_vector, std::vector< unsigned int > &dof_to_boundary_mapping, const Function< spacedim > *const a=0, std::vector< unsigned int > component_mapping=std::vector< unsigned int >())
template<int dim, int spacedim>
void create_boundary_mass_matrix (const hp::MappingCollection< dim, spacedim > &mapping, const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim-1 > &q, SparseMatrix< double > &matrix, const typename FunctionMap< spacedim >::type &boundary_functions, Vector< double > &rhs_vector, std::vector< unsigned int > &dof_to_boundary_mapping, const Function< spacedim > *const a=0, std::vector< unsigned int > component_mapping=std::vector< unsigned int >())
template<int dim, int spacedim>
void create_boundary_mass_matrix (const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim-1 > &q, SparseMatrix< double > &matrix, const typename FunctionMap< spacedim >::type &boundary_functions, Vector< double > &rhs_vector, std::vector< unsigned int > &dof_to_boundary_mapping, const Function< spacedim > *const a=0, std::vector< unsigned int > component_mapping=std::vector< unsigned int >())
template<int dim, int spacedim>
void create_laplace_matrix (const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const hp::MappingCollection< dim, spacedim > &mapping, const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const hp::MappingCollection< dim, spacedim > &mapping, const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
template<int dim, int spacedim>
void create_laplace_matrix (const hp::DoFHandler< dim, spacedim > &dof, const hp::QCollection< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > &rhs, Vector< double > &rhs_vector, const Function< spacedim > *const a=0)
::ExceptionBaseExcComponentMismatch ()

Detailed Description

This namespace provides functions that assemble certain standard matrices for a given triangulation, using a given finite element, a given mapping and a quadrature formula.

Conventions for all functions

There exist two versions of each function. One with a Mapping argument and one without. If a code uses a mapping different from MappingQ1 the functions with mapping argument should be used. Code that uses only MappingQ1 may also use the functions without Mapping argument. Each of these latter functions create a MappingQ1 object and just call the respective functions with that object as mapping argument.

All functions take a sparse matrix object to hold the matrix to be created. The functions assume that the matrix is initialized with a sparsity pattern (SparsityPattern) corresponding to the given degree of freedom handler, i.e. the sparsity structure is already as needed. You can do this by calling the DoFTools::make_sparsity_pattern() function.

Furthermore it is assumed that no relevant data is in the matrix. Some entries will be overwritten and some others will contain invalid data if the matrix wasn't empty before. Therefore you may want to clear the matrix before assemblage.

All created matrices are `raw': they are not condensed, i.e. hanging nodes are not eliminated. The reason is that you may want to add several matrices and could then condense afterwards only once, instead of for every matrix. To actually do computations with these matrices, you have to condense the matrix using the ConstraintMatrix::condense function; you also have to condense the right hand side accordingly and distribute the solution afterwards.

If you want to use boundary conditions with the matrices generated by the functions of this class, you have to use a function like ProblemBase<>::apply_dirichlet_bc to matrix and right hand side.

Supported matrices

At present there are functions to create the following matrices:

Make sure that the order of the Quadrature formula given to these functions is sufficiently high to compute the matrices with the required accuracy. For the choice of this quadrature rule you need to take into account the polynomial degree of the FiniteElement basis functions, the roughness of the coefficient a, as well as the degree of the given Mapping (if any).

Note, that for system elements the mass matrix and the laplace matrix is implemented such that each components couple only with itself, i.e. there is no coupling of shape functions belonging to different components. If the degrees of freedom have been sorted according to their vector component (e.g., using DoFRenumbering::component_wise()), then the resulting matrices will be block diagonal.

If the finite element for which the mass matrix or the laplace matrix is to be built has more than one component, this function accepts a single coefficient as well as a vector valued coefficient function. For the latter case make sure that the number of components coincides with the number of components of the system finite element.

Matrices on the boundary

The create_boundary_mass_matrix() creates the matrix with entries $m_{ij} = \int_{\Gamma} \phi_i \phi_j dx$, where $\Gamma$ is the union of boundary parts with indicators contained in a FunctionMap passed to the function (i.e. if you want to set up the mass matrix for the parts of the boundary with indicators zero and 2, you pass the function a map of unsigned chars as parameter boundary_functions containing the keys zero and 2). The size of the matrix is equal to the number of degrees of freedom that have support on the boundary, i.e. it is not a matrix on all degrees of freedom, but only a subset. (The $\phi_i$ in the formula are this subsect of basis functions which have at least part of their support on $\Gamma$.) In order to determine which shape functions are to be considered, and in order to determine in which order, the function takes a dof_to_boundary_mapping; this object maps global DoF numbers to a numbering of the degrees of freedom located on the boundary, and can be obtained using the function DoFTools::map_dof_to_boundary_indices().

In order to work, the function needs a matrix of the correct size, built on top of a corresponding sparsity pattern. Since we only work on a subset of the degrees of freedom, we can't use the matrices and sparsity patterns that are created for the entire set of degrees of freedom. Rather, you should use the DoFHandler::make_boundary_sparsity_pattern() function to create the correct sparsity pattern, and build a matrix on top of it.

Note that at present there is no function that computes the mass matrix for all shape functions, though such a function would be trivial to implement.

Right hand sides

In many cases, you will not only want to build the matrix, but also a right hand side, which will give a vector with $f_i = \int_\Omega f(x) \phi_i(x) dx$. For this purpose, each function exists in two versions, one only building the matrix and one also building the right hand side vector. If you want to create a right hand side vector without creating a matrix, you can use the VectorTools::create_right_hand_side() function. The use of the latter may be useful if you want to create many right hand side vectors.

Author:
Wolfgang Bangerth, 1998, Ralf Hartmann, 2001

Function Documentation

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const Mapping< dim, spacedim > &  mapping,
const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > *const  a = 0 
)

Assemble the mass matrix. If no coefficient is given, it is assumed to be unity.

If the library is configured to use multithreading, this function works in parallel.

See the general doc of this class for more information.

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > *const  a = 0 
)

Calls the create_mass_matrix() function, see above, with mapping=MappingQ1<dim>().

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const Mapping< dim, spacedim > &  mapping,
const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Assemble the mass matrix and a right hand side vector. If no coefficient is given, it is assumed to be unity.

If the library is configured to use multithreading, this function works in parallel.

See the general doc of this class for more information.

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Calls the create_mass_matrix() function, see above, with mapping=MappingQ1<dim>().

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const hp::MappingCollection< dim, spacedim > &  mapping,
const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > *const  a = 0 
)

Same function as above, but for hp objects.

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > *const  a = 0 
)

Same function as above, but for hp objects.

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const hp::MappingCollection< dim, spacedim > &  mapping,
const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Same function as above, but for hp objects.

template<int dim, typename number , int spacedim>
void MatrixCreator::create_mass_matrix ( const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< number > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Same function as above, but for hp objects.

template<int dim, int spacedim>
void MatrixCreator::create_boundary_mass_matrix ( const Mapping< dim, spacedim > &  mapping,
const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim-1 > &  q,
SparseMatrix< double > &  matrix,
const typename FunctionMap< spacedim >::type &  boundary_functions,
Vector< double > &  rhs_vector,
std::vector< unsigned int > &  dof_to_boundary_mapping,
const Function< spacedim > *const  weight = 0,
std::vector< unsigned int >  component_mapping = std::vector< unsigned int >() 
)

Assemble the mass matrix and a right hand side vector along the boundary.

The matrix is assumed to already be initialized with a suiting sparsity pattern (the DoFHandler provides an appropriate function).

If the library is configured to use multithreading, this function works in parallel.

  • weight: an optional weight for the computation of the mass matrix. If no weight is given, it is set to one.
  • component_mapping: if the components in boundary_functions and dof do not coincide, this vector allows them to be remapped. If the vector is not empty, it has to have one entry for each component in dof. This entry is the component number in boundary_functions that should be used for this component in dof. By default, no remapping is applied.
Todo:
This function does not work for finite elements with cell-dependent shape functions.
template<int dim, int spacedim>
void MatrixCreator::create_boundary_mass_matrix ( const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim-1 > &  q,
SparseMatrix< double > &  matrix,
const typename FunctionMap< spacedim >::type &  boundary_functions,
Vector< double > &  rhs_vector,
std::vector< unsigned int > &  dof_to_boundary_mapping,
const Function< spacedim > *const  a = 0,
std::vector< unsigned int >  component_mapping = std::vector< unsigned int >() 
)

Calls the create_boundary_mass_matrix() function, see above, with mapping=MappingQ1<dim>().

template<int dim, int spacedim>
void MatrixCreator::create_boundary_mass_matrix ( const hp::MappingCollection< dim, spacedim > &  mapping,
const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim-1 > &  q,
SparseMatrix< double > &  matrix,
const typename FunctionMap< spacedim >::type &  boundary_functions,
Vector< double > &  rhs_vector,
std::vector< unsigned int > &  dof_to_boundary_mapping,
const Function< spacedim > *const  a = 0,
std::vector< unsigned int >  component_mapping = std::vector< unsigned int >() 
)

Same function as above, but for hp objects.

template<int dim, int spacedim>
void MatrixCreator::create_boundary_mass_matrix ( const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim-1 > &  q,
SparseMatrix< double > &  matrix,
const typename FunctionMap< spacedim >::type &  boundary_functions,
Vector< double > &  rhs_vector,
std::vector< unsigned int > &  dof_to_boundary_mapping,
const Function< spacedim > *const  a = 0,
std::vector< unsigned int >  component_mapping = std::vector< unsigned int >() 
)

Same function as above, but for hp objects. Same function as above, but for hp objects.

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const Mapping< dim, spacedim > &  mapping,
const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > *const  a = 0 
)

Assemble the Laplace matrix. If no coefficient is given, it is assumed to be constant one.

If the library is configured to use multithreading, this function works in parallel.

See the general doc of this class for more information.

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > *const  a = 0 
)

Calls the create_laplace_matrix() function, see above, with mapping=MappingQ1<dim>().

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const Mapping< dim, spacedim > &  mapping,
const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Assemble the Laplace matrix and a right hand side vector. If no coefficient is given, it is assumed to be constant one.

If the library is configured to use multithreading, this function works in parallel.

See the general doc of this class for more information.

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const DoFHandler< dim, spacedim > &  dof,
const Quadrature< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Calls the create_laplace_matrix() function, see above, with mapping=MappingQ1<dim>().

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const hp::MappingCollection< dim, spacedim > &  mapping,
const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > *const  a = 0 
)

Like the functions above, but for hp dof handlers, mappings, and quadrature collections.

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > *const  a = 0 
)

Like the functions above, but for hp dof handlers, mappings, and quadrature collections.

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const hp::MappingCollection< dim, spacedim > &  mapping,
const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Like the functions above, but for hp dof handlers, mappings, and quadrature collections.

template<int dim, int spacedim>
void MatrixCreator::create_laplace_matrix ( const hp::DoFHandler< dim, spacedim > &  dof,
const hp::QCollection< dim > &  q,
SparseMatrix< double > &  matrix,
const Function< spacedim > &  rhs,
Vector< double > &  rhs_vector,
const Function< spacedim > *const  a = 0 
)

Like the functions above, but for hp dof handlers, mappings, and quadrature collections.

::ExceptionBase& MatrixCreator::ExcComponentMismatch ( ) [static]

Exception

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Wed May 23 2012 06:08:30 by doxygen 1.7.3