From a9d5dde493a97972c5081ca3e06160e597a08946 Mon Sep 17 00:00:00 2001 From: wolf Date: Tue, 3 Nov 1998 18:47:46 +0000 Subject: [PATCH] Implement the interpolation from child to mother cells. git-svn-id: https://svn.dealii.org/trunk@625 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/Todo | 14 ++++++ deal.II/deal.II/include/dofs/dof_accessor.h | 21 ++++++++ deal.II/deal.II/include/fe/fe.h | 6 +-- deal.II/deal.II/source/dofs/dof_accessor.cc | 53 +++++++++++++++++++++ 4 files changed, 91 insertions(+), 3 deletions(-) diff --git a/deal.II/deal.II/Todo b/deal.II/deal.II/Todo index 20797c0285..6b1563a766 100644 --- a/deal.II/deal.II/Todo +++ b/deal.II/deal.II/Todo @@ -156,3 +156,17 @@ Remove all fe& in vectors.h No support points for non-Lagrangian elements? Check for Langrange in interpolation? + + +Implement the random distortion in Triangulation for hanging nodes + also. Hanging nodes need to be reset to the correct mean value + at the end, which is simple for 2D but difficult for 3D. Maybe take + a look at how we get to the original location of the point in the + execute_refinement function and copy the relevant lines. + + +Review DoFCellAccessor::get_interpolated_dof_values: it may be + necessary to write not the non-zero values, but those for which + the corresponding row of the transfer (restriction) matrix + has nonzero entries. This also catches those zeroes which should + really be written but are not because the nodal value was zero. \ No newline at end of file diff --git a/deal.II/deal.II/include/dofs/dof_accessor.h b/deal.II/deal.II/include/dofs/dof_accessor.h index 428062c35c..05daf44445 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.h @@ -599,6 +599,27 @@ class DoFCellAccessor : public DoFSubstructAccessor { void get_dof_values (const dVector &values, dVector &dof_values) const; + /** + * Return the interpolation of the given + * finite element function to the present + * cell. In the simplest case, the cell + * is a terminal one, i.e. has no children; + * then, the returned value is the vector + * of nodal values on that cell. You could + * then as well get the desired values + * through the #get_dof_values# function + * In the other case, when the cell has + * children, we use the restriction + * matrices provided by the finite element + * class to compute the interpolation + * from the children to the present cell. + * + * It is assumed that the vector already + * has the right size beforehand. + */ + void get_interpolated_dof_values (const dVector &values, + dVector &interpolated_values) const; + /** * Exception diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index f66146930a..492b427b81 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -260,7 +260,7 @@ struct FiniteElementBase : * constraints then. */ const dFMatrix & constraints () const; - + /** * Comparison operator. We also check for * equality of the constraint matrix, @@ -330,7 +330,8 @@ struct FiniteElementBase : * * To compute the interpolation of a * finite element field to a cell, you - * may use the #interpolation# function. + * may use the #get_interpolated_dof_values# + * function of the #DoFCellAccessor# class. * See there for more information. * * Upon assembling the transfer matrix @@ -510,7 +511,6 @@ struct FiniteElementBase : * the assumption does not hold, then the distribution has to happen * with all nodes on the small and the large cells. This is not * implemented in the #DoFHandler# class as of now. - * \end{itemize} * * @author Wolfgang Bangerth, 1998 */ diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index dbe373b49a..3e0844099f 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -503,6 +503,59 @@ DoFCellAccessor<2>::get_dof_values (const dVector &values, +template +void +DoFCellAccessor::get_interpolated_dof_values (const dVector &values, + dVector &interpolated_values) const { + const unsigned int total_dofs = dof_handler->get_fe().total_dofs; + + Assert (dof_handler != 0, ExcInvalidObject()); + Assert (&dof_handler->get_fe() != 0, ExcInvalidObject()); + Assert (interpolated_values.size() == total_dofs, + ExcVectorDoesNotMatch()); + Assert (values.size() == dof_handler->n_dofs(), + ExcVectorDoesNotMatch()); + + if (!has_children()) + // if this cell has no children: simply + // return the exact values on this cell + get_dof_values (values, interpolated_values); + else + // otherwise clobber them from the children + { + dVector tmp1(total_dofs); + dVector tmp2(total_dofs); + + interpolated_values.clear (); + + for (unsigned int child=0; child::children_per_cell; + ++child) + { + // get the values from the present + // child, if necessary by + // interpolation itself + child(child)->get_interpolated_dof_values (values, + tmp1); + // interpolate these to the mother + // cell + dof_handler->get_fe()->restrict(child).vmult (tmp2, tmp1); + + // now write those entries in tmp2 + // which are != 0 into the output + // vector. Note that we may not + // add them up, since we would then + // end in adding up the contribution + // from nodes on boundaries of + // children more than once. + for (unsigned int i=0; i