]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement the interpolation from child to mother cells.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 3 Nov 1998 18:47:46 +0000 (18:47 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 3 Nov 1998 18:47:46 +0000 (18:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@625 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/Todo
deal.II/deal.II/include/dofs/dof_accessor.h
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/source/dofs/dof_accessor.cc

index 20797c02851d2033cb392e44389769688f6057f6..6b1563a7668ef04b01cde761eca25f246659df98 100644 (file)
@@ -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<dim>::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
index 428062c35cbf0659ebe55600696bc98312bf8d9a..05daf444451feae87fac8f3186ba6ffc8f94129f 100644 (file)
@@ -599,6 +599,27 @@ class DoFCellAccessor :  public DoFSubstructAccessor<dim> {
     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
index f66146930a2273750b21ee3b90f04ed0ac8c1524..492b427b814f266cd0e9266c6bc72342e98b7f04 100644 (file)
@@ -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
  */
index dbe373b49aeacbe0fa1edf94cfc5304b02fd85a6..3e0844099fd02acf0ef9953e41f3f55258408131 100644 (file)
@@ -503,6 +503,59 @@ DoFCellAccessor<2>::get_dof_values (const dVector &values,
 
 
 
+template <int dim>
+void
+DoFCellAccessor<dim>::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<GeometryInfo<dim>::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<total_dofs; ++i)
+           if (tmp2(i) != 0)
+             interpolated_values(i) = tmp2(i);
+    };
+};
+
+
+
+
 
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.