]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add functions to distribute local to global vectors and matrices.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 16 Aug 1998 12:25:18 +0000 (12:25 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 16 Aug 1998 12:25:18 +0000 (12:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@496 0785d39b-7218-0410-832d-ea1e28bc413d

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

index ae4f381f124800ffac4e205dc9b7410cc3ccc018..115b3e0c1f9323216dfa385ed5014ab91b859256 100644 (file)
@@ -8,6 +8,8 @@
 template <int dim> class Triangulation;
 template <int dim> class DoFHandler;
 class dVector;
+class dFMatrix;
+class dSMatrix;
 
 
 #include <grid/tria_accessor.h>
@@ -93,6 +95,10 @@ class DoFAccessor {
                                      * Exception
                                      */
     DeclException0 (ExcVectorDoesNotMatch);
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcMatrixDoesNotMatch);
 
   protected:
                                     /**
@@ -224,6 +230,32 @@ class DoFLineAccessor :  public DoFAccessor<dim>, public BaseClass {
                                      */
     TriaIterator<dim,DoFLineAccessor<dim,BaseClass> > child (const unsigned int) const;
 
+                                    /**
+                                     * Distribute a local (cell based) vector
+                                     * to a global one by mapping the local
+                                     * numbering of the degrees of freedom
+                                     * to the global one and entering the
+                                     * local values into the global vector.
+                                     *
+                                     * The elements are {\it added} up to
+                                     * the elements in the global vector,
+                                     * rather than just set, since this is
+                                     * usually what one wants.
+                                     */
+    void distribute_local_to_global (const dVector &local_source,
+                                    dVector       &global_destination) const;
+
+                                    /**
+                                     * This function does much the same as the
+                                     * #distribute_local_to_global(dVector,dVector)#
+                                     * function, but operates on matrices
+                                     * instead of vectors. The sparse matrix
+                                     * is supposed to have non-zero entry
+                                     * slots where required.
+                                     */
+    void distribute_local_to_global (const dFMatrix &local_source,
+                                    dSMatrix       &global_destination) const;
+    
                                     /**
                                      * Implement the copy operator needed
                                      * for the iterator classes.
@@ -317,6 +349,32 @@ class DoFQuadAccessor :  public DoFAccessor<dim>, public BaseClass {
                                      */
     TriaIterator<dim,DoFQuadAccessor<dim,BaseClass> > child (const unsigned int) const;
 
+                                    /**
+                                     * Distribute a local (cell based) vector
+                                     * to a global one by mapping the local
+                                     * numbering of the degrees of freedom
+                                     * to the global one and entering the
+                                     * local values into the global vector.
+                                     *
+                                     * The elements are {\it added} up to
+                                     * the elements in the global vector,
+                                     * rather than just set, since this is
+                                     * usually what one wants.
+                                     */
+    void distribute_local_to_global (const dVector &local_source,
+                                    dVector       &global_destination) const;
+
+                                    /**
+                                     * This function does much the same as the
+                                     * #distribute_local_to_global(dVector,dVector)#
+                                     * function, but operates on matrices
+                                     * instead of vectors. The sparse matrix
+                                     * is supposed to have non-zero entry
+                                     * slots where required.
+                                     */
+    void distribute_local_to_global (const dFMatrix &local_source,
+                                    dSMatrix       &global_destination) const;
+    
                                     /**
                                      * Implement the copy operator needed
                                      * for the iterator classes.
index 06c09029348af152f4777130dcd95a038b5bdc5e..49d45b2c6b7407b0f749fd52ef82ee433bbbc86c 100644 (file)
@@ -6,6 +6,11 @@
 #include <grid/tria_iterator.h>
 #include <grid/tria_iterator.templates.h>
 #include <fe/fe.h>
+
+#include <lac/dvector.h>
+#include <lac/dfmatrix.h>
+#include <lac/dsmatrix.h>
+
 #include <vector>
 
 
@@ -122,6 +127,61 @@ DoFLineAccessor<dim,BaseClass>::child (const unsigned int i) const {
 
 
 
+template <int dim, typename BaseClass>
+void DoFLineAccessor<dim,BaseClass>::
+distribute_local_to_global (const dVector &local_source,
+                           dVector       &global_destination) const {
+  Assert (dof_handler != 0, ExcInvalidObject());
+  Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+  Assert (local_source.size() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
+                                 dof_handler->get_selected_fe().dofs_per_line),
+         ExcVectorDoesNotMatch());
+  Assert (dof_handler->n_dofs() == global_destination.size(),
+         ExcVectorDoesNotMatch());
+
+  const unsigned int n_dofs = local_source.size();
+  
+                                  // get indices of dofs
+  vector<int> dofs (n_dofs);
+  get_dof_indices (dofs);
+  
+                                  // distribute cell vector
+  for (unsigned int j=0; j<n_dofs; ++j)
+    global_destination(dofs[j]) += local_source(j);
+};
+
+
+
+template <int dim, typename BaseClass>
+void DoFLineAccessor<dim,BaseClass>::
+distribute_local_to_global (const dFMatrix &local_source,
+                           dSMatrix       &global_destination) const {
+  Assert (dof_handler != 0, ExcInvalidObject());
+  Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+  Assert (local_source.m() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
+                              dof_handler->get_selected_fe().dofs_per_line),
+         ExcVectorDoesNotMatch());
+  Assert (local_source.m() == local_source.n(),
+         ExcMatrixDoesNotMatch());
+  Assert (dof_handler->n_dofs() == global_destination.m(),
+         ExcMatrixDoesNotMatch());
+  Assert (global_destination.m() == global_destination.n(),
+         ExcMatrixDoesNotMatch());
+
+  const unsigned int n_dofs = local_source.m();
+
+                                  // get indices of dofs
+  vector<int> dofs (n_dofs);
+  get_dof_indices (dofs);
+  
+                                  // distribute cell matrix
+  for (unsigned int i=0; i<n_dofs; ++i)
+    for (unsigned int j=0; j<n_dofs; ++j)
+      global_destination.add(dofs[i], dofs[j], local_source(i,j));
+};
+
+
+
 template <int dim, typename BaseClass>
 void DoFLineAccessor<dim,BaseClass>::copy_from (const DoFLineAccessor<dim,BaseClass> &a) {
   BaseClass::copy_from (a);
@@ -258,6 +318,63 @@ DoFQuadAccessor<dim,BaseClass>::child (const unsigned int i) const {
 
 
 
+template <int dim, typename BaseClass>
+void DoFQuadAccessor<dim,BaseClass>::
+distribute_local_to_global (const dVector &local_source,
+                           dVector       &global_destination) const {
+  Assert (dof_handler != 0, ExcInvalidObject());
+  Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+  Assert (local_source.size() == (4*dof_handler->get_selected_fe().dofs_per_vertex +
+                                 4*dof_handler->get_selected_fe().dofs_per_line +
+                                 dof_handler->get_selected_fe().dofs_per_quad),
+         ExcVectorDoesNotMatch());
+  Assert (dof_handler->n_dofs() == global_destination.size(),
+         ExcVectorDoesNotMatch());
+
+  const unsigned int n_dofs = local_source.size();
+
+                                  // get indices of dofs
+  vector<int> dofs (n_dofs);
+  get_dof_indices (dofs);
+  
+                                  // distribute cell vector
+  for (unsigned int j=0; j<n_dofs; ++j)
+    global_destination(dofs[j]) += local_source(j);
+};
+
+
+
+template <int dim, typename BaseClass>
+void DoFQuadAccessor<dim,BaseClass>::
+distribute_local_to_global (const dFMatrix &local_source,
+                           dSMatrix       &global_destination) const {
+  Assert (dof_handler != 0, ExcInvalidObject());
+  Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+  Assert (local_source.m() == (4*dof_handler->get_selected_fe().dofs_per_vertex +
+                              4*dof_handler->get_selected_fe().dofs_per_line +
+                              dof_handler->get_selected_fe().dofs_per_quad),
+         ExcMatrixDoesNotMatch());
+  Assert (local_source.m() == local_source.n(),
+         ExcMatrixDoesNotMatch());
+  Assert (dof_handler->n_dofs() == global_destination.m(),
+         ExcMatrixDoesNotMatch());
+  Assert (global_destination.m() == global_destination.n(),
+         ExcMatrixDoesNotMatch());
+  
+  const unsigned int n_dofs = local_source.m();
+
+                                  // get indices of dofs
+  vector<int> dofs (n_dofs);
+  get_dof_indices (dofs);
+  
+                                  // distribute cell matrix
+  for (unsigned int i=0; i<n_dofs; ++i)
+    for (unsigned int j=0; j<n_dofs; ++j)
+      global_destination.add(dofs[i], dofs[j], local_source(i,j));
+};
+
+
+
 template <int dim, typename BaseClass>
 void DoFQuadAccessor<dim,BaseClass>::copy_from (const DoFQuadAccessor<dim,BaseClass> &a) {
   BaseClass::copy_from (a);

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.