testAcceleratedPowerMethod.cpp
Go to the documentation of this file.
1 /* ----------------------------------------------------------------------------
2 
3  * GTSAM Copyright 2010-2019, Georgia Tech Research Corporation,
4  * Atlanta, Georgia 30332-0415
5  * All Rights Reserved
6  * Authors: Frank Dellaert, et al. (see THANKS for the full author list)
7 
8  * See LICENSE for the license information
9 
10  * -------------------------------------------------------------------------- */
11 
22 #include <gtsam/base/Matrix.h>
23 #include <gtsam/base/VectorSpace.h>
24 #include <gtsam/inference/Symbol.h>
28 
29 #include <Eigen/Core>
30 #include <Eigen/Dense>
31 #include <Eigen/Eigenvalues>
32 #include <iostream>
33 #include <random>
34 
35 using namespace std;
36 using namespace gtsam;
37 
38 /* ************************************************************************* */
39 TEST(AcceleratedPowerMethod, acceleratedPowerIteration) {
40  // test power iteration, beta is set to 0
41  Sparse A(6, 6);
42  A.coeffRef(0, 0) = 6;
43  A.coeffRef(1, 1) = 5;
44  A.coeffRef(2, 2) = 4;
45  A.coeffRef(3, 3) = 3;
46  A.coeffRef(4, 4) = 2;
47  A.coeffRef(5, 5) = 1;
48  Vector initial = (Vector(6) << 0.24434602, 0.22829942, 0.70094486, 0.15463092, 0.55871359,
49  0.2465342).finished();
50  const double ev1 = 6.0;
51 
52  // test accelerated power iteration
54  apf.compute(100, 1e-5);
55  EXPECT_LONGS_EQUAL(6, apf.eigenvector().rows());
56 
57  Vector6 actual1 = apf.eigenvector();
58  const double ritzValue = actual1.dot(A * actual1);
59  const double ritzResidual = (A * actual1 - ritzValue * actual1).norm();
60  EXPECT_DOUBLES_EQUAL(0, ritzResidual, 1e-5);
61 
62  EXPECT_DOUBLES_EQUAL(ev1, apf.eigenvalue(), 1e-5);
63 }
64 
65 /* ************************************************************************* */
66 TEST(AcceleratedPowerMethod, useFactorGraphSparse) {
67  // Let's make a scalar synchronization graph with 4 nodes
69 
70  // Get eigenvalues and eigenvectors with Eigen
71  auto L = fg.hessian();
73 
74  // find the index of the max eigenvalue
75  size_t maxIdx = 0;
76  for (auto i = 0; i < solver.eigenvalues().rows(); ++i) {
77  if (solver.eigenvalues()(i).real() >= solver.eigenvalues()(maxIdx).real())
78  maxIdx = i;
79  }
80  // Store the max eigenvalue and its according eigenvector
81  const auto ev1 = solver.eigenvalues()(maxIdx).real();
82 
83  Vector disturb = Vector4::Random();
84  disturb.normalize();
85  Vector initial = L.first.row(0);
86  double magnitude = initial.norm();
87  initial += 0.03 * magnitude * disturb;
89  apf.compute(100, 1e-5);
90  // Check if the eigenvalue is the maximum eigen value
91  EXPECT_DOUBLES_EQUAL(ev1, apf.eigenvalue(), 1e-8);
92 
93  // Check if the according ritz residual converged to the threshold
94  Vector actual1 = apf.eigenvector();
95  const double ritzValue = actual1.dot(L.first * actual1);
96  const double ritzResidual = (L.first * actual1 - ritzValue * actual1).norm();
97  EXPECT_DOUBLES_EQUAL(0, ritzResidual, 1e-5);
98 }
99 
100 /* ************************************************************************* */
101 TEST(AcceleratedPowerMethod, useFactorGraphDense) {
102  // Let's make a scalar synchronization graph with 10 nodes
104 
105  // Get eigenvalues and eigenvectors with Eigen
106  auto L = fg.hessian();
108 
109  // find the index of the max eigenvalue
110  size_t maxIdx = 0;
111  for (auto i = 0; i < solver.eigenvalues().rows(); ++i) {
112  if (solver.eigenvalues()(i).real() >= solver.eigenvalues()(maxIdx).real())
113  maxIdx = i;
114  }
115  // Store the max eigenvalue and its according eigenvector
116  const auto ev1 = solver.eigenvalues()(maxIdx).real();
117 
118  Vector disturb = Vector10::Random();
119  disturb.normalize();
120  Vector initial = L.first.row(0);
121  double magnitude = initial.norm();
122  initial += 0.03 * magnitude * disturb;
124  apf.compute(100, 1e-5);
125  // Check if the eigenvalue is the maximum eigen value
126  EXPECT_DOUBLES_EQUAL(ev1, apf.eigenvalue(), 1e-8);
127 
128  // Check if the according ritz residual converged to the threshold
129  Vector actual1 = apf.eigenvector();
130  const double ritzValue = actual1.dot(L.first * actual1);
131  const double ritzResidual = (L.first * actual1 - ritzValue * actual1).norm();
132  EXPECT_DOUBLES_EQUAL(0, ritzResidual, 1e-5);
133 }
134 
135 /* ************************************************************************* */
136 int main() {
137  TestResult tr;
138  return TestRegistry::runAllTests(tr);
139 }
140 /* ************************************************************************* */
TestRegistry::runAllTests
static int runAllTests(TestResult &result)
Definition: TestRegistry.cpp:27
gtsam::linear::test::example::createDenseGraph
GaussianFactorGraph createDenseGraph()
Definition: powerMethodExample.h:49
GaussianFactorGraph.h
Linear Factor Graph where all factors are Gaussians.
Eigen::SparseMatrix< double >
e
Array< double, 1, 3 > e(1./3., 0.5, 2.)
EXPECT_LONGS_EQUAL
#define EXPECT_LONGS_EQUAL(expected, actual)
Definition: Test.h:154
TestHarness.h
Matrix.h
typedef and functions to augment Eigen's MatrixXd
gtsam::GaussianFactorGraph::hessian
std::pair< Matrix, Vector > hessian(const Ordering &ordering) const
Definition: GaussianFactorGraph.cpp:264
TEST
TEST(AcceleratedPowerMethod, acceleratedPowerIteration)
Definition: testAcceleratedPowerMethod.cpp:39
real
float real
Definition: datatypes.h:10
gtsam::Vector
Eigen::VectorXd Vector
Definition: Vector.h:38
gtsam::linear::test::example::createSparseGraph
GaussianFactorGraph createSparseGraph()
Definition: powerMethodExample.h:35
rows
int rows
Definition: Tutorial_commainit_02.cpp:1
solver
BiCGSTAB< SparseMatrix< double > > solver
Definition: BiCGSTAB_simple.cpp:5
A
Matrix< SCALARA, Dynamic, Dynamic, opt_A > A
Definition: bench_gemm.cpp:48
gtsam::GaussianFactorGraph
Definition: GaussianFactorGraph.h:73
AcceleratedPowerMethod.h
accelerated power method for fast eigenvalue and eigenvector computation
gtsam::AcceleratedPowerMethod::compute
bool compute(size_t maxIterations, double tol)
Definition: AcceleratedPowerMethod.h:156
gtsam::AcceleratedPowerMethod
Compute maximum Eigenpair with accelerated power method.
Definition: AcceleratedPowerMethod.h:51
gtsam::PowerMethod::eigenvalue
double eigenvalue() const
Return the eigenvalue.
Definition: PowerMethod.h:147
L
MatrixXd L
Definition: LLT_example.cpp:6
Symbol.h
EXPECT_DOUBLES_EQUAL
#define EXPECT_DOUBLES_EQUAL(expected, actual, threshold)
Definition: Test.h:161
TestResult
Definition: TestResult.h:26
gtsam::PowerMethod::eigenvector
Vector eigenvector() const
Return the eigenvector.
Definition: PowerMethod.h:150
gtsam
traits
Definition: SFMdata.h:40
std
Definition: BFloat16.h:88
main
int main()
Definition: testAcceleratedPowerMethod.cpp:136
Eigen::EigenSolver
Computes eigenvalues and eigenvectors of general matrices.
Definition: EigenSolver.h:64
powerMethodExample.h
Create sparse and dense factor graph for PowerMethod/AcceleratedPowerMethod.
initial
Definition: testScenarioRunner.cpp:148
ceres::Vector
Eigen::Matrix< double, Eigen::Dynamic, 1 > Vector
Definition: gtsam/3rdparty/ceres/eigen.h:38
VectorSpace.h
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9


gtsam
Author(s):
autogenerated on Fri Nov 1 2024 03:39:33