chrono::ChSolverMKL Class Reference

Description

Interface to the Intel MKL Pardiso parallel sparse direct solver.

Sparse linear direct solver. Cannot handle VI and complementarity problems, so it cannot be used with NSC formulations.

The solver is equipped with two main features:

  • sparsity pattern lock
  • sparsity pattern learning

See ChDirectSolverLS for more details.

If appropriate and warranted by the problem setup, it is highly recommended to enable the sparsity pattern lock. This can significantly improve performance for more complex problems (larger size and/or problems which include constraints).

Minimal usage example, to be put anywhere in the code, before starting the main simulation loop:

auto mkl_solver = chrono_types::make_shared<ChSolverMKL>();
system.SetSolver(mkl_solver);

See ChSystemDescriptor for more information about the problem formulation and the data structures passed to the solver.

#include <ChSolverMKL.h>

Inheritance diagram for chrono::ChSolverMKL:
Collaboration diagram for chrono::ChSolverMKL:

Public Member Functions

 ChSolverMKL (int num_threads=0)
 Construct an MKL Pardiso sparse direct solver object and specify the number of OpenMP threads. More...
 
virtual Type GetType () const override
 Return type of the solver.
 
Eigen::PardisoLU< ChSparseMatrix > & GetMklEngine ()
 Get a handle to the underlying MKL engine.
 
- Public Member Functions inherited from chrono::ChDirectSolverLS
void LockSparsityPattern (bool val)
 Enable/disable locking the sparsity pattern (default: false). More...
 
void UseSparsityPatternLearner (bool val)
 Enable/disable use of the sparsity pattern learner (default: enabled). More...
 
void ForceSparsityPatternUpdate ()
 Force a call to the sparsity pattern learner to update sparsity pattern on the underlying matrix. More...
 
void SetSparsityEstimate (double sparsity)
 Set estimate for matrix sparsity, a value in [0,1], with 0 indicating a fully dense matrix (default: 0.9). More...
 
virtual void SetMatrixSymmetryType (MatrixSymmetryType symmetry)
 Set the matrix symmetry type (default: GENERAL).
 
void UsePermutationVector (bool val)
 Enable/disable use of permutation vector (default: false). More...
 
void LeverageRhsSparsity (bool val)
 Enable/disable leveraging sparsity in right-hand side vector (default: false). More...
 
virtual void EnableNullPivotDetection (bool val, double threshold=0)
 Enable detection of null pivots. More...
 
void ResetTimers ()
 Reset timers for internal phases in Solve and Setup.
 
double GetTimeSolve_Assembly () const
 Get cumulative time for assembly operations in Solve phase.
 
double GetTimeSolve_SolverCall () const
 Get cumulative time for Pardiso calls in Solve phase.
 
double GetTimeSetup_Assembly () const
 Get cumulative time for assembly operations in Setup phase.
 
double GetTimeSetup_SolverCall () const
 Get cumulative time for Pardiso calls in Setup phase.
 
int GetNumSetupCalls () const
 Return the number of calls to the solver's Setup function.
 
int GetNumSolveCalls () const
 Return the number of calls to the solver's Setup function.
 
ChSparseMatrixGetMatrix ()
 Get a handle to the underlying matrix.
 
virtual bool Setup (ChSystemDescriptor &sysd) override
 Perform the solver setup operations. More...
 
virtual double Solve (ChSystemDescriptor &sysd) override
 Solve linear system. More...
 
virtual void ArchiveOUT (ChArchiveOut &marchive) override
 Method to allow serialization of transient data to archives.
 
virtual void ArchiveIN (ChArchiveIn &marchive) override
 Method to allow de serialization of transient data from archives.
 
- Public Member Functions inherited from chrono::ChSolver
void SetVerbose (bool mv)
 Set verbose output from solver.
 

Additional Inherited Members

- Public Types inherited from chrono::ChDirectSolverLS
enum  MatrixSymmetryType { MatrixSymmetryType::GENERAL, MatrixSymmetryType::SYMMETRIC_POSDEF, MatrixSymmetryType::SYMMETRIC_INDEF, MatrixSymmetryType::STRUCTURAL_SYMMETRIC }
 
- Public Types inherited from chrono::ChSolver
enum  Type {
  Type::PSOR = 0, Type::PSSOR, Type::PJACOBI, Type::PMINRES,
  Type::BARZILAIBORWEIN, Type::APGD, Type::SPARSE_LU, Type::SPARSE_QR,
  Type::PARDISO, Type::MUMPS, Type::GMRES, Type::MINRES,
  Type::BICGSTAB, CUSTOM
}
 Available types of solvers. More...
 
- Protected Member Functions inherited from chrono::ChDirectSolverLS
virtual bool SolveRequiresMatrix () const override
 Indicate whether or not the Solve() phase requires an up-to-date problem matrix. More...
 
- Protected Attributes inherited from chrono::ChDirectSolverLS
ChSparseMatrix m_mat
 problem matrix
 
int m_dim
 problem size
 
MatrixSymmetryType m_symmetry
 symmetry of problem matrix
 
double m_sparsity
 user-supplied estimate of matrix sparsity
 
ChVectorDynamic< double > m_rhs
 right-hand side vector
 
ChVectorDynamic< double > m_sol
 solution vector
 
int m_solve_call
 counter for calls to Solve
 
int m_setup_call
 counter for calls to Setup
 
bool m_lock
 is the matrix sparsity pattern locked?
 
bool m_use_learner
 use the sparsity pattern learner?
 
bool m_force_update
 force a call to the sparsity pattern learner?
 
bool m_use_perm
 use of the permutation vector?
 
bool m_use_rhs_sparsity
 leverage right-hand side sparsity?
 
bool m_null_pivot_detection
 enable detection of zero pivots?
 
ChTimer m_timer_setup_assembly
 timer for matrix assembly
 
ChTimer m_timer_setup_solvercall
 timer for factorization
 
ChTimer m_timer_solve_assembly
 timer for RHS assembly
 
ChTimer m_timer_solve_solvercall
 timer for solution
 
- Protected Attributes inherited from chrono::ChSolver
bool verbose
 

Constructor & Destructor Documentation

◆ ChSolverMKL()

chrono::ChSolverMKL::ChSolverMKL ( int  num_threads = 0)

Construct an MKL Pardiso sparse direct solver object and specify the number of OpenMP threads.

Passing the default value num_threads=0 results in using a number of threads equal to the number of available processors (as returned by the function omp_get_num_procs)


The documentation for this class was generated from the following files:
  • /builds/uwsbel/chrono/src/chrono_mkl/ChSolverMKL.h
  • /builds/uwsbel/chrono/src/chrono_mkl/ChSolverMKL.cpp