From 09f401be18d0bbcb7a9c67b838f1df793142afe9 Mon Sep 17 00:00:00 2001 From: young Date: Tue, 6 Sep 2011 08:23:59 +0000 Subject: [PATCH] Upgrade PETSc matrices to petsc-dev git-svn-id: https://svn.dealii.org/trunk@24258 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/deal.II/lac/petsc_matrix_base.h | 14 ++++- deal.II/source/lac/petsc_matrix_base.cc | 59 ++++++++++++++++++- deal.II/source/lac/petsc_sparse_matrix.cc | 14 ++++- 3 files changed, 82 insertions(+), 5 deletions(-) diff --git a/deal.II/include/deal.II/lac/petsc_matrix_base.h b/deal.II/include/deal.II/lac/petsc_matrix_base.h index 4e16e9f1f4..5939f5c5ac 100644 --- a/deal.II/include/deal.II/lac/petsc_matrix_base.h +++ b/deal.II/include/deal.II/lac/petsc_matrix_base.h @@ -1108,7 +1108,12 @@ namespace PETScWrappers * $1000\times32$-bit machine * precision. */ - PetscTruth is_symmetric (const double tolerance = 1.e-12); +#if DEAL_II_PETSC_VERSION_DEV() + PetscBool +#else + PetscTruth +#endif + is_symmetric (const double tolerance = 1.e-12); #if DEAL_II_PETSC_VERSION_GTE(2,3,0) /** @@ -1120,7 +1125,12 @@ namespace PETScWrappers * $1000\times32$-bit machine * precision. */ - PetscTruth is_hermitian (const double tolerance = 1.e-12); +#if DEAL_II_PETSC_VERSION_DEV() + PetscBool +#else + PetscTruth +#endif + is_hermitian (const double tolerance = 1.e-12); #endif /* diff --git a/deal.II/source/lac/petsc_matrix_base.cc b/deal.II/source/lac/petsc_matrix_base.cc index d5ca400588..ebabbdfd7a 100644 --- a/deal.II/source/lac/petsc_matrix_base.cc +++ b/deal.II/source/lac/petsc_matrix_base.cc @@ -90,7 +90,11 @@ namespace PETScWrappers MatrixBase::~MatrixBase () { +#if DEAL_II_PETSC_VERSION_DEV() + const int ierr = MatDestroy (&matrix); +#else const int ierr = MatDestroy (matrix); +#endif AssertThrow (ierr == 0, ExcPETScError(ierr)); } @@ -100,7 +104,11 @@ namespace PETScWrappers MatrixBase::clear () { // destroy the matrix... +#if DEAL_II_PETSC_VERSION_DEV() + int ierr = MatDestroy (&matrix); +#else int ierr = MatDestroy (matrix); +#endif AssertThrow (ierr == 0, ExcPETScError(ierr)); // ...and replace it by an empty // sequential matrix @@ -146,19 +154,33 @@ namespace PETScWrappers #endif IS index_set; +#if DEAL_II_PETSC_VERSION_DEV() + ISCreateGeneral (get_mpi_communicator(), 1, &petsc_row, PETSC_COPY_VALUES, &index_set); +#else ISCreateGeneral (get_mpi_communicator(), 1, &petsc_row, &index_set); +#endif + #if DEAL_II_PETSC_VERSION_LT(2,3,0) const int ierr = MatZeroRows(matrix, index_set, &new_diag_value); +#else +#if DEAL_II_PETSC_VERSION_DEV() + const int ierr + = MatZeroRowsIS(matrix, index_set, new_diag_value, PETSC_NULL, PETSC_NULL); #else const int ierr = MatZeroRowsIS(matrix, index_set, new_diag_value); +#endif #endif AssertThrow (ierr == 0, ExcPETScError(ierr)); +#if DEAL_II_PETSC_VERSION_DEV() + ISDestroy (&index_set); +#else ISDestroy (index_set); +#endif compress (); } @@ -183,20 +205,34 @@ namespace PETScWrappers // to call them even if #rows is empty, // since this is a collective operation IS index_set; +#if DEAL_II_PETSC_VERSION_DEV() + ISCreateGeneral (get_mpi_communicator(), rows.size(), + &petsc_rows[0], PETSC_COPY_VALUES, &index_set); +#else ISCreateGeneral (get_mpi_communicator(), rows.size(), &petsc_rows[0], &index_set); +#endif #if DEAL_II_PETSC_VERSION_LT(2,3,0) const int ierr = MatZeroRows(matrix, index_set, &new_diag_value); +#else +#if DEAL_II_PETSC_VERSION_DEV() + const int ierr + = MatZeroRowsIS(matrix, index_set, new_diag_value, PETSC_NULL, PETSC_NULL); #else const int ierr = MatZeroRowsIS(matrix, index_set, new_diag_value); +#endif #endif AssertThrow (ierr == 0, ExcPETScError(ierr)); +#if DEAL_II_PETSC_VERSION_DEV() + ISDestroy (&index_set); +#else ISDestroy (index_set); +#endif compress (); } @@ -570,20 +606,39 @@ namespace PETScWrappers AssertThrow (ierr == 0, ExcPETScError(ierr)); } +#if DEAL_II_PETSC_VERSION_DEV() + PetscBool +#else PetscTruth +#endif MatrixBase::is_symmetric (const double tolerance) { - PetscTruth truth; +#if DEAL_II_PETSC_VERSION_DEV() + PetscBool +#else + PetscTruth +#endif + truth; // First flush PETSc caches compress (); MatIsSymmetric (matrix, tolerance, &truth); return truth; } +#if DEAL_II_PETSC_VERSION_DEV() + PetscBool +#else PetscTruth +#endif MatrixBase::is_hermitian (const double tolerance) { - PetscTruth truth; +#if DEAL_II_PETSC_VERSION_DEV() + PetscBool +#else + PetscTruth +#endif + truth; + // First flush PETSc caches compress (); diff --git a/deal.II/source/lac/petsc_sparse_matrix.cc b/deal.II/source/lac/petsc_sparse_matrix.cc index 83dbffbd92..0fb4b70d5d 100644 --- a/deal.II/source/lac/petsc_sparse_matrix.cc +++ b/deal.II/source/lac/petsc_sparse_matrix.cc @@ -74,7 +74,7 @@ namespace PETScWrappers } - + void SparseMatrix::reinit (const unsigned int m, const unsigned int n, @@ -83,7 +83,11 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one +#if DEAL_II_PETSC_VERSION_DEV() + const int ierr = MatDestroy (&matrix); +#else const int ierr = MatDestroy (matrix); +#endif AssertThrow (ierr == 0, ExcPETScError(ierr)); do_reinit (m, n, n_nonzero_per_row, is_symmetric); @@ -99,7 +103,11 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one +#if DEAL_II_PETSC_VERSION_DEV() + const int ierr = MatDestroy (&matrix); +#else const int ierr = MatDestroy (matrix); +#endif AssertThrow (ierr == 0, ExcPETScError(ierr)); do_reinit (m, n, row_lengths, is_symmetric); @@ -115,7 +123,11 @@ namespace PETScWrappers { // get rid of old matrix and generate a // new one +#if DEAL_II_PETSC_VERSION_DEV() + const int ierr = MatDestroy (&matrix); +#else const int ierr = MatDestroy (matrix); +#endif AssertThrow (ierr == 0, ExcPETScError(ierr)); do_reinit (sparsity_pattern, preset_nonzero_locations); -- 2.39.5