/****************************************************************************** * * AMDiS - Adaptive multidimensional simulations * * Copyright (C) 2013 Dresden University of Technology. All Rights Reserved. * Web: https://fusionforge.zih.tu-dresden.de/projects/amdis * * Authors: * Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al. * * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. * * * This file is part of AMDiS * * See also license.opensource.txt in the distribution. * ******************************************************************************/ /** \file PetscSolver.h */ #ifndef AMDIS_PETSC_SOLVER_SEQ_H #define AMDIS_PETSC_SOLVER_SEQ_H #ifdef HAVE_SEQ_PETSC #include "solver/LinearSolver.h" #include "solver/MatrixStreams.h" #include "solver/PetscTypes.h" #include "Timer.h" #include #include #include #include #include #include #include #include #include namespace AMDiS { class PetscRunner : public OEMRunner { public: /// Constructor of standard PETSc runner. Reads ksp and pc parameters from initfile. PetscRunner(LinearSolver* oem_); typedef SolverMatrix > BlockMatrix; /// initialize the solver \ref ksp and preconditioner \ref pc void init(const BlockMatrix& A, const Mat& fullMatrix); /// solve the linear equation \f$A\cdot x = b \f$ by applying the PETSc solver \ref ksp int solve(const Mat& A, Vec& x, const Vec& b); /// destroy solver \ref ksp and preconditioner \ref pc virtual void exit() { KSPDestroy(&ksp); } /// get the PETSc solver \ref ksp KSP getSolver() { return ksp; } /// get the PETSc preconditioner \ref pc PC getPc() { return pc; } ~PetscRunner() { } protected: static void createSubSolver(KSP &ksp_, Mat m, std::string kspPrefix_); static void setSolver(KSP ksp_, std::string kspPrefix_, KSPType kspType, PCType pcType, PetscReal rtol = PETSC_DEFAULT, PetscReal atol = PETSC_DEFAULT, PetscInt maxIt = PETSC_DEFAULT); LinearSolver& oem; private: /// PETSc solver object KSP ksp; /// PETSc preconditioner object PC pc; /// KSP database prefix std::string kspPrefix; bool zeroStartVector; }; /** * \ingroup Solver * * \brief * Wrapper for the external PETSc solver: * http://www.mcs.anl.gov/petsc/ * * This is a suite of data structures and routines for the * scalable (parallel) solution of scientific applications. */ template< typename Runner = PetscRunner > class PetscSolver : public LinearSolver { public: /// Creator class used in the LinearSolverMap. class Creator : public LinearSolverCreator { public: virtual ~Creator() {} /// Returns a new PetscSolver object. LinearSolver* create() { return new PetscSolver(this->name); } }; PetscSolver(std::string n) : LinearSolver(n), runner(this) {} ~PetscSolver() {} virtual OEMRunner* getRunner() { return &runner; } void setNestedVectors(bool createNestedVector = true) { vecSol.createNestedVector = createNestedVector; vecRhs.createNestedVector = createNestedVector; } protected: Runner runner; int solveLinearSystem(const SolverMatrix >& A, SystemVector& x, SystemVector& b, bool createMatrixData, bool storeMatrixData) { FUNCNAME("PetscSolver::solveLinearSystem()"); Timer t; // transfer matrix and rhs-vector to PETSc data-structures if (createMatrixData) { petscMat << A; if (x.getSize() == 1) runner.init(A, petscMat.nestMat[0]); else runner.init(A, petscMat.matrix); } vecSol << x; vecRhs << b; INFO(info, 8)("fill PETSc matrix needed %.5f seconds\n", t.elapsed()); // solve the linear system using PETSc solvers t.reset(); if (vecSol.nestVec.size() == 1) error = runner.solve(petscMat.nestMat[0], vecSol.nestVec[0], vecRhs.nestVec[0]); else error = runner.solve(petscMat.matrix, vecSol.vector, vecRhs.vector); // transfer solution from PETSc vector to SystemVector vecSol.vector >> x; vecSol.destroy(); vecRhs.destroy(); if (!storeMatrixData) petscMat.destroy(); runner.exit(); return error; } private: /// PETSc System-Matrix PetscMatrix petscMat; /// Solution and RHS vectors. PetscVector vecSol, vecRhs; }; } // end namespace AMDiS #endif // HAVE_SEQ_PETSC #endif // AMDIS_PETSC_SOLVER_SEQ_H