From 4a2702a7ea8ec49ad1ee50d16ecb19353bcc087f Mon Sep 17 00:00:00 2001 From: wolf Date: Tue, 16 Mar 2004 20:27:37 +0000 Subject: [PATCH] Implement local_apply_boundary_values. git-svn-id: https://svn.dealii.org/trunk@8792 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/numerics/matrices.h | 63 +++++++++++++++-- .../numerics/matrices.all_dimensions.cc | 69 +++++++++++++++++++ 2 files changed, 128 insertions(+), 4 deletions(-) diff --git a/deal.II/deal.II/include/numerics/matrices.h b/deal.II/deal.II/include/numerics/matrices.h index e898ed0ab8..5b218d01d6 100644 --- a/deal.II/deal.II/include/numerics/matrices.h +++ b/deal.II/deal.II/include/numerics/matrices.h @@ -714,10 +714,35 @@ class MatrixTools : public MatrixCreator * the classes that are used to wrap * PETSc objects. * - * For a replacement function, - * see the documentation of the - * @ref{FilteredMatrix} class in - * the @p{LAC} sublibrary. + * Note that this function is not very + * efficient: it needs to alternatingly + * read and write into the matrix, a + * situation that PETSc does not handle + * too well. In addition, we only get rid + * of rows corresponding to boundary + * nodes, but the corresponding case of + * deleting the respective columns + * (i.e. if @arg eliminate_columns is @p + * true) is not presently implemented, + * and probably will never because it is + * too expensive without direct access to + * the PETSc data structures. A third + * reason against this function is that + * it doesn't handle the case where the + * matrix is distributed across an MPI + * system. + * + * In order to still be able to eliminate + * boundary values, it is better to get + * rid of them before the local matrices + * and vectors are distributed to the + * global ones, because then we don't + * have to mess with the sparse data + * structures. The + * local_apply_boundary_values() function + * does that, and is recommended for use + * instead of the global one for PETSc + * matrices and vectors. */ #ifdef DEAL_II_USE_PETSC static void @@ -727,6 +752,36 @@ class MatrixTools : public MatrixCreator PETScWrappers::VectorBase &right_hand_side, const bool eliminate_columns = true); #endif + + /** + * Rather than applying boundary values + * to the global matrix and vector after + * assembly, this function does so @em + * before assembling, by modifying the + * local matrix and vector + * contributions. If you call this + * function on all local contributions, + * the resulting matrix will have the + * same entries, and the final call to + * apply_boundary_values() on the global + * system will not be necessary. + * + * Since this function does not have to + * work on the complicated data + * structures of sparse matrices, it is + * relatively cheap. It may therefore be + * a win if you have many fixed degrees + * of freedom (e.g. boundary nodes), or + * if access to the sparse matrix is + * expensive (e.g. for block sparse + * matrices, or for PETSc matrices). + */ + static void + local_apply_boundary_values (const std::map &boundary_values, + const std::vector &local_dof_indices, + FullMatrix &local_matrix, + Vector &local_rhs, + const bool eliminate_columns); /** * Exception diff --git a/deal.II/deal.II/source/numerics/matrices.all_dimensions.cc b/deal.II/deal.II/source/numerics/matrices.all_dimensions.cc index 7834575749..3301c70663 100644 --- a/deal.II/deal.II/source/numerics/matrices.all_dimensions.cc +++ b/deal.II/deal.II/source/numerics/matrices.all_dimensions.cc @@ -13,6 +13,7 @@ #include +#include #include #include #include @@ -642,6 +643,74 @@ apply_boundary_values (const std::map &boundary_values, #endif +void +MatrixTools:: +local_apply_boundary_values (const std::map &boundary_values, + const std::vector &local_dof_indices, + FullMatrix &local_matrix, + Vector &local_rhs, + const bool eliminate_columns) +{ + Assert (local_dof_indices.size() == local_matrix.m(), + ExcDimensionMismatch(local_dof_indices.size(), + local_matrix.m())); + Assert (local_dof_indices.size() == local_matrix.n(), + ExcDimensionMismatch(local_dof_indices.size(), + local_matrix.n())); + Assert (local_dof_indices.size() == local_rhs.size(), + ExcDimensionMismatch(local_dof_indices.size(), + local_rhs.size())); + + // if there is nothing to do, then exit + // right away + if (boundary_values.size() == 0) + return; + + // otherwise traverse all the dofs used in + // the local matrices and vectors and see + // what's there to do + const unsigned int n_local_dofs = local_dof_indices.size(); + for (unsigned int i=0; i::const_iterator + boundary_value = boundary_values.find (local_dof_indices[i]); + if (boundary_value != boundary_values.end()) + { + // remove this row, except for the + // diagonal element + for (unsigned j=0; jsecond; + + // finally do the elimination step + // if requested + if (eliminate_columns == true) + { + for (unsigned int row=0; rowsecond; + local_matrix(row,i) = 0; + } + } + } + } +} + + -- 2.39.5