From: hartmann Date: Mon, 17 Jun 2002 14:21:54 +0000 (+0000) Subject: New extrapolate function taking into account hanging node constraints. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5138c479081e3dc27ccc3476475a961d4ea7ff81;p=dealii-svn.git New extrapolate function taking into account hanging node constraints. git-svn-id: https://svn.dealii.org/trunk@6151 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h index c4cfebc61f..d3b0b534c1 100644 --- a/deal.II/deal.II/include/fe/fe_tools.h +++ b/deal.II/deal.II/include/fe/fe_tools.h @@ -341,18 +341,49 @@ class FETools * Note that the resulting field * does not satisfy continuity * requirements of the given - * finite elements. You have to - * apply - * @ref{ConstraintMatrix}@p{::distribute} - * with the hanging node - * constraints of the second DoF - * handler object to make the - * output continuous again. + * finite elements. + * + * When you use continuous + * elements on grids with hanging + * nodes, please use the + * @p{extrapolate} function with + * an additional + * @p{ConstraintMatrix} argument, + * see below. + */ + template + static void extrapolate(const DoFHandler &dof1, + const Vector &z1, + const DoFHandler &dof2, + Vector &z2); + + /** + * Gives the patchwise + * extrapolation of a @p{dof1} + * function @p{z1} to a @p{dof2} + * function @p{z2}. @p{dof1} and + * @p{dof2} need to be @ref{DoFHandler} + * based on the same triangulation. + * @p{constraints} is a hanging + * node constraints object + * corresponding to + * @p{dof2}. This object is + * particular important when + * interpolating onto continuous + * elements on grids with hanging + * nodes (locally refined grids). + * + * This function is interesting + * for e.g. extrapolating + * patchwise a piecewise linear + * solution to a piecewise + * quadratic solution. */ template static void extrapolate(const DoFHandler &dof1, const Vector &z1, const DoFHandler &dof2, + const ConstraintMatrix &constraints, Vector &z2); /** diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 5ad56f1cdb..3bbfb5b61a 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -375,6 +375,20 @@ void FETools::extrapolate(const DoFHandler &dof1, const Vector &u1, const DoFHandler &dof2, Vector &u2) +{ + ConstraintMatrix dummy; + dummy.close(); + extrapolate(dof1, u1, dof2, dummy, u2); +} + + + +template +void FETools::extrapolate(const DoFHandler &dof1, + const Vector &u1, + const DoFHandler &dof2, + const ConstraintMatrix &constraints, + Vector &u2) { Assert(dof1.get_fe().n_components() == dof2.get_fe().n_components(), ExcDimensionMismatch(dof1.get_fe().n_components(), dof2.get_fe().n_components())); @@ -383,7 +397,7 @@ void FETools::extrapolate(const DoFHandler &dof1, Assert(u2.size()==dof2.n_dofs(), ExcDimensionMismatch(u2.size(), dof2.n_dofs())); Vector u3(dof2.n_dofs()); - interpolate(dof1, u1, dof2, u3); + interpolate(dof1, u1, dof2, constraints, u3); const unsigned int dofs_per_cell = dof2.get_fe().dofs_per_cell; const Triangulation &tria=dof1.get_tria(); @@ -423,6 +437,9 @@ void FETools::extrapolate(const DoFHandler &dof1, } } } + + // Apply hanging node constraints. + constraints.distribute(u2); } diff --git a/deal.II/doc/news/2002/c-3-4.html b/deal.II/doc/news/2002/c-3-4.html index 04447c0734..c1d0f0c879 100644 --- a/deal.II/doc/news/2002/c-3-4.html +++ b/deal.II/doc/news/2002/c-3-4.html @@ -105,6 +105,20 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK

deal.II

    +
  1. New: The existing FETools::extrapolate functions does not + take into account hanging node constraints. Therefore, it works + for continuous elements on globally refined grids and on + discontinuous elements, only. Now, there is a new FETools::extrapolate function with an + additional ConstraintMatrix argument + for the hanging node constraints. This new function works for + continuous as well as for discontinuous elements. But, the old + function is still supported. +
    + (RH 2002/06/17) +

    +
  2. New: There are now Triangulation::save/load_user_pointers(vector<void *> &) functions similar to the respective