]> https://gitweb.dealii.org/ - dealii.git/commitdiff
PETScWrappers::MatrixBase::add(): assert indices are in range.
authorDavid Wells <drwells@email.unc.edu>
Sun, 11 Aug 2024 19:12:16 +0000 (13:12 -0600)
committerDavid Wells <drwells@email.unc.edu>
Mon, 12 Aug 2024 02:42:37 +0000 (20:42 -0600)
include/deal.II/lac/petsc_matrix_base.h

index 660d3aa682054285555561a1c65df57228c36c38..599255ad2992a81be3f690f257c61356be1f229c 100644 (file)
@@ -1480,30 +1480,29 @@ namespace PETScWrappers
 
     prepare_action(VectorOperation::add);
 
-    const PetscInt  petsc_i = row;
-    const PetscInt *col_index_ptr;
-
-    const PetscScalar *col_value_ptr;
-    int                n_columns;
+    const auto petsc_row = row;
+    AssertIntegerConversion(petsc_row, row);
 
-    // If we don't elide zeros, the pointers are already available...
+    PetscInt n_columns = 0;
+    column_indices.resize(n_cols);
+    column_values.resize(n_cols);
     if (elide_zero_values == false)
       {
-        col_index_ptr = reinterpret_cast<const PetscInt *>(col_indices);
-        col_value_ptr = values;
-        n_columns     = n_cols;
+        n_columns = static_cast<PetscInt>(n_cols);
+        AssertIntegerConversion(n_columns, n_cols);
+
+        for (size_type j = 0; j < n_cols; ++j)
+          {
+            AssertIsFinite(values[j]);
+            column_indices[j] = static_cast<PetscInt>(col_indices[j]);
+            AssertIntegerConversion(column_indices[j], col_indices[j]);
+            column_values[j] = values[j];
+          }
       }
     else
       {
         // Otherwise, extract nonzero values in each row and get the
         // respective index.
-        if (column_indices.size() < n_cols)
-          {
-            column_indices.resize(n_cols);
-            column_values.resize(n_cols);
-          }
-
-        n_columns = 0;
         for (size_type j = 0; j < n_cols; ++j)
           {
             const PetscScalar value = values[j];
@@ -1516,13 +1515,15 @@ namespace PETScWrappers
               }
           }
         AssertIndexRange(n_columns, n_cols + 1);
-
-        col_index_ptr = column_indices.data();
-        col_value_ptr = column_values.data();
       }
 
-    const PetscErrorCode ierr = MatSetValues(
-      matrix, 1, &petsc_i, n_columns, col_index_ptr, col_value_ptr, ADD_VALUES);
+    const PetscErrorCode ierr = MatSetValues(matrix,
+                                             1,
+                                             &petsc_row,
+                                             n_columns,
+                                             column_indices.data(),
+                                             column_values.data(),
+                                             ADD_VALUES);
     AssertThrow(ierr == 0, ExcPETScError(ierr));
   }
 
@@ -1611,14 +1612,20 @@ namespace PETScWrappers
   inline bool
   MatrixBase::in_local_range(const size_type index) const
   {
-    PetscInt begin, end;
+    PetscInt petsc_begin, petsc_end;
 
     const PetscErrorCode ierr =
-      MatGetOwnershipRange(static_cast<const Mat &>(matrix), &begin, &end);
+      MatGetOwnershipRange(static_cast<const Mat &>(matrix),
+                           &petsc_begin,
+                           &petsc_end);
     AssertThrow(ierr == 0, ExcPETScError(ierr));
 
-    return ((index >= static_cast<size_type>(begin)) &&
-            (index < static_cast<size_type>(end)));
+    const auto begin = static_cast<size_type>(petsc_begin);
+    AssertIntegerConversion(begin, petsc_begin);
+    const auto end = static_cast<size_type>(petsc_end);
+    AssertIntegerConversion(end, petsc_end);
+
+    return ((index >= begin) && (index < end));
   }
 
 

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.