]> https://gitweb.dealii.org/ - dealii.git/commitdiff
PETScWrappers fix improper usage of WORLD and GetArray
authorStefano Zampini <stefano.zampini@gmail.com>
Thu, 10 Nov 2022 16:52:39 +0000 (17:52 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Tue, 20 Dec 2022 09:52:30 +0000 (10:52 +0100)
add a method  to get the Vec

include/deal.II/lac/la_parallel_block_vector.templates.h
include/deal.II/lac/petsc_vector_base.h
include/deal.II/lac/read_write_vector.templates.h
include/deal.II/lac/vector.templates.h
source/lac/petsc_parallel_vector.cc
source/lac/petsc_vector_base.cc

index 330a82adee3acc6d6120e6015261842c991059e6..d530d1819bc0f9fafe0c1fdbbeaff92ee5b6b5ab 100644 (file)
@@ -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<const Vec &>(petsc_vec.block(i)),
-                        &start_ptr);
+          const PetscScalar *start_ptr;
+          PetscErrorCode     ierr =
+            VecGetArrayRead(static_cast<const Vec &>(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<const Vec &>(petsc_vec.block(i)),
-                                 &start_ptr);
+          ierr =
+            VecRestoreArrayRead(static_cast<const Vec &>(petsc_vec.block(i)),
+                                &start_ptr);
           AssertThrow(ierr == 0, ExcPETScError(ierr));
 
           // spread ghost values between processes?
index 9fe7ba6e8b08b171f6ff4c33386a908ad61ba753..5e84c1e42e68585191dd3eeae4c1ac4543acd461 100644 (file)
@@ -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));
       }
   }
index d1b9d6c070666daf6b5d31af1fe114ee7b066489..ac63050ac94710c4278379c7300c8e6e9469f789 100644 (file)
@@ -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<const Vec &>(petsc_vec), &start_ptr);
+    const PetscScalar *start_ptr;
+    PetscErrorCode     ierr =
+      VecGetArrayRead(static_cast<const Vec &>(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<const Vec &>(petsc_vec), &start_ptr);
+    ierr = VecRestoreArrayRead(static_cast<const Vec &>(petsc_vec), &start_ptr);
     AssertThrow(ierr == 0, ExcPETScError(ierr));
   }
 #endif
index 7f5beb5f5bf73d20e9afb43974f5069f156596ae..79920ed76e060a2b76a801b38265c2654e4e685f 100644 (file)
@@ -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);
index 794c687542be2b24b8d2e74c6d13dc0de74d0061..61a5cdd68e5bf3a73b6b3363efade6e641b02586 100644 (file)
@@ -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());
index 71cf16a3daaa4eeb07392649fffb517756d42fe9..671a198973549de1e136db13955f8419e7a71b2c 100644 (file)
@@ -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
   {

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.