From: Sebastian Kinnewig Date: Thu, 18 Jan 2024 10:10:27 +0000 (+0100) Subject: Check if the global vector in distribute_local_to_global has ghost elements. X-Git-Tag: relicensing~146^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e24c00c88ae5dd37b2c4665ec84cdee767079567;p=dealii.git Check if the global vector in distribute_local_to_global has ghost elements. --- diff --git a/include/deal.II/lac/affine_constraints.h b/include/deal.II/lac/affine_constraints.h index 8dde2feb57..305a39828c 100644 --- a/include/deal.II/lac/affine_constraints.h +++ b/include/deal.II/lac/affine_constraints.h @@ -2629,6 +2629,7 @@ AffineConstraints::distribute_local_to_global( const std::vector &local_dof_indices, OutVector &global_vector) const { + Assert(global_vector.has_ghost_elements() == false, ExcGhostsPresent()); Assert(local_vector.size() == local_dof_indices.size(), ExcDimensionMismatch(local_vector.size(), local_dof_indices.size())); distribute_local_to_global(local_vector.begin(), diff --git a/include/deal.II/lac/affine_constraints.templates.h b/include/deal.II/lac/affine_constraints.templates.h index cb51e1acf9..e027a8d58e 100644 --- a/include/deal.II/lac/affine_constraints.templates.h +++ b/include/deal.II/lac/affine_constraints.templates.h @@ -2531,6 +2531,7 @@ AffineConstraints::distribute_local_to_global( AssertDimension(local_vector.size(), local_dof_indices_row.size()); AssertDimension(local_matrix.m(), local_dof_indices_row.size()); AssertDimension(local_matrix.n(), local_dof_indices_col.size()); + Assert(global_vector.has_ghost_elements() == false, ExcGhostsPresent()); // diagonal checks if we have only one index set (if both are equal // diagonal should be set to true). @@ -4210,6 +4211,7 @@ AffineConstraints::distribute_local_to_global( AssertDimension(local_matrix.n(), local_dof_indices.size()); AssertDimension(local_matrix.m(), local_dof_indices.size()); + Assert(global_vector.has_ghost_elements() == false, ExcGhostsPresent()); Assert(global_matrix.m() == global_matrix.n(), ExcNotQuadratic()); if (use_vectors == true) { @@ -4367,6 +4369,7 @@ AffineConstraints::distribute_local_to_global( AssertDimension(local_matrix.n(), local_dof_indices.size()); AssertDimension(local_matrix.m(), local_dof_indices.size()); + Assert(global_vector.has_ghost_elements() == false, ExcGhostsPresent()); Assert(global_matrix.m() == global_matrix.n(), ExcNotQuadratic()); Assert(global_matrix.n_block_rows() == global_matrix.n_block_cols(), ExcNotQuadratic());