]> https://gitweb.dealii.org/ - dealii.git/commitdiff
First attempt to make the PETScWrappers also work when PETSC is configured with 64...
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 21 Sep 2009 08:36:44 +0000 (08:36 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Mon, 21 Sep 2009 08:36:44 +0000 (08:36 +0000)
git-svn-id: https://svn.dealii.org/trunk@19485 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/petsc_matrix_base.h
deal.II/lac/include/lac/petsc_solver.h
deal.II/lac/include/lac/petsc_vector_base.h
deal.II/lac/source/petsc_matrix_base.cc
deal.II/lac/source/petsc_parallel_sparse_matrix.cc
deal.II/lac/source/petsc_solver.cc
deal.II/lac/source/petsc_sparse_matrix.cc
deal.II/lac/source/petsc_vector_base.cc
deal.II/lac/source/slepc_solver.cc

index 17121547a0aeb7a8b3402c48b3ea19c7a9c75327..442929fd06b2b1322a1bed88bc559fe055f53769 100644 (file)
@@ -1210,7 +1210,11 @@ namespace PETScWrappers
                                        * adding/inserting local data into
                                        * the (large) sparse matrix.
                                        */
-      std::vector<unsigned int> column_indices;
+#ifdef PETSC_USE_64BIT_INDICES
+      std::vector<PetscInt> column_indices;
+#else
+      std::vector<int> column_indices;
+#endif
 
                                       /**
                                        * An internal array of double values
@@ -1471,13 +1475,19 @@ namespace PETScWrappers
   {
     prepare_action(LastAction::insert);
 
-    const signed int petsc_i = row;
+#ifdef PETSC_USE_64BIT_INDICES
+    const PetscInt petsc_i = row;
+    PetscInt * col_index_ptr;
+#else
+    const int petsc_i = row;
     int * col_index_ptr;
+#endif
     PetscScalar const* col_value_ptr;
     int n_columns;
 
                                   // If we don't elide zeros, the pointers
                                   // are already available...
+#ifndef PETSC_USE_64BIT_INDICES
     if (elide_zero_values == false)
       {
        col_index_ptr = (int*)col_indices;
@@ -1485,6 +1495,7 @@ namespace PETScWrappers
        n_columns = n_cols;
       }
     else
+#endif
       {
                                   // Otherwise, extract nonzero values in
                                   // each row and get the respective index.
@@ -1510,7 +1521,7 @@ namespace PETScWrappers
          }
        Assert(n_columns <= (int)n_cols, ExcInternalError());
 
-       col_index_ptr = (int*)&column_indices[0];
+       col_index_ptr = &column_indices[0];
        col_value_ptr = &column_values[0];
       }
 
@@ -1618,13 +1629,19 @@ namespace PETScWrappers
   {
     prepare_action(LastAction::add);
 
-    const signed int petsc_i = row;
+#ifdef PETSC_USE_64BIT_INDICES
+    const PetscInt petsc_i = row;
+    PetscInt * col_index_ptr;
+#else
+    const int petsc_i = row;
     int * col_index_ptr;
+#endif
     PetscScalar const* col_value_ptr;
     int n_columns;
 
                                   // If we don't elide zeros, the pointers
                                   // are already available...
+#ifndef PETSC_USE_64BIT_INDICES
     if (elide_zero_values == false)
       {
        col_index_ptr = (int*)col_indices;
@@ -1632,6 +1649,7 @@ namespace PETScWrappers
        n_columns = n_cols;
       }
     else
+#endif
       {
                                   // Otherwise, extract nonzero values in
                                   // each row and get the respective index.
@@ -1657,7 +1675,7 @@ namespace PETScWrappers
          }
        Assert(n_columns <= (int)n_cols, ExcInternalError());
 
-       col_index_ptr = (int*)&column_indices[0];
+       col_index_ptr = &column_indices[0];
        col_value_ptr = &column_values[0];
       }
 
@@ -1734,7 +1752,11 @@ namespace PETScWrappers
   bool
   MatrixBase::in_local_range (const unsigned int index) const
   {
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt begin, end;
+#else
     int begin, end;
+#endif
     const int ierr = MatGetOwnershipRange (static_cast<const Mat &>(matrix),
                                           &begin, &end);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
index 45deb8ea3be54ca23bdde7c19bf5fa58f77e4c30..b552fd6670d74925c47ac85aec287b8f83c5d39e 100644 (file)
@@ -162,9 +162,17 @@ namespace PETScWrappers
                                         * convergence has been reached.
                                         */
       static
-      int
+#ifdef PETSC_USE_64BIT_INDICES
+      PetscErrorCode 
+#else
+      int 
+#endif
       convergence_test (KSP                 ksp,
+#ifdef PETSC_USE_64BIT_INDICES
+                        const PetscInt      iteration,
+#else
                         const int           iteration,
+#endif
                         const PetscReal     residual_norm,
                         KSPConvergedReason *reason,
                         void               *solver_control);
index e4be758cda8206ebf718d7ca6c5449944c085414..8762b5f9a677df95d6d5d09b36d9ab97e38f5c75 100644 (file)
@@ -823,7 +823,11 @@ namespace PETScWrappers
           AssertThrow (ierr == 0, ExcPETScError(ierr));
         }
 
