]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add new function to distribute cell to dof data.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 6 Apr 1998 17:31:31 +0000 (17:31 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 6 Apr 1998 17:31:31 +0000 (17:31 +0000)
git-svn-id: https://svn.dealii.org/trunk@147 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_handler.h
deal.II/deal.II/source/dofs/dof_handler.cc

index 1d39fb5de81bd7086a1f3c8bceeb4c564d7677aa..efc298267d4f0596f7aa9b0be8db9799a4796e23 100644 (file)
@@ -27,6 +27,7 @@ template <int dim, class Accessor> class TriaActiveIterator;
 
 template <int dim> class Triangulation;
 
+class dVector;
 class dSMatrix;
 class dSMatrixStruct;
 class ConstraintMatrix;
@@ -601,7 +602,24 @@ class DoFHandler : public DoFDimensionInfo<dim> {
                                      * the transfer process.
                                      */
     unsigned int max_transfer_entries (const unsigned int max_level_diff) const;
-    
+
+                                    /**
+                                     * Take a vector of values which live on
+                                     * cells (e.g. an error per cell) and
+                                     * distribute it to the dofs in such a
+                                     * way that a finite element field results,
+                                     * which can then be further processed,
+                                     * e.g. for output.
+                                     *
+                                     * It is assumed that the number of
+                                     * elements in #cell_data# equals the
+                                     * number of active cells. The size of
+                                     * #dof_data# is adjusted to the right
+                                     * size.
+                                     */
+    void distribute_cell_to_dof_vector (const vector<double> &cell_data,
+                                       dVector              &dof_data) const;
+
                                     /**
                                      *  @name Cell iterator functions
                                      */
@@ -1034,6 +1052,13 @@ class DoFHandler : public DoFDimensionInfo<dim> {
                                      *  Exception
                                      */
     DeclException0 (ExcFunctionNotUseful);
+                                    /**
+                                     * Exception
+                                     */
+    DeclException2 (ExcWrongSize,
+                   int, int,
+                   << "The dimension " << arg1 << " of the vector is wrong. "
+                   << "It should be " << arg2);
     
   protected:
                                     /**
index 9bcedc4f060cbf225ee803f9011d72f491b35dbc..fa51a436f33a0b3438c090a3e3f19c213b6a4e3c 100644 (file)
@@ -1465,6 +1465,53 @@ unsigned int DoFHandler<2>::max_transfer_entries (const unsigned int max_level_d
 
 
 
+template <int dim>
+void DoFHandler<dim>::distribute_cell_to_dof_vector (const vector<double> &cell_data,
+                                                    dVector              &dof_data) const {
+  Assert (cell_data.size()==tria->n_active_cells(),
+         ExcWrongSize (cell_data.size(), tria->n_active_cells()));
+
+                                  // assign the right size to the output vector
+  dof_data.reinit (n_dofs());
+                                  // count how often we have added a value
+                                  // in the sum for each dof
+  vector<unsigned char> touch_count (n_dofs(), 0);
+
+  active_cell_iterator cell = begin_active(),
+                      endc = end();
+  unsigned int present_cell = 0;
+  vector<int> dof_indices;
+  const unsigned int total_dofs = selected_fe->total_dofs;
+  
+  for (; cell!=endc; ++cell, ++present_cell) 
+    {
+      dof_indices.resize (0);
+      cell->get_dof_indices (dof_indices);
+      for (unsigned int i=0; i<total_dofs; ++i)
+       {
+                                          // sum up contribution of the
+                                          // present_cell to this dof
+         dof_data(dof_indices[i]) += cell_data[present_cell];
+                                          // note that we added another
+                                          // summand
+         ++touch_count[dof_indices[i]];
+       };
+    };
+
+                                  // compute the mean value on all the
+                                  // dofs by dividing with the number
+                                  // of summands.
+  for (unsigned int i=0; i<n_dofs(); ++i)
+    {
+      Assert (touch_count[i]!=0, ExcInternalError());
+      dof_data(i) /=  touch_count[i];
+    };
+};
+
+
+
+
+
 void DoFHandler<1>::reserve_space (const FiniteElement<1> &fe) {
   Assert (tria->n_levels() > 0, ExcInvalidTriangulation());
                                    // delete all levels and set them up

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.