KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
eigen_sparse_lu_solver.h
Go to the documentation of this file.
1 /* KRATOS _ _ ____ _
2 // | | (_)_ __ ___ __ _ _ __/ ___| ___ | |_ _____ _ __ ___
3 // | | | | '_ \ / _ \/ _` | '__\___ \ / _ \| \ \ / / _ \ '__/ __|
4 // | |___| | | | | __/ (_| | | ___) | (_) | |\ V / __/ | \__ |
5 // |_____|_|_| |_|\___|\__,_|_| |____/ \___/|_| \_/ \___|_| |___/ Application
6 //
7 // Author: Thomas Oberbichler
8 */
9 
10 #if !defined(KRATOS_EIGEN_SPARSE_LU_SOLVER_H_INCLUDED)
11 #define KRATOS_EIGEN_SPARSE_LU_SOLVER_H_INCLUDED
12 
13 // External includes
14 #include <Eigen/Sparse>
15 
16 // Project includes
17 #include "includes/define.h"
18 #include "linear_solvers_define.h"
20 #include "spaces/ublas_space.h"
23 
24 namespace Kratos {
25 
26 template <typename TScalar = double>
28 {
29 public:
30  using Scalar = TScalar;
33 
34 private:
35  Eigen::SparseLU<SparseMatrix> m_solver;
36 
37 public:
38  static std::string Name()
39  {
40  return "eigen_sparse_lu";
41  }
42 
43  void Initialize(Parameters settings)
44  {
45  }
46 
47  bool Compute(Eigen::Map<const SparseMatrix> a)
48  {
49  m_solver.compute(a);
50 
51  const bool success = m_solver.info() == Eigen::Success;
52 
53  KRATOS_ERROR_IF(!success) << m_solver.lastErrorMessage() << std::endl;
54 
55  return success;
56  }
57 
58  bool Solve(Eigen::Ref<const Vector> b, Eigen::Ref<Vector> x) const
59  {
60  x = m_solver.solve(b);
61 
62  const bool success = m_solver.info() == Eigen::Success;
63 
64  return success;
65  }
66 
67  void PrintInfo(std::ostream &rOStream) const
68  {
69  rOStream << "EigenDirectSolver <" << Name() << "> finished.";
70  }
71 
72  std::string GetSolverErrorMessages() const
73  {
74  return m_solver.lastErrorMessage();
75  }
76 };
77 
78 } // namespace Kratos
79 
80 #endif // defined(KRATOS_EIGEN_SPARSE_LU_SOLVER_H_INCLUDED)
Definition: eigen_sparse_lu_solver.h:28
bool Compute(Eigen::Map< const SparseMatrix > a)
Definition: eigen_sparse_lu_solver.h:47
bool Solve(Eigen::Ref< const Vector > b, Eigen::Ref< Vector > x) const
Definition: eigen_sparse_lu_solver.h:58
TScalar Scalar
Definition: eigen_sparse_lu_solver.h:30
Kratos::EigenDynamicVector< Scalar > Vector
Definition: eigen_sparse_lu_solver.h:32
void Initialize(Parameters settings)
Definition: eigen_sparse_lu_solver.h:43
std::string GetSolverErrorMessages() const
Definition: eigen_sparse_lu_solver.h:72
Kratos::EigenSparseMatrix< Scalar > SparseMatrix
Definition: eigen_sparse_lu_solver.h:31
static std::string Name()
Definition: eigen_sparse_lu_solver.h:38
void PrintInfo(std::ostream &rOStream) const
Definition: eigen_sparse_lu_solver.h:67
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
#define KRATOS_ERROR_IF(conditional)
Definition: exception.h:162
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
Eigen::SparseMatrix< _Scalar, Eigen::RowMajor, int > EigenSparseMatrix
Definition: linear_solvers_define.h:35
Eigen::Matrix< _Scalar, Eigen::Dynamic, 1 > EigenDynamicVector
Definition: linear_solvers_define.h:33
a
Definition: generate_stokes_twofluid_element.py:77
b
Definition: generate_total_lagrangian_mixed_volumetric_strain_element.py:31
x
Definition: sensitivityMatrix.py:49