+#ifdef PETSC_USE_64BIT_INDICES
+      const PetscInt petsc_i = index;
+#else
       const signed int petsc_i = index;
+#endif
 
       const int ierr
         = VecSetValues (vector, 1, &petsc_i, &value, INSERT_VALUES);
@@ -863,7 +867,11 @@ namespace PETScWrappers
         return *this;
 
                                        // use the PETSc function to add something
+#ifdef PETSC_USE_64BIT_INDICES
+      const PetscInt petsc_i = index;
+#else
       const signed int petsc_i = index;
+#endif
       const int ierr
         = VecSetValues (vector, 1, &petsc_i, &value, ADD_VALUES);
       AssertThrow (ierr == 0, ExcPETScError(ierr));
@@ -901,7 +909,11 @@ namespace PETScWrappers
         return *this;
 
                                        // use the PETSc function to add something
+#ifdef PETSC_USE_64BIT_INDICES
+      const PetscInt petsc_i = index;
+#else
       const signed int petsc_i = index;
+#endif
       const PetscScalar subtractand = -value;
       const int ierr
         = VecSetValues (vector, 1, &petsc_i, &subtractand, ADD_VALUES);
@@ -938,7 +950,11 @@ namespace PETScWrappers
       if (value == 1.)
         return *this;
 
+#ifdef PETSC_USE_64BIT_INDICES
+      const PetscInt petsc_i = index;
+#else
       const signed int petsc_i = index;
+#endif
 
       const PetscScalar new_value
         = static_cast<PetscScalar>(*this) * value;
@@ -978,7 +994,11 @@ namespace PETScWrappers
       if (value == 1.)
         return *this;
 
+ #ifdef PETSC_USE_64BIT_INDICES
+      const PetscInt petsc_i = index;
+#else
       const signed int petsc_i = index;
+#endif
 
       const PetscScalar new_value
         = static_cast<PetscScalar>(*this) / value;
