diff --git a/Solver/eigenSolver.cpp b/Solver/eigenSolver.cpp index d2585ebcf9b446f084e3f82db3a436f5b7370026..f594831c0b98769699de5414c2021f9fc8c4fcea 100644 --- a/Solver/eigenSolver.cpp +++ b/Solver/eigenSolver.cpp @@ -105,8 +105,10 @@ bool eigenSolver::solve(int numEigenValues, std::string which) Msg::Error("SLEPc diverged after %d iterations", its); else if(reason == EPS_DIVERGED_BREAKDOWN) Msg::Error("SLEPc generic breakdown in method"); +#if !(PETSC_VERSION_RELEASE == 0) // petsc-dev else if(reason == EPS_DIVERGED_NONSYMMETRIC) Msg::Error("The operator is nonsymmetric"); +#endif // get number of converged approximate eigenpairs PetscInt nconv; @@ -152,11 +154,20 @@ bool eigenSolver::solve(int numEigenValues, std::string which) } _eigenVectors.push_back(ev); } +#if (PETSC_VERSION_RELEASE == 0) // petsc-dev + _try(VecDestroy(&xr)); + _try(VecDestroy(&xi)); +#else _try(VecDestroy(xr)); _try(VecDestroy(xi)); +#endif } +#if (PETSC_VERSION_RELEASE == 0) // petsc-dev + _try(EPSDestroy(&eps)); +#else _try(EPSDestroy(eps)); +#endif if(reason == EPS_CONVERGED_TOL){ Msg::Debug("SLEPc done"); diff --git a/Solver/linearSystemPETSc.cpp b/Solver/linearSystemPETSc.cpp index 520e0c13cdfd6f86f760aa057927b8d8446f4da8..809d532f9fc2b4a6ae1b39bf845004600f38cdf6 100644 --- a/Solver/linearSystemPETSc.cpp +++ b/Solver/linearSystemPETSc.cpp @@ -145,9 +145,15 @@ bool linearSystemPETScBlockDouble::isAllocated() const void linearSystemPETScBlockDouble::clear() { if(_isAllocated){ +#if (PETSC_VERSION_RELEASE == 0) // petsc-dev + MatDestroy(&_a); + VecDestroy(&_x); + VecDestroy(&_b); +#else MatDestroy(_a); VecDestroy(_x); VecDestroy(_b); +#endif } _isAllocated = false; } diff --git a/Solver/linearSystemPETSc.hpp b/Solver/linearSystemPETSc.hpp index 67c35945efbd6b69341078cdc5e33b52868c2067..783b16142d6bb8437f8246cb1bd0927202fc945f 100644 --- a/Solver/linearSystemPETSc.hpp +++ b/Solver/linearSystemPETSc.hpp @@ -2,6 +2,12 @@ #include <petsc.h> #include <petscksp.h> #include "linearSystemPETSc.h" + +#if (PETSC_VERSION_RELEASE == 0) // petsc-dev +#define PetscTruth PetscBool +#define PetscOptionsGetTruth PetscOptionsGetBool +#endif + static void _try(int ierr) { CHKERRABORT(PETSC_COMM_WORLD, ierr); @@ -39,8 +45,12 @@ template <class scalar> linearSystemPETSc<scalar>::~linearSystemPETSc() { clear(); - if (_kspAllocated) - _try (KSPDestroy (_ksp)); + if(_kspAllocated) +#if (PETSC_VERSION_RELEASE == 0) // petsc-dev + _try(KSPDestroy(&_ksp)); +#else + _try(KSPDestroy(_ksp)); +#endif } template <class scalar> @@ -141,9 +151,15 @@ template <class scalar> void linearSystemPETSc<scalar>::clear() { if(_isAllocated){ +#if (PETSC_VERSION_RELEASE == 0) // petsc-dev + _try(MatDestroy(&_a)); + _try(VecDestroy(&_x)); + _try(VecDestroy(&_b)); +#else _try(MatDestroy(_a)); _try(VecDestroy(_x)); _try(VecDestroy(_b)); +#endif } _isAllocated = false; }