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_dense_col_piv_householder_qr_solver.h
Go to the documentation of this file.
1 /* KRATOS _ _ ____ _
2 // | | (_)_ __ ___ __ _ _ __/ ___| ___ | |_ _____ _ __ ___
3 // | | | | '_ \ / _ \/ _` | '__\___ \ / _ \| \ \ / / _ \ '__/ __|
4 // | |___| | | | | __/ (_| | | ___) | (_) | |\ V / __/ | \__ |
5 // |_____|_|_| |_|\___|\__,_|_| |____/ \___/|_| \_/ \___|_| |___/ Application
6 //
7 // Author: Quirin Aumann
8 */
9 
10 #if !defined(KRATOS_EIGEN_DENSE_COL_PIV_HOUSEHOLDER_QR_SOLVER_H_INCLUDED)
11 #define KRATOS_EIGEN_DENSE_COL_PIV_HOUSEHOLDER_QR_SOLVER_H_INCLUDED
12 
13 // External includes
14 #include <Eigen/QR>
15 
16 // Project includes
17 #include "includes/define.h"
18 #include "linear_solvers_define.h"
19 
20 namespace Kratos {
21 
22 template <typename TScalar = double>
24 {
25 public:
26  using Scalar = TScalar;
29 
30 private:
31  Eigen::ColPivHouseholderQR<Matrix> m_solver;
32 
33 public:
34  static std::string Name()
35  {
36  return "complex_dense_col_piv_householder_qr";
37  }
38 
39  void Initialize(Parameters settings)
40  {
41  }
42 
43  bool Compute(Eigen::Map<Matrix> a)
44  {
45  m_solver.compute(a);
46 
47  const bool success = m_solver.info() == Eigen::Success;
48 
49  return success;
50  }
51 
52  bool Solve(Eigen::Ref<const Vector> b, Eigen::Ref<Vector> x) const
53  {
54  x = m_solver.solve(b);
55 
56  const bool success = m_solver.info() == Eigen::Success;
57 
58  return success;
59  }
60 
61  bool SolveMultiple(Eigen::Ref<const Matrix> b, Eigen::Ref<Matrix> x) const
62  {
63  x = m_solver.solve(b);
64 
65  const bool success = m_solver.info() == Eigen::Success;
66 
67  return success;
68  }
69 
70  void PrintInfo(std::ostream &rOStream) const
71  {
72  rOStream << "EigenDirectSolver <" << Name() << "> finished.";
73  }
74 
75  std::string GetSolverErrorMessages() const
76  {
77  return "No additional information";
78  }
79 };
80 
81 } // namespace Kratos
82 
83 #endif // defined(KRATOS_EIGEN_DENSE_COL_PIV_HOUSEHOLDER_QR_SOLVER_H_INCLUDED)
Definition: eigen_dense_col_piv_householder_qr_solver.h:24
bool SolveMultiple(Eigen::Ref< const Matrix > b, Eigen::Ref< Matrix > x) const
Definition: eigen_dense_col_piv_householder_qr_solver.h:61
void Initialize(Parameters settings)
Definition: eigen_dense_col_piv_householder_qr_solver.h:39
std::string GetSolverErrorMessages() const
Definition: eigen_dense_col_piv_householder_qr_solver.h:75
bool Compute(Eigen::Map< Matrix > a)
Definition: eigen_dense_col_piv_householder_qr_solver.h:43
bool Solve(Eigen::Ref< const Vector > b, Eigen::Ref< Vector > x) const
Definition: eigen_dense_col_piv_householder_qr_solver.h:52
Kratos::EigenDynamicMatrix< Scalar > Matrix
Definition: eigen_dense_col_piv_householder_qr_solver.h:27
void PrintInfo(std::ostream &rOStream) const
Definition: eigen_dense_col_piv_householder_qr_solver.h:70
Kratos::EigenDynamicVector< Scalar > Vector
Definition: eigen_dense_col_piv_householder_qr_solver.h:28
static std::string Name()
Definition: eigen_dense_col_piv_householder_qr_solver.h:34
TScalar Scalar
Definition: eigen_dense_col_piv_householder_qr_solver.h:26
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
Eigen::Matrix< _Scalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > EigenDynamicMatrix
Definition: linear_solvers_define.h:32
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