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: Dario Mangoni, Radu Serban
13 // =============================================================================
14 
15 #ifndef CHSOLVERMKL_H
16 #define CHSOLVERMKL_H
17 
18 #include "chrono_pardisomkl/ChApiPardisoMKL.h"
19 #include "chrono/solver/ChDirectSolverLS.h"
20 
21 #include <Eigen/PardisoSupport>
22 
23 namespace chrono {
24 
25 /// @addtogroup pardisomkl_module
26 /// @{
27 
28 /** \class ChSolverPardisoMKL
29 \brief Interface to the Intel MKL Pardiso parallel sparse direct solver.
30 
31 Sparse linear direct solver.
32 Cannot handle VI and complementarity problems, so it cannot be used with NSC formulations.
33 
34 The solver is equipped with two main features:
35 - sparsity pattern lock
36 - sparsity pattern learning
37 
38 See ChDirectSolverLS for more details.
39 
40 <div class="ce-warning">
41 If appropriate and warranted by the problem setup, it is \e highly recommended to enable the sparsity pattern \e lock.
42 This can significantly improve performance for more complex problems (larger size and/or problems which include
43 constraints).
44 </div>
45 
46 Minimal usage example, to be put anywhere in the code, before starting the main simulation loop:
47 \code{.cpp}
48 auto mkl_solver = chrono_types::make_shared<ChSolverPardisoMKL>();
49 system.SetSolver(mkl_solver);
50 \endcode
51 
52 See ChSystemDescriptor for more information about the problem formulation and the data structures passed to the solver.
53 */
54 class ChApiPardisoMKL ChSolverPardisoMKL : public ChDirectSolverLS {
55   public:
56     /// Construct an MKL Pardiso sparse direct solver object and specify the number of OpenMP threads.
57     /// Passing the default value num_threads=0 results in using a number of threads equal to the number
58     /// of available processors (as returned by the function omp_get_num_procs)
59     ChSolverPardisoMKL(int num_threads = 0);
60 
~ChSolverPardisoMKL()61     ~ChSolverPardisoMKL() {}
62 
GetType()63     virtual Type GetType() const override { return Type::PARDISO_MKL; }
64 
65     /// Get a handle to the underlying MKL engine.
GetMklEngine()66     Eigen::PardisoLU<ChSparseMatrix>& GetMklEngine() { return m_engine; }
67 
68   private:
69     /// Factorize the current sparse matrix and return true if successful.
70     virtual bool FactorizeMatrix() override;
71 
72     /// Solve the linear system using the current factorization and right-hand side vector.
73     /// Load the solution vector (already of appropriate size) and return true if succesful.
74     virtual bool SolveSystem() override;
75 
76     /// Display an error message corresponding to the last failure.
77     /// This function is only called if Factorize or Solve returned false.
78     virtual void PrintErrorMessage() override;
79 
80     Eigen::PardisoLU<ChSparseMatrix> m_engine;  ///< underlying Eigen Pardiso interface
81 };
82 
83 /// @} mkl_module
84 
85 }  // end namespace chrono
86 
87 #endif
88