@@ -997,7 +1017,12 @@ namespace PETScWrappers
   bool
   VectorBase::in_local_range (const unsigned int index) const
   {
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt begin, end;
+#else
     int begin, end;
+#endif
+
     const int ierr = VecGetOwnershipRange (static_cast<const Vec &>(vector),
                                           &begin, &end);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
index 2508eaf432f849e3161e646437dc246f8af1d3c3..9bddf2731ce633973d66f13b1cba837529ddfc12 100644 (file)
@@ -47,15 +47,15 @@ namespace PETScWrappers
 
                                        // get a representation of the present
                                        // row
-      int ncols;
-
 #if (PETSC_VERSION_MAJOR <= 2) && \
     ((PETSC_VERSION_MINOR < 2) ||  \
      ((PETSC_VERSION_MINOR == 2) && (PETSC_VERSION_SUBMINOR == 0)))
+      int          ncols;
       int         *colnums;
       PetscScalar *values;
 #else
-      const int         *colnums;
+      PetscInt           ncols;
+      const PetscInt    *colnums;
       const PetscScalar *values;
 #endif
 
@@ -212,8 +212,12 @@ namespace PETScWrappers
   MatrixBase::el (const unsigned int i,
                   const unsigned int j) const
   {
-    const signed int petsc_i = i;
-    const signed int petsc_j = j;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt 
+#else
+    int 
+#endif
+      petsc_i = i, petsc_j = j;
     PetscScalar value;
 
     const int ierr
@@ -259,7 +263,12 @@ namespace PETScWrappers
   unsigned int
   MatrixBase::m () const
   {
-    int n_rows, n_cols;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      n_rows, n_cols;
     int ierr = MatGetSize (matrix, &n_rows, &n_cols);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
@@ -271,7 +280,12 @@ namespace PETScWrappers
   unsigned int
   MatrixBase::n () const
   {
-    int n_rows, n_cols;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      n_rows, n_cols;
     int ierr = MatGetSize (matrix, &n_rows, &n_cols);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
@@ -283,7 +297,12 @@ namespace PETScWrappers
   unsigned int
   MatrixBase::local_size () const
   {
-    int n_rows, n_cols;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      n_rows, n_cols;
     int ierr = MatGetLocalSize (matrix, &n_rows, &n_cols);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
@@ -295,7 +314,12 @@ namespace PETScWrappers
   std::pair<unsigned int, unsigned int>
   MatrixBase::local_range () const
   {
-    int begin, end;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      begin, end;
     const int ierr = MatGetOwnershipRange (static_cast<const Mat &>(matrix),
                                           &begin, &end);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
@@ -332,15 +356,15 @@ namespace PETScWrappers
 
                                      // get a representation of the present
                                      // row
-    int ncols;
-
-#if (PETSC_VERSION_MAJOR <= 2) && \
+ #if (PETSC_VERSION_MAJOR <= 2) && \
     ((PETSC_VERSION_MINOR < 2) ||  \
      ((PETSC_VERSION_MINOR == 2) && (PETSC_VERSION_SUBMINOR == 0)))
+    int ncols;
     int         *colnums;
     PetscScalar *values;
 #else
-    const int         *colnums;
+    PetscInt ncols;
+    const PetscInt    *colnums;
     const PetscScalar *values;
 #endif
 
index 0d804eb214602d049a459578fafcbd964d384e19..b3125fac3d845aef8b0c1674a3f6a0142dd3edf8 100644 (file)
@@ -237,10 +237,19 @@ namespace PETScWrappers
                                        // signed integers. so we have to
                                        // convert, unless we want to play dirty
                                        // tricks with conversions of pointers
+#ifdef PETSC_USE_64BIT_INDICES
+      const std::vector<PetscInt> int_row_lengths (row_lengths.begin(),
+                                                  row_lengths.end());
+      PetscInt
+#else
       const std::vector<signed int> int_row_lengths (row_lengths.begin(),
                                                      row_lengths.end());
+      int
+#endif
+       petsc_m = m, petsc_n = n;
 
 //TODO: There must be a significantly better way to provide information about the off-diagonal blocks of the matrix. this way, petsc keeps allocating tiny chunks of memory, and gets completely hung up over this
+
       const int ierr
         = MatCreateMPIAIJ(communicator,
                          local_rows, local_columns,
@@ -301,8 +310,13 @@ namespace PETScWrappers
 
                                        // then count the elements in- and
                                        // out-of-window for the rows we own
-      std::vector<int> row_lengths_in_window (local_row_end - local_row_start);
-      std::vector<int> row_lengths_out_of_window (local_row_end - local_row_start);
+#ifdef PETSC_USE_64BIT_INDICES
+      std::vector<PetscInt> 
+#else
+      std::vector<int> 
+#endif
+       row_lengths_in_window (local_row_end - local_row_start),
+       row_lengths_out_of_window (local_row_end - local_row_start);
       for (unsigned int row = local_row_start; row<local_row_end; ++row)
         for (unsigned int c=0; c<sparsity_pattern.row_length(row); ++c)
           {
@@ -361,7 +375,12 @@ namespace PETScWrappers
     ((PETSC_VERSION_MINOR < 2) ||  \
      ((PETSC_VERSION_MINOR == 2) && (PETSC_VERSION_SUBMINOR == 0)))
 
-          std::vector<int> row_entries;
+#ifdef PETSC_USE_64BIT_INDICES
+          std::vector<PetscInt> 
+#else
+          std::vector<int> 
+#endif
+            row_entries;
           std::vector<PetscScalar> row_values;
           for (unsigned int i=0; i<sparsity_pattern.n_rows(); ++i)
             {
@@ -386,10 +405,13 @@ namespace PETScWrappers
                                            // dummy entry at the end to make
                                            // sure petsc doesn't read past the
                                            // end
-          std::vector<int> rowstart_in_window (local_row_end -
-                                               local_row_start + 1,
-                                               0);
-          std::vector<int> colnums_in_window;
+#ifdef PETSC_USE_64BIT_INDICES
+          std::vector<PetscInt> 
+#else
+          std::vector<int> 
+#endif
+           rowstart_in_window (local_row_end - local_row_start + 1, 0),
+           colnums_in_window;
           {
             unsigned int n_cols = 0;
             for (unsigned int i=local_row_start; i<local_row_end; ++i)
@@ -444,7 +466,12 @@ namespace PETScWrappers
 
             for (unsigned int i=local_row_start; i<local_row_end; ++i)
               {
-                const int petsc_i = i;
+#ifdef PETSC_USE_64BIT_INDICES
+               PetscInt 
+#else
+               int 
+#endif
+                 petsc_i = i;
                 MatSetValues (matrix, 1, &petsc_i,
                               sparsity_pattern.row_length(i),
                               &colnums_in_window[rowstart_in_window[i-local_row_start]],
index 04cc47144ba337d30af57cdf734207be522d52dd..2eae9aecfaad6d175caf140a54148ff6c74760d1 100644 (file)
@@ -221,7 +221,11 @@ namespace PETScWrappers
 
   int
   SolverBase::convergence_test (KSP                 /*ksp*/,
-                                const int           iteration,
+#ifdef PETSC_USE_64BIT_INDICES
+                               const PetscInt      iteration,
+#else
+                               const int           iteration,
+#endif
                                 const PetscReal     residual_norm,
                                 KSPConvergedReason *reason,
                                 void               *solver_control_x)
index 183127fd41c6113caf3198d81a3fdbb894bec22d..7bee6c7f1e2e3d5d98871f2594b8c377e1ff7fca 100644 (file)
@@ -180,8 +180,12 @@ namespace PETScWrappers
                                      // signed integers. so we have to
                                      // convert, unless we want to play dirty
                                      // tricks with conversions of pointers
-    const std::vector<signed int> int_row_lengths (row_lengths.begin(),
-                                                   row_lengths.end());
+#ifdef PETSC_USE_64BIT_INDICES
+    const std::vector<PetscInt>
+#else
+    const std::vector<int>
+#endif
+      int_row_lengths (row_lengths.begin(), row_lengths.end());
     const int ierr
       = MatCreateSeqAIJ(PETSC_COMM_SELF, m, n, 0,
                         &int_row_lengths[0], &matrix);
@@ -232,7 +236,12 @@ namespace PETScWrappers
                                      // class.
     if (preset_nonzero_locations == true)
       {
-        std::vector<int> row_entries;
+#ifdef PETSC_USE_64BIT_INDICES
+       std::vector<PetscInt>
+#else
+       std::vector<int>
+#endif
+         row_entries;
         std::vector<PetscScalar> row_values;
         for (unsigned int i=0; i<sparsity_pattern.n_rows(); ++i)
           {
@@ -241,7 +250,12 @@ namespace PETScWrappers
             for (unsigned int j=0; j<row_lengths[i]; ++j)
               row_entries[j] = sparsity_pattern.column_number (i,j);
               
-            const int int_row = i;
+#ifdef PETSC_USE_64BIT_INDICES
+           const PetscInt
+#else
+           const int
+#endif
+             int_row = i;
             MatSetValues (matrix, 1, &int_row,
                           row_lengths[i], &row_entries[0],
                           &row_values[0], INSERT_VALUES);
index ea4d9575b86e735c617a1c90afb8f212320e3064..5298f4b1940af4994b3f90dccd343e76cd834f0a 100644 (file)
@@ -63,7 +63,12 @@ namespace PETScWrappers
                                            // element is actually locally
                                            // available
           int ierr;
-          int begin, end;
+#ifdef PETSC_USE_64BIT_INDICES
+         PetscInt
+#else
+         int
+#endif
+           begin, end;
           ierr = VecGetOwnershipRange (static_cast<const Vec &>(vector),
                                        &begin, &end);
           AssertThrow (ierr == 0, ExcPETScError(ierr));
@@ -180,7 +185,12 @@ namespace PETScWrappers
   unsigned int
   VectorBase::size () const
   {
-    int sz;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      sz;
     const int ierr = VecGetSize (vector, &sz);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
@@ -192,7 +202,12 @@ namespace PETScWrappers
   unsigned int
   VectorBase::local_size () const
   {
-    int sz;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      sz;
     const int ierr = VecGetLocalSize (vector, &sz);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
 
@@ -204,7 +219,12 @@ namespace PETScWrappers
   std::pair<unsigned int, unsigned int>
   VectorBase::local_range () const
   {
-    int begin, end;
+#ifdef PETSC_USE_64BIT_INDICES
+    PetscInt
+#else
+    int
+#endif
+      begin, end;
     const int ierr = VecGetOwnershipRange (static_cast<const Vec &>(vector),
                                           &begin, &end);
     AssertThrow (ierr == 0, ExcPETScError(ierr));
@@ -990,15 +1010,13 @@ namespace PETScWrappers
                                     // (unlike the above calls)
     if (n_elements != 0)
       {
-#if (PETSC_VERSION_MAJOR <= 2) &&                                      \
-     ((PETSC_VERSION_MINOR < 2) ||                                             \
-      ((PETSC_VERSION_MINOR == 2) && (PETSC_VERSION_SUBMINOR == 0)))
-       const int * petsc_indices = indices;
-#else
+#ifdef PETSC_USE_64BIT_INDICES
        std::vector<PetscInt> petsc_ind (n_elements);
        for (unsigned int i=0; i<n_elements; ++i)
          petsc_ind[i] = indices[i];
-        const PetscInt * petsc_indices = &petsc_ind[0];
+       const PetscInt *petsc_indices = &petsc_ind[0];
+#else
+       const int * petsc_indices = (const int*)indices;
 #endif
 
        InsertMode mode = ADD_VALUES;
index 943cdc73b9dac97d7da982fa154cc6f9a2e893da..4fcc15590365cf65aada815c6153cd7cf06ffd09 100644 (file)
@@ -129,7 +129,11 @@ namespace SLEPcWrappers
                                     // get number of converged
                                     // eigenstates
     ierr = EPSGetConverged (solver_data->eps, 
+#ifdef PETSC_USE_64BIT_INDICES
+                           reinterpret_cast<PetscInt *>(n_converged));
+#else
                            reinterpret_cast<int *>(n_converged));
+#endif
     AssertThrow (ierr == 0, ExcSLEPcError(ierr));
   }
 

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.