]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Make distribute_cell_to_dof_vector aware of finite elements made up of several compon...
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 7 May 1999 12:31:53 +0000 (12:31 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 7 May 1999 12:31:53 +0000 (12:31 +0000)
git-svn-id: https://svn.dealii.org/trunk@1299 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 6a1fd77c0d5d9826ad7b93a7fe44f07e46166581..f5c6a12d658e64721bd6e6e8af790c766b5e913a 100644 (file)
@@ -722,10 +722,27 @@ class DoFHandler
                                      * The output vector, being a data
                                      * vector on the grid, always consists
                                      * of elements of type #double#.
+                                     *
+                                     * In case the finite element used by
+                                     * this DoFHandler consists of more than
+                                     * one component, you should give which
+                                     * component in the output vector should
+                                     * be used to store the finite element
+                                     * field in; the default is zero (no other
+                                     * value is allowed if the finite element
+                                     * consists only of one component). All
+                                     * other components of the vector remain
+                                     * untouched, i.e. their contents are
+                                     * not changed.
+                                     *
+                                     * It is assumed that the output vector
+                                     * #dof_data# already has the right size,
+                                     * i.e. #n_dofs()# elements.
                                      */
     template <typename Number>
     void distribute_cell_to_dof_vector (const Vector<Number> &cell_data,
-                                       Vector<double>       &dof_data) const;
+                                       Vector<double>       &dof_data,
+                                       const unsigned int    component = 0) const;
 
                                     /**
                                      * Create a mapping from degree of freedom
@@ -1425,7 +1442,16 @@ class DoFHandler
     DeclException2 (ExcInvalidMaskDimension,
                    int, int,
                    << "The dimension of the mask " << arg1 << " does not match"
-                   << " the number of components in the finite element object.");
+                   << " the number of components in the finite element object.");    
+                                    /**
+                                     * Exception
+                                     */
+    DeclException2 (ExcInvalidComponent,
+                   int, int,
+                   << "The component you gave (" << arg1 << ") "
+                   << "is invalid with respect to the number "
+                   << "of components in the finite element "
+                   << "(" << arg2 << ")");
     
   protected:
     
index fdb33a5dff89c6572615f652e5c0d406dc8117b4..295d96993f0b39251b2a5a52b8983fcad08b086e 100644 (file)
@@ -2553,12 +2553,21 @@ unsigned int DoFHandler<3>::max_couplings_between_boundary_dofs () const {
 template <int dim>
 template <typename Number>
 void DoFHandler<dim>::distribute_cell_to_dof_vector (const Vector<Number> &cell_data,
-                                                    Vector<double>       &dof_data) const {
+                                                    Vector<double>       &dof_data,
+                                                    const unsigned int    component) 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());
+  Assert (dof_data.size()==n_dofs(), ExcWrongSize (dof_data.size(), n_dofs()));
+  Assert (component < selected_fe->n_components,
+         ExcInvalidComponent(component, selected_fe->n_components));
+
+                                  // store a flag whether we should care
+                                  // about different components. this is
+                                  // just a simplification, we could ask
+                                  // for this at every single place
+                                  // equally well
+  const bool consider_components = (selected_fe->n_components != 1);
+  
                                   // count how often we have added a value
                                   // in the sum for each dof
   vector<unsigned char> touch_count (n_dofs(), 0);
@@ -2568,28 +2577,40 @@ void DoFHandler<dim>::distribute_cell_to_dof_vector (const Vector<Number> &cell_
   unsigned int present_cell = 0;
   const unsigned int total_dofs = selected_fe->total_dofs;
   vector<int> dof_indices (total_dofs);
-  
+
   for (; cell!=endc; ++cell, ++present_cell) 
     {
       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]];
-       };
+                                        // consider this dof only if it
+                                        // is the right component. if there
+                                        // is only one component, short cut
+                                        // the test
+       if (!consider_components ||
+           (selected_fe->system_to_component_index(i).first == component))
+         {
+                                            // 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];
+                                      // assert that each dof was used
+                                      // at least once. this needs not be
+                                      // the case if the vector has more than
+                                      // one component
+      Assert (consider_components || (touch_count[i]!=0),
+             ExcInternalError());
+      if (touch_count[i] != 0)
+       dof_data(i) /=  touch_count[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.