Don't access non-local rows of a matrix.

git-svn-id: https://svn.dealii.org/trunk@8989 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/lac/include/lac/petsc_matrix_base.h b/deal.II/lac/include/lac/petsc_matrix_base.h
index 558d2bab2ae4ba23f08dafc675884def0ae095cb..2503272619be2d7bfb6f03a83fed8b2113497650 100644 (file)
--- a/deal.II/lac/include/lac/petsc_matrix_base.h
+++ b/deal.II/lac/include/lac/petsc_matrix_base.h
@@ -95,6 +95,22 @@ namespace PETScWrappers
                                               * Exception
                                               */
             DeclException0 (ExcBeyondEndOfMatrix);
+                                             /**
+                                              * Exception
+                                              */
+            DeclException1 (ExcPETScError,
+                            int,
+                            << "An error with error number " << arg1
+                            << " occured while calling a PETSc function");
+                                             /**
+                                              * Exception
+                                              */
+            DeclException3 (ExcAccessToNonlocalRow,
+                            int, int, int,
+                            << "You tried to access row " << arg1
+                            << " of a distributed matrix, but only rows "
+                            << arg2 << " through " << arg3
+                            << " are stored locally and can be accessed.");
             
           private:
                                              /**

diff --git a/deal.II/lac/source/petsc_matrix_base.cc b/deal.II/lac/source/petsc_matrix_base.cc
index 82517195629e15c53926bd188a0a00b60a4e9de1..4819b82aa92087b334bb6b7f102ba98efbd31606 100644 (file)
--- a/deal.II/lac/source/petsc_matrix_base.cc
+++ b/deal.II/lac/source/petsc_matrix_base.cc
@@ -13,6 +13,8 @@
 
 
 #include <lac/petsc_matrix_base.h>
+#include <lac/petsc_sparse_matrix.h>
+#include <lac/petsc_parallel_sparse_matrix.h>
 #include <lac/petsc_vector.h>
 
 #ifdef DEAL_II_USE_PETSC
@@ -40,8 +42,39 @@ namespace PETScWrappers
       
                                        // otherwise first flush PETSc caches
       matrix->compress ();
-      
-      const signed int a_row_int = this->a_row;
+
+                                       // next get the number of the row
+                                       // relative to the rows this process
+                                       // actually stores (remember that
+                                       // matrices can be distributed)
+      signed int local_row;
+      if (dynamic_cast<const PETScWrappers::SparseMatrix *>(matrix) != 0)
+        {
+          local_row = this->a_row;
+        }
+      else if (dynamic_cast<const PETScWrappers::MPI::SparseMatrix *>(matrix)
+               != 0)
+        {
+                                           // first verify that the requested
+                                           // element is actually locally
+                                           // available
+          int ierr;
+          int begin, end;
+          ierr = MatGetOwnershipRange (static_cast<const Mat &>(*matrix),
+                                       &begin, &end);
+          AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+          AssertThrow ((this->a_row >= static_cast<unsigned int>(begin)) &&
+                       (this->a_row < static_cast<unsigned int>(end)),
+                       ExcAccessToNonlocalRow (this->a_row, begin, end-1));
+
+                                           // alright, so we seem to have this
+                                           // particular row on the present
+                                           // processor
+          local_row = this->a_row - begin;
+        }
+      else
+        Assert (false, ExcInternalError());
 
                                        // get a representation of the present
                                        // row
@@ -50,7 +83,7 @@ namespace PETScWrappers
       PetscScalar *values;
 
       int ierr;
-      ierr = MatGetRow(*matrix, a_row_int, &ncols, &colnums, &values);
+      ierr = MatGetRow(*matrix, local_row, &ncols, &colnums, &values);
       AssertThrow (ierr == 0, MatrixBase::ExcPETScError(ierr));
 
                                        // copy it into our caches
@@ -59,7 +92,7 @@ namespace PETScWrappers
       value_cache.reset (new std::vector<PetscScalar> (values, values+ncols));
 
                                        // and finally restore the matrix
-      ierr = MatRestoreRow(*matrix, a_row_int, &ncols, &colnums, &values);
+      ierr = MatRestoreRow(*matrix, local_row, &ncols, &colnums, &values);
       AssertThrow (ierr == 0, MatrixBase::ExcPETScError(ierr));
     }
   }