base/Matrix.h
Go to the documentation of this file.
1 /* ----------------------------------------------------------------------------
2 
3  * GTSAM Copyright 2010, 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 
23 // \callgraph
24 
25 #pragma once
26 
28 #include <gtsam/base/Vector.h>
29 
30 #include <vector>
31 
37 namespace gtsam {
38 
39 typedef Eigen::MatrixXd Matrix;
41 
42 // Create handy typedefs and constants for square-size matrices
43 // MatrixMN, MatrixN = MatrixNN, I_NxN, and Z_NxN, for M,N=1..9
44 #define GTSAM_MAKE_MATRIX_DEFS(N) \
45 using Matrix##N = Eigen::Matrix<double, N, N>; \
46 using Matrix1##N = Eigen::Matrix<double, 1, N>; \
47 using Matrix2##N = Eigen::Matrix<double, 2, N>; \
48 using Matrix3##N = Eigen::Matrix<double, 3, N>; \
49 using Matrix4##N = Eigen::Matrix<double, 4, N>; \
50 using Matrix5##N = Eigen::Matrix<double, 5, N>; \
51 using Matrix6##N = Eigen::Matrix<double, 6, N>; \
52 using Matrix7##N = Eigen::Matrix<double, 7, N>; \
53 using Matrix8##N = Eigen::Matrix<double, 8, N>; \
54 using Matrix9##N = Eigen::Matrix<double, 9, N>; \
55 static const Eigen::MatrixBase<Matrix##N>::IdentityReturnType I_##N##x##N = Matrix##N::Identity(); \
56 static const Eigen::MatrixBase<Matrix##N>::ConstantReturnType Z_##N##x##N = Matrix##N::Zero();
57 
67 
68 // Matrix expressions for accessing parts of matrices
69 typedef Eigen::Block<Matrix> SubMatrix;
70 typedef Eigen::Block<const Matrix> ConstSubMatrix;
71 
72 // Matrix formatting arguments when printing.
73 // Akin to Matlab style.
74 const Eigen::IOFormat& matlabFormat();
75 
79 template <class MATRIX>
80 bool equal_with_abs_tol(const Eigen::DenseBase<MATRIX>& A, const Eigen::DenseBase<MATRIX>& B, double tol = 1e-9) {
81 
82  const size_t n1 = A.cols(), m1 = A.rows();
83  const size_t n2 = B.cols(), m2 = B.rows();
84 
85  if(m1!=m2 || n1!=n2) return false;
86 
87  for(size_t i=0; i<m1; i++)
88  for(size_t j=0; j<n1; j++) {
89  if(!fpEqual(A(i,j), B(i,j), tol, false)) {
90  return false;
91  }
92  }
93  return true;
94 }
95 
99 inline bool operator==(const Matrix& A, const Matrix& B) {
100  return equal_with_abs_tol(A,B,1e-9);
101 }
102 
106 inline bool operator!=(const Matrix& A, const Matrix& B) {
107  return !(A==B);
108  }
109 
113 GTSAM_EXPORT bool assert_equal(const Matrix& A, const Matrix& B, double tol = 1e-9);
114 
118 GTSAM_EXPORT bool assert_inequal(const Matrix& A, const Matrix& B, double tol = 1e-9);
119 
123 GTSAM_EXPORT bool assert_equal(const std::list<Matrix>& As, const std::list<Matrix>& Bs, double tol = 1e-9);
124 
128 GTSAM_EXPORT bool linear_independent(const Matrix& A, const Matrix& B, double tol = 1e-9);
129 
133 GTSAM_EXPORT bool linear_dependent(const Matrix& A, const Matrix& B, double tol = 1e-9);
134 
136 template<class MATRIX>
137 inline MATRIX prod(const MATRIX& A, const MATRIX&B) {
138  MATRIX result = A * B;
139  return result;
140 }
141 
145 GTSAM_EXPORT void print(const Matrix& A, const std::string& s, std::ostream& stream);
146 
150 GTSAM_EXPORT void print(const Matrix& A, const std::string& s = "");
151 
155 GTSAM_EXPORT void save(const Matrix& A, const std::string &s, const std::string& filename);
156 
162 GTSAM_EXPORT std::istream& operator>>(std::istream& inputStream, Matrix& destinationMatrix);
163 
173 template<class MATRIX>
174 Eigen::Block<const MATRIX> sub(const MATRIX& A, size_t i1, size_t i2, size_t j1, size_t j2) {
175  size_t m=i2-i1, n=j2-j1;
176  return A.block(i1,j1,m,n);
177 }
178 
187 template <typename Derived1, typename Derived2>
188 void insertSub(Eigen::MatrixBase<Derived1>& fullMatrix, const Eigen::MatrixBase<Derived2>& subMatrix, size_t i, size_t j) {
189  fullMatrix.block(i, j, subMatrix.rows(), subMatrix.cols()) = subMatrix;
190 }
191 
195 GTSAM_EXPORT Matrix diag(const std::vector<Matrix>& Hs);
196 
203 template<class MATRIX>
204 const typename MATRIX::ConstColXpr column(const MATRIX& A, size_t j) {
205  return A.col(j);
206 }
207 
214 template<class MATRIX>
215 const typename MATRIX::ConstRowXpr row(const MATRIX& A, size_t j) {
216  return A.row(j);
217 }
218 
224 template<class MATRIX>
225 void zeroBelowDiagonal(MATRIX& A, size_t cols=0) {
226  const size_t m = A.rows(), n = A.cols();
227  const size_t k = (cols) ? std::min(cols, std::min(m,n)) : std::min(m,n);
228  for (size_t j=0; j<k; ++j)
229  A.col(j).segment(j+1, m-(j+1)).setZero();
230 }
231 
235 inline Matrix trans(const Matrix& A) { return A.transpose(); }
236 
238 template <int OutM, int OutN, int OutOptions, int InM, int InN, int InOptions>
239 struct Reshape {
240  //TODO replace this with Eigen's reshape function as soon as available. (There is a PR already pending : https://bitbucket.org/eigen/eigen/pull-request/41/reshape/diff)
243  return in.data();
244  }
245 };
246 
248 template <int M, int InOptions>
249 struct Reshape<M, M, InOptions, M, M, InOptions> {
252  return in;
253  }
254 };
255 
257 template <int M, int N, int InOptions>
258 struct Reshape<M, N, InOptions, M, N, InOptions> {
261  return in;
262  }
263 };
264 
266 template <int M, int N, int InOptions>
267 struct Reshape<N, M, InOptions, M, N, InOptions> {
270  return in.transpose();
271  }
272 };
273 
274 template <int OutM, int OutN, int OutOptions, int InM, int InN, int InOptions>
276  static_assert(InM * InN == OutM * OutN);
278 }
279 
286 GTSAM_EXPORT std::pair<Matrix,Matrix> qr(const Matrix& A);
287 
293 GTSAM_EXPORT void inplace_QR(Matrix& A);
294 
303 GTSAM_EXPORT std::list<std::tuple<Vector, double, double> >
305 
313 GTSAM_EXPORT void householder_(Matrix& A, size_t k, bool copy_vectors=true);
314 
321 GTSAM_EXPORT void householder(Matrix& A, size_t k);
322 
330 GTSAM_EXPORT Vector backSubstituteUpper(const Matrix& U, const Vector& b, bool unit=false);
331 
339 //TODO: is this function necessary? it isn't used
340 GTSAM_EXPORT Vector backSubstituteUpper(const Vector& b, const Matrix& U, bool unit=false);
341 
349 GTSAM_EXPORT Vector backSubstituteLower(const Matrix& L, const Vector& b, bool unit=false);
350 
357 GTSAM_EXPORT Matrix stack(size_t nrMatrices, ...);
358 GTSAM_EXPORT Matrix stack(const std::vector<Matrix>& blocks);
359 
370 GTSAM_EXPORT Matrix collect(const std::vector<const Matrix *>& matrices, size_t m = 0, size_t n = 0);
371 GTSAM_EXPORT Matrix collect(size_t nrMatrices, ...);
372 
379 GTSAM_EXPORT void vector_scale_inplace(const Vector& v, Matrix& A, bool inf_mask = false); // row
380 GTSAM_EXPORT Matrix vector_scale(const Vector& v, const Matrix& A, bool inf_mask = false); // row
381 GTSAM_EXPORT Matrix vector_scale(const Matrix& A, const Vector& v, bool inf_mask = false); // column
382 
394 inline Matrix3 skewSymmetric(double wx, double wy, double wz) {
395  return (Matrix3() << 0.0, -wz, +wy, +wz, 0.0, -wx, -wy, +wx, 0.0).finished();
396 }
397 
398 template <class Derived>
400  return skewSymmetric(w(0), w(1), w(2));
401 }
402 
404 GTSAM_EXPORT Matrix inverse_square_root(const Matrix& A);
405 
407 GTSAM_EXPORT Matrix cholesky_inverse(const Matrix &A);
408 
421 GTSAM_EXPORT void svd(const Matrix& A, Matrix& U, Vector& S, Matrix& V);
422 
430 GTSAM_EXPORT std::tuple<int, double, Vector>
431 DLT(const Matrix& A, double rank_tol = 1e-9);
432 
438 GTSAM_EXPORT Matrix expm(const Matrix& A, size_t K=7);
439 
440 std::string formatMatrixIndented(const std::string& label, const Matrix& matrix, bool makeVectorHorizontal = false);
441 
448 template <int N>
452 
456  OptionalJacobian<N, N> H2 = {}) const {
457  const MatrixN invA = A.inverse();
458  const VectorN c = invA * b;
459  // The derivative in A is just -[c[0]*invA c[1]*invA ... c[N-1]*invA]
460  if (H1)
461  for (size_t j = 0; j < N; j++)
462  H1->template middleCols<N>(N * j) = -c[j] * invA;
463  // The derivative in b is easy, as invA*b is just a linear map:
464  if (H2) *H2 = invA;
465  return c;
466  }
467 };
468 
474 template <typename T, int N>
476  inline constexpr static auto M = traits<T>::dimension;
479 
480  // The function phi should calculate f(a)*b, with derivatives in a and b.
481  // Naturally, the derivative in b is f(a).
482  typedef std::function<VectorN(
485 
488 
490  VectorN operator()(const T& a, const VectorN& b,
491  OptionalJacobian<N, M> H1 = {},
492  OptionalJacobian<N, N> H2 = {}) const {
493  MatrixN A;
494  phi_(a, b, {}, A); // get A = f(a) by calling f once
495  const MatrixN invA = A.inverse();
496  const VectorN c = invA * b;
497 
498  if (H1) {
500  phi_(a, c, H, {}); // get derivative H of forward mapping
501  *H1 = -invA* H;
502  }
503  if (H2) *H2 = invA;
504  return c;
505  }
506 
507  private:
508  const Operator phi_;
509 };
510 
511 GTSAM_EXPORT Matrix LLt(const Matrix& A);
512 
513 GTSAM_EXPORT Matrix RtR(const Matrix& A);
514 
515 GTSAM_EXPORT Vector columnNormSquare(const Matrix &A);
516 } // namespace gtsam
w
RowVector3d w
Definition: Matrix_resize_int.cpp:3
gtsam::weighted_eliminate
list< std::tuple< Vector, double, double > > weighted_eliminate(Matrix &A, Vector &b, const Vector &sigmas)
Definition: Matrix.cpp:261
H
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 y set format x g set format y g set format x2 g set format y2 g set format z g set angles radians set nogrid set key title set key left top Right noreverse box linetype linewidth samplen spacing width set nolabel set noarrow set nologscale set logscale x set set pointsize set encoding default set nopolar set noparametric set set set set surface set nocontour set clabel set mapping cartesian set nohidden3d set cntrparam order set cntrparam linear set cntrparam levels auto set cntrparam points set size set set xzeroaxis lt lw set x2zeroaxis lt lw set yzeroaxis lt lw set y2zeroaxis lt lw set tics in set ticslevel set tics set mxtics default set mytics default set mx2tics default set my2tics default set xtics border mirror norotate autofreq set ytics border mirror norotate autofreq set ztics border nomirror norotate autofreq set nox2tics set noy2tics set timestamp bottom norotate set rrange[*:*] noreverse nowriteback set trange[*:*] noreverse nowriteback set urange[*:*] noreverse nowriteback set vrange[*:*] noreverse nowriteback set xlabel matrix size set x2label set timefmt d m y n H
Definition: gnuplot_common_settings.hh:74
gtsam::stream
static void stream(std::ostream &os, const DiscreteValues &x, const KeyFormatter &keyFormatter)
Definition: DiscreteValues.cpp:30
gtsam::MatrixRowMajor
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > MatrixRowMajor
Definition: base/Matrix.h:40
gtsam::DLT
std::tuple< int, double, Vector > DLT(const Matrix &A, double rank_tol)
Definition: Matrix.cpp:556
Eigen
Namespace containing all symbols from the Eigen library.
Definition: jet.h:637
gtsam::operator>>
istream & operator>>(istream &inputStream, Matrix &destinationMatrix)
Definition: Matrix.cpp:163
gtsam::column
const MATRIX::ConstColXpr column(const MATRIX &A, size_t j)
Definition: base/Matrix.h:204
gtsam::MultiplyWithInverseFunction
Definition: base/Matrix.h:475
gtsam::zeroBelowDiagonal
void zeroBelowDiagonal(MATRIX &A, size_t cols=0)
Definition: base/Matrix.h:225
Eigen::Block
Expression of a fixed-size or dynamic-size block.
Definition: Block.h:103
Vector.h
typedef and functions to augment Eigen's VectorXd
benchmark.n1
n1
Definition: benchmark.py:79
gtsam::MultiplyWithInverse::operator()
VectorN operator()(const MatrixN &A, const VectorN &b, OptionalJacobian< N, N *N > H1={}, OptionalJacobian< N, N > H2={}) const
A.inverse() * b, with optional derivatives.
Definition: base/Matrix.h:454
s
RealScalar s
Definition: level1_cplx_impl.h:126
e
Array< double, 1, 3 > e(1./3., 0.5, 2.)
test_constructor::sigmas
Vector1 sigmas
Definition: testHybridNonlinearFactor.cpp:52
gtsam::householder
void householder(Matrix &A, size_t k)
Definition: Matrix.cpp:342
c
Scalar Scalar * c
Definition: benchVecAdd.cpp:17
gtsam::inverse_square_root
Matrix inverse_square_root(const Matrix &A)
Definition: Matrix.cpp:540
gtsam::diag
Matrix diag(const std::vector< Matrix > &Hs)
Definition: Matrix.cpp:196
gtsam::expm
T expm(const Vector &x, int K=7)
Definition: Lie.h:317
gtsam::MultiplyWithInverseFunction::Operator
std::function< VectorN(const T &, const VectorN &, OptionalJacobian< N, M >, OptionalJacobian< N, N >)> Operator
Definition: base/Matrix.h:484
gtsam::Reshape::reshape
static ReshapedType reshape(const Eigen::Matrix< double, InM, InN, InOptions > &in)
Definition: base/Matrix.h:242
B
Definition: test_numpy_dtypes.cpp:299
m1
Matrix3d m1
Definition: IOFormat.cpp:2
gtsam::householder_
void householder_(Matrix &A, size_t k, bool copy_vectors)
Definition: Matrix.cpp:315
gtsam::equal_with_abs_tol
bool equal_with_abs_tol(const Eigen::DenseBase< MATRIX > &A, const Eigen::DenseBase< MATRIX > &B, double tol=1e-9)
Definition: base/Matrix.h:80
gtsam::MultiplyWithInverseFunction::VectorN
Eigen::Matrix< double, N, 1 > VectorN
Definition: base/Matrix.h:477
gtsam::Matrix
Eigen::MatrixXd Matrix
Definition: base/Matrix.h:39
gtsam::matlabFormat
const Eigen::IOFormat & matlabFormat()
Definition: Matrix.cpp:129
gtsam::qr
pair< Matrix, Matrix > qr(const Matrix &A)
Definition: Matrix.cpp:224
gtsam::operator==
bool operator==(const Matrix &A, const Matrix &B)
Definition: base/Matrix.h:99
gtsam::reshape
Reshape< OutM, OutN, OutOptions, InM, InN, InOptions >::ReshapedType reshape(const Eigen::Matrix< double, InM, InN, InOptions > &m)
Definition: base/Matrix.h:275
gtsam::MultiplyWithInverseFunction::operator()
VectorN operator()(const T &a, const VectorN &b, OptionalJacobian< N, M > H1={}, OptionalJacobian< N, N > H2={}) const
f(a).inverse() * b, with optional derivatives
Definition: base/Matrix.h:490
gtsam::Vector
Eigen::VectorXd Vector
Definition: Vector.h:39
gtsam::skewSymmetric
Matrix3 skewSymmetric(double wx, double wy, double wz)
Definition: base/Matrix.h:394
result
Values result
Definition: OdometryOptimize.cpp:8
gtsam::MultiplyWithInverseFunction::phi_
const Operator phi_
Definition: base/Matrix.h:508
n
int n
Definition: BiCGSTAB_simple.cpp:1
m2
MatrixType m2(n_dims)
gtsam::print
void print(const Matrix &A, const string &s, ostream &stream)
Definition: Matrix.cpp:145
benchmark.n2
n2
Definition: benchmark.py:85
A
Definition: test_numpy_dtypes.cpp:298
gtsam::sub
Eigen::Block< const MATRIX > sub(const MATRIX &A, size_t i1, size_t i2, size_t j1, size_t j2)
Definition: base/Matrix.h:174
j
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2
gtsam::vector_scale_inplace
void vector_scale_inplace(const Vector &v, Matrix &A, bool inf_mask)
Definition: Matrix.cpp:470
gtsam::row
const MATRIX::ConstRowXpr row(const MATRIX &A, size_t j)
Definition: base/Matrix.h:215
relicense.filename
filename
Definition: relicense.py:57
GTSAM_MAKE_MATRIX_DEFS
#define GTSAM_MAKE_MATRIX_DEFS(N)
Definition: base/Matrix.h:44
gtsam::stack
Matrix stack(size_t nrMatrices,...)
Definition: Matrix.cpp:385
gtsam::MultiplyWithInverseFunction::MatrixN
Eigen::Matrix< double, N, N > MatrixN
Definition: base/Matrix.h:478
gtsam::RtR
Matrix RtR(const Matrix &A)
Definition: Matrix.cpp:518
gtsam::cholesky_inverse
Matrix cholesky_inverse(const Matrix &A)
Definition: Matrix.cpp:527
gtsam::Reshape
Reshape functor.
Definition: base/Matrix.h:239
gtsam::fpEqual
bool fpEqual(double a, double b, double tol, bool check_relative_also)
Definition: Vector.cpp:42
L
MatrixXd L
Definition: LLT_example.cpp:6
OptionalJacobian.h
Special class for optional Jacobian arguments.
m
Matrix3f m
Definition: AngleAxis_mimic_euler.cpp:1
gtsam::columnNormSquare
Vector columnNormSquare(const Matrix &A)
Definition: Matrix.cpp:213
Eigen::Triplet< double >
gtsam::Reshape< M, M, InOptions, M, M, InOptions >::ReshapedType
const typedef Eigen::Matrix< double, M, M, InOptions > & ReshapedType
Definition: base/Matrix.h:250
Eigen::Map
A matrix or vector expression mapping an existing array of data.
Definition: Map.h:94
gtsam::Reshape::ReshapedType
Eigen::Map< const Eigen::Matrix< double, OutM, OutN, OutOptions > > ReshapedType
Definition: base/Matrix.h:241
matrix
Map< Matrix< T, Dynamic, Dynamic, ColMajor >, 0, OuterStride<> > matrix(T *data, int rows, int cols, int stride)
Definition: gtsam/3rdparty/Eigen/blas/common.h:110
gtsam::svd
void svd(const Matrix &A, Matrix &U, Vector &S, Matrix &V)
Definition: Matrix.cpp:548
ConstColXpr
const typedef Block< const Derived, internal::traits< Derived >::RowsAtCompileTime, 1, !IsRowMajor > ConstColXpr
Definition: BlockMethods.h:15
gtsam::MultiplyWithInverseFunction::M
constexpr static auto M
Definition: base/Matrix.h:476
gtsam::b
const G & b
Definition: Group.h:79
gtsam::save
void save(const Matrix &A, const string &s, const string &filename)
Definition: Matrix.cpp:156
gtsam::MultiplyWithInverse
Definition: base/Matrix.h:449
a
ArrayXXi a
Definition: Array_initializer_list_23_cxx11.cpp:1
gtsam::inplace_QR
void inplace_QR(Matrix &A)
Definition: Matrix.cpp:626
gtsam::Reshape< N, M, InOptions, M, N, InOptions >::reshape
static ReshapedType reshape(const Eigen::Matrix< double, M, N, InOptions > &in)
Definition: base/Matrix.h:269
gtsam::backSubstituteUpper
Vector backSubstituteUpper(const Matrix &U, const Vector &b, bool unit)
Definition: Matrix.cpp:365
gtsam::Reshape< N, M, InOptions, M, N, InOptions >::ReshapedType
Eigen::Matrix< double, M, N, InOptions >::ConstTransposeReturnType ReshapedType
Definition: base/Matrix.h:268
gtsam::operator!=
bool operator!=(const Matrix &A, const Matrix &B)
Definition: base/Matrix.h:106
gtsam
traits
Definition: SFMdata.h:40
gtsam::traits
Definition: Group.h:36
i1
double i1(double x)
Definition: i1.c:150
K
#define K
Definition: igam.h:8
gtsam::OptionalJacobian
Definition: OptionalJacobian.h:38
gtsam::vector_scale
Matrix vector_scale(const Vector &v, const Matrix &A, bool inf_mask)
Definition: Matrix.cpp:486
gtsam::trans
Matrix trans(const Matrix &A)
Definition: base/Matrix.h:235
gtsam::MultiplyWithInverse::MatrixN
Eigen::Matrix< double, N, N > MatrixN
Definition: base/Matrix.h:451
gtsam::backSubstituteLower
Vector backSubstituteLower(const Matrix &L, const Vector &b, bool unit)
Definition: Matrix.cpp:355
gtsam::insertSub
void insertSub(Eigen::MatrixBase< Derived1 > &fullMatrix, const Eigen::MatrixBase< Derived2 > &subMatrix, size_t i, size_t j)
Definition: base/Matrix.h:188
v
Array< int, Dynamic, 1 > v
Definition: Array_initializer_list_vector_cxx11.cpp:1
gtsam::Reshape< M, M, InOptions, M, M, InOptions >::reshape
static ReshapedType reshape(const Eigen::Matrix< double, M, M, InOptions > &in)
Definition: base/Matrix.h:251
gtsam::assert_equal
bool assert_equal(const Matrix &expected, const Matrix &actual, double tol)
Definition: Matrix.cpp:41
Eigen::PlainObjectBase< Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > >::data
EIGEN_DEVICE_FUNC const EIGEN_STRONG_INLINE Scalar * data() const
Definition: PlainObjectBase.h:247
min
#define min(a, b)
Definition: datatypes.h:19
gtsam::tol
const G double tol
Definition: Group.h:79
U
@ U
Definition: testDecisionTree.cpp:342
V
MatrixXcd V
Definition: EigenSolver_EigenSolver_MatrixType.cpp:15
gtsam::Reshape< M, N, InOptions, M, N, InOptions >::ReshapedType
const typedef Eigen::Matrix< double, M, N, InOptions > & ReshapedType
Definition: base/Matrix.h:259
gtsam::assert_inequal
bool assert_inequal(const Matrix &A, const Matrix &B, double tol)
Definition: Matrix.cpp:61
Eigen::Matrix
The matrix class, also used for vectors and row-vectors.
Definition: 3rdparty/Eigen/Eigen/src/Core/Matrix.h:178
setZero
v setZero(3)
N
#define N
Definition: igam.h:9
Eigen::MatrixBase
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
gtsam::linear_independent
bool linear_independent(const Matrix &A, const Matrix &B, double tol)
Definition: Matrix.cpp:101
gtsam::collect
Matrix collect(const std::vector< const Matrix * > &matrices, size_t m, size_t n)
Definition: Matrix.cpp:431
cols
int cols
Definition: Tutorial_commainit_02.cpp:1
gtsam::Reshape< M, N, InOptions, M, N, InOptions >::reshape
static ReshapedType reshape(const Eigen::Matrix< double, M, N, InOptions > &in)
Definition: base/Matrix.h:260
ConstRowXpr
const typedef Block< const Derived, 1, internal::traits< Derived >::ColsAtCompileTime, IsRowMajor > ConstRowXpr
Definition: BlockMethods.h:18
gtsam::prod
MATRIX prod(const MATRIX &A, const MATRIX &B)
Definition: base/Matrix.h:137
gtsam::MultiplyWithInverse::VectorN
Eigen::Matrix< double, N, 1 > VectorN
Definition: base/Matrix.h:450
j1
double j1(double x)
Definition: j1.c:174
gtsam::MultiplyWithInverseFunction::MultiplyWithInverseFunction
MultiplyWithInverseFunction(const Operator &phi)
Construct with function as explained above.
Definition: base/Matrix.h:487
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
gtsam::linear_dependent
bool linear_dependent(const Matrix &A, const Matrix &B, double tol)
Definition: Matrix.cpp:115
S
DiscreteKey S(1, 2)
gtsam::formatMatrixIndented
std::string formatMatrixIndented(const std::string &label, const Matrix &matrix, bool makeVectorHorizontal)
Definition: Matrix.cpp:587
gtsam::LLt
Matrix LLt(const Matrix &A)
Definition: Matrix.cpp:511
M
Matrix< RealScalar, Dynamic, Dynamic > M
Definition: bench_gemm.cpp:51


gtsam
Author(s):
autogenerated on Sun Feb 16 2025 04:02:04