From 54d326cfd1e59de8fbb645b51a05548ff9cea1fd Mon Sep 17 00:00:00 2001 From: Stefano Zampini Date: Thu, 10 Nov 2022 17:52:39 +0100 Subject: [PATCH] PETScWrappers fix improper usage of WORLD and GetArray add a method to get the Vec --- .../lac/la_parallel_block_vector.templates.h | 13 ++--- include/deal.II/lac/petsc_vector_base.h | 20 +++++--- .../deal.II/lac/read_write_vector.templates.h | 8 +-- include/deal.II/lac/vector.templates.h | 6 +-- source/lac/petsc_parallel_vector.cc | 8 +-- source/lac/petsc_vector_base.cc | 49 +++++++++++-------- 6 files changed, 60 insertions(+), 44 deletions(-) diff --git a/include/deal.II/lac/la_parallel_block_vector.templates.h b/include/deal.II/lac/la_parallel_block_vector.templates.h index 330a82adee..d530d1819b 100644 --- a/include/deal.II/lac/la_parallel_block_vector.templates.h +++ b/include/deal.II/lac/la_parallel_block_vector.templates.h @@ -311,10 +311,10 @@ namespace LinearAlgebra StandardExceptions::ExcInvalidState()); // get a representation of the vector and copy it - PetscScalar * start_ptr; - PetscErrorCode ierr = - VecGetArray(static_cast(petsc_vec.block(i)), - &start_ptr); + const PetscScalar *start_ptr; + PetscErrorCode ierr = + VecGetArrayRead(static_cast(petsc_vec.block(i)), + &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); const size_type vec_size = this->block(i).locally_owned_size(); @@ -323,8 +323,9 @@ namespace LinearAlgebra this->block(i).begin()); // restore the representation of the vector - ierr = VecRestoreArray(static_cast(petsc_vec.block(i)), - &start_ptr); + ierr = + VecRestoreArrayRead(static_cast(petsc_vec.block(i)), + &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); // spread ghost values between processes? diff --git a/include/deal.II/lac/petsc_vector_base.h b/include/deal.II/lac/petsc_vector_base.h index 9fe7ba6e8b..5e84c1e42e 100644 --- a/include/deal.II/lac/petsc_vector_base.h +++ b/include/deal.II/lac/petsc_vector_base.h @@ -777,6 +777,14 @@ namespace PETScWrappers */ operator const Vec &() const; + /** + * Return a reference to the underlying PETSc type. It can be used to + * modify the underlying data, so use it only when you know what you + * are doing. + */ + Vec & + petsc_vector(); + /** * Estimate for the memory consumption (not implemented for this class). */ @@ -1197,8 +1205,8 @@ namespace PETScWrappers ierr = VecGetSize(locally_stored_elements, &lsize); AssertThrow(ierr == 0, ExcPETScError(ierr)); - PetscScalar *ptr; - ierr = VecGetArray(locally_stored_elements, &ptr); + const PetscScalar *ptr; + ierr = VecGetArrayRead(locally_stored_elements, &ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); for (PetscInt i = 0; i < n_idx; ++i) @@ -1221,7 +1229,7 @@ namespace PETScWrappers } } - ierr = VecRestoreArray(locally_stored_elements, &ptr); + ierr = VecRestoreArrayRead(locally_stored_elements, &ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); ierr = VecGhostRestoreLocalForm(vector, &locally_stored_elements); @@ -1236,8 +1244,8 @@ namespace PETScWrappers PetscErrorCode ierr = VecGetOwnershipRange(vector, &begin, &end); AssertThrow(ierr == 0, ExcPETScError(ierr)); - PetscScalar *ptr; - ierr = VecGetArray(vector, &ptr); + const PetscScalar *ptr; + ierr = VecGetArrayRead(vector, &ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); for (PetscInt i = 0; i < n_idx; ++i) @@ -1258,7 +1266,7 @@ namespace PETScWrappers *(values_begin + i) = *(ptr + index - begin); } - ierr = VecRestoreArray(vector, &ptr); + ierr = VecRestoreArrayRead(vector, &ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); } } diff --git a/include/deal.II/lac/read_write_vector.templates.h b/include/deal.II/lac/read_write_vector.templates.h index d1b9d6c070..ac63050ac9 100644 --- a/include/deal.II/lac/read_write_vector.templates.h +++ b/include/deal.II/lac/read_write_vector.templates.h @@ -557,16 +557,16 @@ namespace LinearAlgebra StandardExceptions::ExcInvalidState()); // get a representation of the vector and copy it - PetscScalar * start_ptr; - PetscErrorCode ierr = - VecGetArray(static_cast(petsc_vec), &start_ptr); + const PetscScalar *start_ptr; + PetscErrorCode ierr = + VecGetArrayRead(static_cast(petsc_vec), &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); const size_type vec_size = petsc_vec.locally_owned_size(); internal::copy_petsc_vector(start_ptr, start_ptr + vec_size, begin()); // restore the representation of the vector - ierr = VecRestoreArray(static_cast(petsc_vec), &start_ptr); + ierr = VecRestoreArrayRead(static_cast(petsc_vec), &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); } #endif diff --git a/include/deal.II/lac/vector.templates.h b/include/deal.II/lac/vector.templates.h index 7f5beb5f5b..79920ed76e 100644 --- a/include/deal.II/lac/vector.templates.h +++ b/include/deal.II/lac/vector.templates.h @@ -108,8 +108,8 @@ namespace internal scatter_context, v, sequential_vector, INSERT_VALUES, SCATTER_FORWARD); AssertThrow(ierr == 0, ExcPETScError(ierr)); - PetscScalar *start_ptr; - ierr = VecGetArray(sequential_vector, &start_ptr); + const PetscScalar *start_ptr; + ierr = VecGetArrayRead(sequential_vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); const PETScWrappers::VectorBase::size_type v_size = v.size(); @@ -119,7 +119,7 @@ namespace internal internal::VectorOperations::copy(start_ptr, start_ptr + out.size(), out.begin()); - ierr = VecRestoreArray(sequential_vector, &start_ptr); + ierr = VecRestoreArrayRead(sequential_vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); ierr = VecScatterDestroy(&scatter_context); diff --git a/source/lac/petsc_parallel_vector.cc b/source/lac/petsc_parallel_vector.cc index 794c687542..61a5cdd68e 100644 --- a/source/lac/petsc_parallel_vector.cc +++ b/source/lac/petsc_parallel_vector.cc @@ -345,10 +345,10 @@ namespace PETScWrappers // get a representation of the vector and // loop over all the elements - PetscScalar *val; - PetscInt nlocal, istart, iend; + const PetscScalar *val; + PetscInt nlocal, istart, iend; - PetscErrorCode ierr = VecGetArray(vector, &val); + PetscErrorCode ierr = VecGetArrayRead(vector, &val); AssertThrow(ierr == 0, ExcPETScError(ierr)); ierr = VecGetLocalSize(vector, &nlocal); @@ -404,7 +404,7 @@ namespace PETScWrappers // restore the representation of the // vector - ierr = VecRestoreArray(vector, &val); + ierr = VecRestoreArrayRead(vector, &val); AssertThrow(ierr == 0, ExcPETScError(ierr)); AssertThrow(out.fail() == false, ExcIO()); diff --git a/source/lac/petsc_vector_base.cc b/source/lac/petsc_vector_base.cc index 71cf16a3da..671a198973 100644 --- a/source/lac/petsc_vector_base.cc +++ b/source/lac/petsc_vector_base.cc @@ -58,8 +58,8 @@ namespace PETScWrappers ierr = VecGetSize(locally_stored_elements, &lsize); AssertThrow(ierr == 0, ExcPETScError(ierr)); - PetscScalar *ptr; - ierr = VecGetArray(locally_stored_elements, &ptr); + const PetscScalar *ptr; + ierr = VecGetArrayRead(locally_stored_elements, &ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); PetscScalar value; @@ -85,7 +85,7 @@ namespace PETScWrappers value = *(ptr + ghostidx + end - begin); } - ierr = VecRestoreArray(locally_stored_elements, &ptr); + ierr = VecRestoreArrayRead(locally_stored_elements, &ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); ierr = @@ -454,8 +454,8 @@ namespace PETScWrappers // get a representation of the vector and // loop over all the elements - PetscScalar * start_ptr; - PetscErrorCode ierr = VecGetArray(vector, &start_ptr); + const PetscScalar *start_ptr; + PetscErrorCode ierr = VecGetArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); PetscScalar mean = 0; @@ -483,7 +483,7 @@ namespace PETScWrappers // restore the representation of the // vector - ierr = VecRestoreArray(vector, &start_ptr); + ierr = VecRestoreArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); return mean; @@ -521,8 +521,8 @@ namespace PETScWrappers { // get a representation of the vector and // loop over all the elements - PetscScalar * start_ptr; - PetscErrorCode ierr = VecGetArray(vector, &start_ptr); + const PetscScalar *start_ptr; + PetscErrorCode ierr = VecGetArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); real_type norm = 0; @@ -550,7 +550,7 @@ namespace PETScWrappers // restore the representation of the // vector - ierr = VecRestoreArray(vector, &start_ptr); + ierr = VecRestoreArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); return norm; @@ -603,8 +603,8 @@ namespace PETScWrappers { // get a representation of the vector and // loop over all the elements - PetscScalar * start_ptr; - PetscErrorCode ierr = VecGetArray(vector, &start_ptr); + const PetscScalar *start_ptr; + PetscErrorCode ierr = VecGetArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); const PetscScalar *ptr = start_ptr, @@ -622,7 +622,7 @@ namespace PETScWrappers // restore the representation of the // vector - ierr = VecRestoreArray(vector, &start_ptr); + ierr = VecRestoreArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); return flag; @@ -657,8 +657,8 @@ namespace PETScWrappers { // get a representation of the vector and // loop over all the elements - PetscScalar * start_ptr; - PetscErrorCode ierr = VecGetArray(vector, &start_ptr); + const PetscScalar *start_ptr; + PetscErrorCode ierr = VecGetArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); const PetscScalar *ptr = start_ptr, @@ -676,7 +676,7 @@ namespace PETScWrappers // restore the representation of the // vector - ierr = VecRestoreArray(vector, &start_ptr); + ierr = VecRestoreArrayRead(vector, &start_ptr); AssertThrow(ierr == 0, ExcPETScError(ierr)); return flag; @@ -845,14 +845,15 @@ namespace PETScWrappers VectorBase::write_ascii(const PetscViewerFormat format) { // TODO[TH]:assert(is_compressed()) + MPI_Comm comm = PetscObjectComm((PetscObject)vector); // Set options PetscErrorCode ierr = - PetscViewerSetFormat(PETSC_VIEWER_STDOUT_WORLD, format); + PetscViewerSetFormat(PETSC_VIEWER_STDOUT_(comm), format); AssertThrow(ierr == 0, ExcPETScError(ierr)); // Write to screen - ierr = VecView(vector, PETSC_VIEWER_STDOUT_WORLD); + ierr = VecView(vector, PETSC_VIEWER_STDOUT_(comm)); AssertThrow(ierr == 0, ExcPETScError(ierr)); } @@ -868,8 +869,8 @@ namespace PETScWrappers // get a representation of the vector and // loop over all the elements - PetscScalar * val; - PetscErrorCode ierr = VecGetArray(vector, &val); + const PetscScalar *val; + PetscErrorCode ierr = VecGetArrayRead(vector, &val); AssertThrow(ierr == 0, ExcPETScError(ierr)); @@ -897,7 +898,7 @@ namespace PETScWrappers // restore the representation of the // vector - ierr = VecRestoreArray(vector, &val); + ierr = VecRestoreArrayRead(vector, &val); AssertThrow(ierr == 0, ExcPETScError(ierr)); AssertThrow(out.fail() == false, ExcIO()); @@ -913,13 +914,19 @@ namespace PETScWrappers } - VectorBase::operator const Vec &() const { return vector; } + Vec & + VectorBase::petsc_vector() + { + return vector; + } + + std::size_t VectorBase::memory_consumption() const { -- 2.39.5