From: Wolfgang Bangerth Date: Sun, 16 Aug 1998 12:25:18 +0000 (+0000) Subject: Add functions to distribute local to global vectors and matrices. X-Git-Tag: v8.0.0~22758 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c732ee68dd51c210929f8a4758349b8476447f98;p=dealii.git Add functions to distribute local to global vectors and matrices. git-svn-id: https://svn.dealii.org/trunk@496 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_accessor.h b/deal.II/deal.II/include/dofs/dof_accessor.h index ae4f381f12..115b3e0c1f 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.h @@ -8,6 +8,8 @@ template class Triangulation; template class DoFHandler; class dVector; +class dFMatrix; +class dSMatrix; #include @@ -93,6 +95,10 @@ class DoFAccessor { * Exception */ DeclException0 (ExcVectorDoesNotMatch); + /** + * Exception + */ + DeclException0 (ExcMatrixDoesNotMatch); protected: /** @@ -224,6 +230,32 @@ class DoFLineAccessor : public DoFAccessor, public BaseClass { */ TriaIterator > 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, public BaseClass { */ TriaIterator > 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. diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index 06c0902934..49d45b2c6b 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -6,6 +6,11 @@ #include #include #include + +#include +#include +#include + #include @@ -122,6 +127,61 @@ DoFLineAccessor::child (const unsigned int i) const { +template +void DoFLineAccessor:: +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 dofs (n_dofs); + get_dof_indices (dofs); + + // distribute cell vector + for (unsigned int j=0; j +void DoFLineAccessor:: +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 dofs (n_dofs); + get_dof_indices (dofs); + + // distribute cell matrix + for (unsigned int i=0; i void DoFLineAccessor::copy_from (const DoFLineAccessor &a) { BaseClass::copy_from (a); @@ -258,6 +318,63 @@ DoFQuadAccessor::child (const unsigned int i) const { +template +void DoFQuadAccessor:: +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 dofs (n_dofs); + get_dof_indices (dofs); + + // distribute cell vector + for (unsigned int j=0; j +void DoFQuadAccessor:: +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 dofs (n_dofs); + get_dof_indices (dofs); + + // distribute cell matrix + for (unsigned int i=0; i void DoFQuadAccessor::copy_from (const DoFQuadAccessor &a) { BaseClass::copy_from (a);