From: kronbichler Date: Mon, 21 Sep 2009 08:36:44 +0000 (+0000) Subject: First attempt to make the PETScWrappers also work when PETSC is configured with 64... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=677c2ec616deee9f65de7de1f58b53d33a76402f;p=dealii-svn.git First attempt to make the PETScWrappers also work when PETSC is configured with 64-bit integers. git-svn-id: https://svn.dealii.org/trunk@19485 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 17121547a0..442929fd06 100644 --- a/deal.II/lac/include/lac/petsc_matrix_base.h +++ b/deal.II/lac/include/lac/petsc_matrix_base.h @@ -1210,7 +1210,11 @@ namespace PETScWrappers * adding/inserting local data into * the (large) sparse matrix. */ - std::vector column_indices; +#ifdef PETSC_USE_64BIT_INDICES + std::vector column_indices; +#else + std::vector 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(matrix), &begin, &end); AssertThrow (ierr == 0, ExcPETScError(ierr)); diff --git a/deal.II/lac/include/lac/petsc_solver.h b/deal.II/lac/include/lac/petsc_solver.h index 45deb8ea3b..b552fd6670 100644 --- a/deal.II/lac/include/lac/petsc_solver.h +++ b/deal.II/lac/include/lac/petsc_solver.h @@ -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); diff --git a/deal.II/lac/include/lac/petsc_vector_base.h b/deal.II/lac/include/lac/petsc_vector_base.h index e4be758cda..8762b5f9a6 100644 --- a/deal.II/lac/include/lac/petsc_vector_base.h +++ b/deal.II/lac/include/lac/petsc_vector_base.h @@ -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(*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(*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(vector), &begin, &end); AssertThrow (ierr == 0, ExcPETScError(ierr)); diff --git a/deal.II/lac/source/petsc_matrix_base.cc b/deal.II/lac/source/petsc_matrix_base.cc index 2508eaf432..9bddf2731c 100644 --- a/deal.II/lac/source/petsc_matrix_base.cc +++ b/deal.II/lac/source/petsc_matrix_base.cc @@ -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 MatrixBase::local_range () const { - int begin, end; +#ifdef PETSC_USE_64BIT_INDICES + PetscInt +#else + int +#endif + begin, end; const int ierr = MatGetOwnershipRange (static_cast(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 diff --git a/deal.II/lac/source/petsc_parallel_sparse_matrix.cc b/deal.II/lac/source/petsc_parallel_sparse_matrix.cc index 0d804eb214..b3125fac3d 100644 --- a/deal.II/lac/source/petsc_parallel_sparse_matrix.cc +++ b/deal.II/lac/source/petsc_parallel_sparse_matrix.cc @@ -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 int_row_lengths (row_lengths.begin(), + row_lengths.end()); + PetscInt +#else const std::vector 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 row_lengths_in_window (local_row_end - local_row_start); - std::vector row_lengths_out_of_window (local_row_end - local_row_start); +#ifdef PETSC_USE_64BIT_INDICES + std::vector +#else + std::vector +#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 row_entries; +#ifdef PETSC_USE_64BIT_INDICES + std::vector +#else + std::vector +#endif + row_entries; std::vector row_values; for (unsigned int i=0; i rowstart_in_window (local_row_end - - local_row_start + 1, - 0); - std::vector colnums_in_window; +#ifdef PETSC_USE_64BIT_INDICES + std::vector +#else + std::vector +#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 int_row_lengths (row_lengths.begin(), - row_lengths.end()); +#ifdef PETSC_USE_64BIT_INDICES + const std::vector +#else + const std::vector +#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 row_entries; +#ifdef PETSC_USE_64BIT_INDICES + std::vector +#else + std::vector +#endif + row_entries; std::vector row_values; for (unsigned int i=0; i(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 VectorBase::local_range () const { - int begin, end; +#ifdef PETSC_USE_64BIT_INDICES + PetscInt +#else + int +#endif + begin, end; const int ierr = VecGetOwnershipRange (static_cast(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 petsc_ind (n_elements); for (unsigned int i=0; ieps, +#ifdef PETSC_USE_64BIT_INDICES + reinterpret_cast(n_converged)); +#else reinterpret_cast(n_converged)); +#endif AssertThrow (ierr == 0, ExcSLEPcError(ierr)); }