From 0d5d785594a2051fd3154e5beb6cbf6c86794240 Mon Sep 17 00:00:00 2001 From: guido Date: Fri, 28 Nov 2003 15:49:37 +0000 Subject: [PATCH] L2-projection between degrees git-svn-id: https://svn.dealii.org/trunk@8213 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_tools.h | 29 +++++++++++++- deal.II/deal.II/source/fe/fe_tools.cc | 58 +++++++++++++++++++++++++++ deal.II/doc/news/c-4-0.html | 9 +++++ 3 files changed, 94 insertions(+), 2 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h index 5e87a72e21..1ab3d595b3 100644 --- a/deal.II/deal.II/include/fe/fe_tools.h +++ b/deal.II/deal.II/include/fe/fe_tools.h @@ -112,6 +112,11 @@ class FETools const FiniteElement &fe2, FullMatrix &difference_matrix); + /** + * Compute the local + * $L^2$-projection matrix from + * fe1 to fe2. + */ template static void get_projection_matrix(const FiniteElement &fe1, const FiniteElement &fe2, @@ -334,8 +339,28 @@ class FETools const InVector& z1, const DoFHandler& dof2, const ConstraintMatrix& constraints2, - OutVector& z1_difference); - + OutVector& z1_difference); + + /** + * $L^2$ projection for + * discontinuous + * elements. Operates the same + * direction as interpolate. + * + * The global projection can be + * computed by local matrices if + * the finite element spaces are + * discontinuous. With continuous + * elements, this is impossible, + * since a global mass matrix + * must be inverted. + */ + template + static void project_dg (const DoFHandler& dof1, + const InVector& u1, + const DoFHandler& dof2, + OutVector& u2); + /** * Gives the patchwise * extrapolation of a @p{dof1} diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 2ac088be83..f491bcb470 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -650,6 +650,48 @@ void FETools::interpolation_difference(const DoFHandler &dof1, } +template +void FETools::project_dg(const DoFHandler &dof1, + const InVector &u1, + const DoFHandler &dof2, + OutVector &u2) +{ + Assert(&dof1.get_tria()==&dof2.get_tria(), ExcTriangulationMismatch()); + Assert(dof1.get_fe().n_components() == dof2.get_fe().n_components(), + ExcDimensionMismatch(dof1.get_fe().n_components(), dof2.get_fe().n_components())); + Assert(u1.size()==dof1.n_dofs(), ExcDimensionMismatch(u1.size(), dof1.n_dofs())); + Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs())); + + typename DoFHandler::active_cell_iterator cell1 = dof1.begin_active(); + typename DoFHandler::active_cell_iterator cell2 = dof2.begin_active(); + typename DoFHandler::active_cell_iterator end = dof2.end(); + + const unsigned int n1 = dof1.get_fe().dofs_per_cell; + const unsigned int n2 = dof2.get_fe().dofs_per_cell; + + Vector u1_local(n1); + Vector u2_local(n2); + std::vector dofs(n2); + + FullMatrix matrix(n2,n1); + get_projection_matrix(dof1.get_fe(), dof2.get_fe(), matrix); + + while (cell2 != end) + { + cell1->get_dof_values(u1, u1_local); + matrix.vmult(u2_local, u1_local); + cell2->get_dof_indices(dofs); + for (unsigned int i=0; i void FETools::extrapolate(const DoFHandler &dof1, const InVector &u1, @@ -1564,6 +1606,10 @@ void FETools::interpolation_difference const DoFHandler &, const ConstraintMatrix &, Vector &); template +void FETools::project_dg +(const DoFHandler &, const Vector &, + const DoFHandler &, Vector &); +template void FETools::extrapolate (const DoFHandler &, const Vector &, const DoFHandler &, Vector &); @@ -1604,6 +1650,10 @@ void FETools::interpolation_difference const DoFHandler &, const ConstraintMatrix &, Vector &); template +void FETools::project_dg +(const DoFHandler &, const Vector &, + const DoFHandler &, Vector &); +template void FETools::extrapolate (const DoFHandler &, const Vector &, const DoFHandler &, Vector &); @@ -1644,6 +1694,10 @@ void FETools::interpolation_difference const DoFHandler &, const ConstraintMatrix &, BlockVector &); template +void FETools::project_dg +(const DoFHandler &, const BlockVector &, + const DoFHandler &, BlockVector &); +template void FETools::extrapolate (const DoFHandler &, const BlockVector &, const DoFHandler &, BlockVector &); @@ -1693,6 +1747,10 @@ void FETools::interpolation_difference const DoFHandler &, const ConstraintMatrix &, BlockVector &); template +void FETools::project_dg +(const DoFHandler &, const BlockVector &, + const DoFHandler &, BlockVector &); +template void FETools::extrapolate (const DoFHandler &, const BlockVector &, const DoFHandler &, BlockVector &); diff --git a/deal.II/doc/news/c-4-0.html b/deal.II/doc/news/c-4-0.html index 22a5c5f8b9..012c4cb7cc 100644 --- a/deal.II/doc/news/c-4-0.html +++ b/deal.II/doc/news/c-4-0.html @@ -304,6 +304,15 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
    +
  1. + New: The function + FETools::project_dg performs + L2-projections between finite element spaces of different + degrees on the same mesh. +
    + (GK 2003/11/28) +

    +
  2. Improved: FiniteElementData has a function tensor_degree(), returning the degree of the -- 2.39.5