Preconditioner.h
Go to the documentation of this file.
1 /*
2  * Preconditioner.h
3  *
4  * Created on: Jun 2, 2014
5  * Author: Yong-Dian Jian
6  * Author: Sungtae An
7  */
8 
9 #pragma once
10 
11 #include <gtsam/base/Vector.h>
12 #include <memory>
13 #include <iosfwd>
14 #include <map>
15 #include <string>
16 
17 namespace gtsam {
18 
19 class GaussianFactorGraph;
20 class KeyInfo;
21 class VectorValues;
22 
23 /* parameters for the preconditioner */
24 struct GTSAM_EXPORT PreconditionerParameters {
25 
26  typedef std::shared_ptr<PreconditionerParameters> shared_ptr;
27 
28  enum Kernel { /* Preconditioner Kernel */
29  GTSAM = 0,
30  CHOLMOD /* experimental */
31  } kernel_ ;
32 
33  enum Verbosity {
34  SILENT = 0,
35  COMPLEXITY = 1,
36  ERROR = 2
37  } verbosity_ ;
38 
39  PreconditionerParameters(): kernel_(GTSAM), verbosity_(SILENT) {}
40  PreconditionerParameters(const PreconditionerParameters &p) : kernel_(p.kernel_), verbosity_(p.verbosity_) {}
42 
43  /* general interface */
44  inline Kernel kernel() const { return kernel_; }
45  inline Verbosity verbosity() const { return verbosity_; }
46 
47  void print() const;
48 
49  virtual void print(std::ostream &os) const;
50 
51  static Kernel kernelTranslator(const std::string &s);
52  static Verbosity verbosityTranslator(const std::string &s);
53  static std::string kernelTranslator(Kernel k);
54  static std::string verbosityTranslator(Verbosity v);
55 
56  /* for serialization */
57  friend std::ostream& operator<<(std::ostream &os, const PreconditionerParameters &p);
58  };
59 
60 /* PCG aims to solve the problem: A x = b by reparametrizing it as
61  * L^{-1} A L^{-T} y = L^{-1} b or M^{-1} A x = M^{-1} b,
62  * where A \approx L L^{T}, or A \approx M
63  * The goal of this class is to provide a general interface to all preconditioners */
64 class GTSAM_EXPORT Preconditioner {
65 public:
66  typedef std::shared_ptr<Preconditioner> shared_ptr;
67  typedef std::vector<size_t> Dimensions;
68 
69  /* Generic Constructor and Destructor */
71  virtual ~Preconditioner() {}
72 
73  /*
74  * Abstract interface for raw vectors. VectorValues is a speed bottleneck
75  * and Yong-Dian has profiled preconditioners (outside GTSAM) with the the
76  * three methods below. In GTSAM, unfortunately, we are still using the
77  * VectorValues methods called in iterative-inl.h
78  */
79 
81  virtual void solve(const Vector& y, Vector &x) const = 0;
82 
84  virtual void transposeSolve(const Vector& y, Vector& x) const = 0;
85 
87  virtual void build(
88  const GaussianFactorGraph &gfg,
89  const KeyInfo &info,
90  const std::map<Key,Vector> &lambda
91  ) = 0;
92 };
93 
94 /*******************************************************************************************/
97  typedef std::shared_ptr<DummyPreconditionerParameters> shared_ptr;
100 };
101 
102 /*******************************************************************************************/
103 class GTSAM_EXPORT DummyPreconditioner : public Preconditioner {
104 public:
106  typedef std::shared_ptr<DummyPreconditioner> shared_ptr;
107 
108 public:
109 
110  DummyPreconditioner() : Base() {}
111  ~DummyPreconditioner() override {}
112 
113  /* Computation Interfaces for raw vector */
114  void solve(const Vector& y, Vector &x) const override { x = y; }
115  void transposeSolve(const Vector& y, Vector& x) const override { x = y; }
116  void build(
117  const GaussianFactorGraph &gfg,
118  const KeyInfo &info,
119  const std::map<Key,Vector> &lambda
120  ) override {}
121 };
122 
123 /*******************************************************************************************/
128 };
129 
130 /*******************************************************************************************/
131 class GTSAM_EXPORT BlockJacobiPreconditioner : public Preconditioner {
132 public:
135  ~BlockJacobiPreconditioner() override ;
136 
137  /* Computation Interfaces for raw vector */
138  void solve(const Vector& y, Vector &x) const override;
139  void transposeSolve(const Vector& y, Vector& x) const override;
140  void build(
141  const GaussianFactorGraph &gfg,
142  const KeyInfo &info,
143  const std::map<Key,Vector> &lambda
144  ) override;
145 
146 protected:
147 
148  void clean() ;
149 
150  std::vector<size_t> dims_;
151  double *buffer_;
152  size_t bufferSize_;
153  size_t nnz_;
154 };
155 
156 /*********************************************************************************************/
157 /* factory method to create preconditioners */
158 std::shared_ptr<Preconditioner> createPreconditioner(const std::shared_ptr<PreconditionerParameters> parameters);
159 
160 }
161 
162 
void print(const Matrix &A, const string &s, ostream &stream)
Definition: Matrix.cpp:155
void build(const GaussianFactorGraph &gfg, const KeyInfo &info, const std::map< Key, Vector > &lambda) override
build/factorize the preconditioner
std::shared_ptr< DummyPreconditioner > shared_ptr
std::ostream & operator<<(std::ostream &os, const Dih6 &m)
Definition: testGroup.cpp:109
Scalar * y
def build
Definition: noxfile.py:86
std::shared_ptr< Preconditioner > shared_ptr
void transposeSolve(const Vector &y, Vector &x) const override
implement x = L^{-T} y
else if n * info
Eigen::VectorXd Vector
Definition: Vector.h:38
Array< int, Dynamic, 1 > v
RealScalar s
static ConjugateGradientParameters parameters
cout<< "The eigenvalues of A are:"<< endl<< ces.eigenvalues()<< endl;cout<< "The matrix of eigenvectors, V, is:"<< endl<< ces.eigenvectors()<< endl<< endl;complex< float > lambda
traits
Definition: chartTesting.h:28
typedef and functions to augment Eigen&#39;s VectorXd
ofstream os("timeSchurFactors.csv")
float * p
void solve(const Vector &y, Vector &x) const override
implement x = L^{-1} y
PreconditionerParameters Base
PreconditionerParameters(const PreconditionerParameters &p)
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
std::shared_ptr< PreconditionerParameters > shared_ptr
std::shared_ptr< Preconditioner > createPreconditioner(const std::shared_ptr< PreconditionerParameters > params)
std::shared_ptr< DummyPreconditionerParameters > shared_ptr
std::vector< size_t > Dimensions


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:35:15