1 // ============================================================================= 2 // PROJECT CHRONO - http://projectchrono.org 3 // 4 // Copyright (c) 2014 projectchrono.org 5 // All rights reserved. 6 // 7 // Use of this source code is governed by a BSD-style license that can be found 8 // in the LICENSE file at the top level of the distribution and at 9 // http://projectchrono.org/license-chrono.txt. 10 // 11 // ============================================================================= 12 // Authors: Alessandro Tasora, Radu Serban 13 // ============================================================================= 14 15 #ifndef CH_SOLVER_H 16 #define CH_SOLVER_H 17 18 #include <vector> 19 20 #include "chrono/solver/ChConstraint.h" 21 #include "chrono/solver/ChSystemDescriptor.h" 22 #include "chrono/solver/ChVariables.h" 23 24 namespace chrono { 25 26 /// @addtogroup chrono_solver 27 /// @{ 28 29 /// Base class for all Chrono solvers (for linear problems or complementarity problems). \n 30 /// See ChSystemDescriptor for more information about the problem formulation and the data structures passed to the 31 /// solver. 32 class ChApi ChSolver { 33 34 public: 35 /// Available types of solvers. 36 enum class Type { 37 // Iterative VI solvers 38 PSOR = 0, ///< Projected SOR (Successive Over-Relaxation) 39 PSSOR, ///< Projected symmetric SOR 40 PJACOBI, ///< Projected Jacobi 41 PMINRES, ///< Projected MINRES 42 BARZILAIBORWEIN, ///< Barzilai-Borwein 43 APGD, ///< Accelerated Projected Gradient Descent 44 ADDM, ///< Alternating Direction Method of Multipliers 45 // Direct linear solvers 46 SPARSE_LU, ///< Sparse supernodal LU factorization 47 SPARSE_QR, ///< Sparse left-looking rank-revealing QR factorization 48 PARDISO_MKL, ///< Pardiso MKL (super-nodal sparse direct solver) 49 PARDISO_PROJECT, ///< Pardiso (from PardisoProject) (super-nodal sparse direct solver) 50 MUMPS, ///< Mumps (MUltifrontal Massively Parallel sparse direct Solver) 51 // Iterative linear solvers 52 GMRES, ///< Generalized Minimal RESidual Algorithm 53 MINRES, ///< MINimum RESidual method 54 BICGSTAB, ///< Bi-conjugate gradient stabilized 55 // Other 56 CUSTOM, 57 }; 58 ~ChSolver()59 virtual ~ChSolver() {} 60 61 /// Return type of the solver. GetType()62 virtual Type GetType() const { return Type::CUSTOM; } 63 64 /// Indicate whether or not the Solve() phase requires an up-to-date problem matrix. 65 /// Typically, direct solvers only need the matrix for the Setup() phase. However, iterative solvers likely require 66 /// the matrix to perform the necessary matrix-vector operations. 67 virtual bool SolveRequiresMatrix() const = 0; 68 69 /// Performs the solution of the problem. 70 /// This function MUST be implemented in children classes, with specialized 71 /// methods such as iterative or direct solvers. 72 /// Returns true if it successfully solves the problem and false otherwise. 73 virtual double Solve(ChSystemDescriptor& sysd ///< system description with constraints and variables 74 ) = 0; 75 76 /// This function does the setup operations for the solver. 77 /// The purpose of this function is to prepare the solver for subsequent calls to the 78 /// solve function. This function is called only as frequently it is determined that 79 /// it is appropriate to perform the setup phase. Setup(ChSystemDescriptor & sysd)80 virtual bool Setup(ChSystemDescriptor& sysd ///< system description with constraints and variables 81 ) { 82 return true; 83 } 84 85 /// Set verbose output from solver. SetVerbose(bool mv)86 void SetVerbose(bool mv) { verbose = mv; } 87 88 /// Method to allow serialization of transient data to archives. 89 virtual void ArchiveOUT(ChArchiveOut& marchive); 90 91 /// Method to allow de-serialization of transient data from archives. 92 virtual void ArchiveIN(ChArchiveIn& marchive); 93 94 protected: ChSolver()95 ChSolver() : verbose(false) {} 96 97 bool verbose; 98 }; 99 100 /// @} chrono_solver 101 102 } // end namespace chrono 103 104 #endif