Go to the documentation of this file.
12 #include <Eigen/SparseCore>
13 #include <Eigen/Eigenvalues>
15 #include <Eigen/SparseCholesky>
17 #include "../SymGEigsSolver.h"
82 template <
typename Scalar =
long double>
115 if (has_true_BM ==
false)
137 return chol_MBM.info();
143 D_MBM_vec = D_MBM_vec.cwiseSqrt();
145 for (
int i = 0;
i < D_MBM_vec.
rows();
i++)
147 D_MBM_vec(
i) =
Complex(1.0, 0.0) / D_MBM_vec(
i);
154 Upper_MBM.template triangularView<Eigen::Upper>().solveInPlace(U_inv);
157 M =
M * right_product.real();
199 result.topLeftCorner(
A.rows(),
A.cols()) =
A;
200 result.topRightCorner(
B.rows(),
B.cols()) =
B;
202 result.bottomRightCorner(
D.rows(),
D.cols()) =
D;
211 result.block(0, 0,
A.rows(),
A.cols()) =
A;
212 result.block(0,
A.cols(),
B.rows(),
B.cols()) =
B;
214 result.block(
A.rows(), 0,
D.rows(),
D.cols()) =
D;
215 result.block(
A.rows(),
A.cols(),
E.rows(),
E.cols()) =
E;
216 result.block(
A.rows(),
A.cols() +
B.cols(),
F.rows(),
F.cols()) =
F;
217 result.block(
A.rows() +
D.rows(), 0,
G.rows(),
G.cols()) =
G;
218 result.block(
A.rows() +
D.rows(),
A.cols(),
H.rows(),
H.cols()) =
H;
219 result.block(
A.rows() +
D.rows(),
A.cols() +
B.cols(),
I.rows(),
I.cols()) =
I;
228 cmp = std::less<Scalar>{};
230 cmp = std::greater<Scalar>{};
234 epairs.insert(std::make_pair(evalues(
i), evectors.col(
i)));
237 for (
auto& epair : epairs)
239 evectors.col(
i) = epair.second;
240 evalues(
i) = epair.first;
250 std::vector<Eigen::Triplet<Scalar>> tripletList;
251 tripletList.reserve(
matrix.cols() -
int(colToRemove.size()));
253 for (
int iRow = 0; iRow <
matrix.cols(); iRow++)
255 if (std::find(colToRemove.begin(), colToRemove.end(), iRow) == colToRemove.end())
269 int BlockSize =
m_nev;
272 for (
int iCol = 0; iCol <
m_nev; iCol++)
275 for (
int iRow = 0; iRow <
m_n; iRow++)
281 if (
sqrt(sum) < tolerance_L2)
284 columnsToDelete.push_back(iCol);
301 if (
A.rows() !=
X.
rows() ||
A.rows() !=
A.cols())
302 throw std::invalid_argument(
"Wrong size");
316 if (
B.rows() !=
A.rows() ||
B.cols() !=
A.cols())
317 throw std::invalid_argument(
"Wrong size");
334 SparseMatrix directions, AX, AR, BX,
AD, ADD, DD, BDD,
BD, XAD, RAD, DAD, XBD, RBD,
BR, sparse_eVecX, sparse_eVecR, sparse_eVecD, inverse_matrix;
335 Matrix XAR, RAR, XBR, gramA, gramB, eVecX, eVecR, eVecD;
336 std::vector<int> columnsToDelete;
368 X =
X * sparse_eVecX;
369 AX = AX * sparse_eVecX;
370 BX = BX * sparse_eVecX;
373 for (
int iter_num = 0; iter_num < max_iter; iter_num++)
389 if (columnsToDelete.size() > 0)
398 columnsToDelete.clear();
430 XAR =
Matrix(
X.transpose() * AR);
436 XAD =
X.transpose() *
AD;
438 DAD = directions.transpose() *
AD;
439 XBD =
X.transpose() *
BD;
442 gramA =
stack_9_matricies(
m_evalues.asDiagonal(), XAR, XAD, XAR.transpose(), RAR, RAD, XAD.transpose(), RAD.transpose(), DAD.transpose());
443 gramB =
stack_9_matricies(Matrix::Identity(
m_nev,
m_nev), XBR, XBD, XBR.transpose(), Matrix::Identity(BlockSize, BlockSize), RBD, XBD.transpose(), RBD.transpose(), Matrix::Identity(BlockSize, BlockSize));
460 int nconv = geigs.compute();
483 sparse_eVecX = eVecX.sparseView();
484 sparse_eVecR = eVecR.sparseView();
485 sparse_eVecD = eVecD.sparseView();
488 ADD = AR * sparse_eVecR;
489 BDD =
BR * sparse_eVecR;
491 DD = DD + directions * sparse_eVecD;
492 ADD = ADD +
AD * sparse_eVecD;
493 BDD = BDD +
BD * sparse_eVecD;
500 sparse_eVecX = eVecX.sparseView();
501 sparse_eVecR = eVecR.sparseView();
504 ADD = AR * sparse_eVecR;
505 BDD =
BR * sparse_eVecR;
508 X =
X * sparse_eVecX + DD;
509 AX = AX * sparse_eVecX + ADD;
510 BX = BX * sparse_eVecX + BDD;
552 #endif // LOBPCG_SOLVER
void resize(Index rows, Index cols)
Eigen::Matrix< Complex, Eigen::Dynamic, Eigen::Dynamic > ComplexMatrix
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
bool flag_with_constraints
LOBPCGSolver(const SparseMatrix &A, const SparseMatrix X)
Array< double, 1, 3 > e(1./3., 0.5, 2.)
const MatrixU matrixU() const
void removeColumns(SparseMatrix &matrix, std::vector< int > &colToRemove)
A direct sparse LDLT Cholesky factorizations without square root.
GaussianFactorGraphValuePair Y
Eigen::Matrix< Complex, Eigen::Dynamic, 1 > ComplexVector
void setPreconditioner(const SparseMatrix &preconditioner)
std::complex< Scalar > Complex
int orthogonalizeInPlace(SparseMatrix &M, SparseMatrix &B, SparseMatrix &true_BM, bool has_true_BM=false)
const VectorType vectorD() const
Matrix stack_9_matricies(Matrix A, Matrix B, Matrix C, Matrix D, Matrix E, Matrix F, Matrix G, Matrix H, Matrix I)
bool flag_with_preconditioner
void setB(const SparseMatrix &B)
BearingRangeFactor< Pose, Point2 > BR
Eigen::Matrix< Scalar, Eigen::Dynamic, Eigen::Dynamic > Matrix
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
void compute(int maxit=10, Scalar tol_div_n=1e-7)
A small structure to hold a non zero as a triplet (i,j,value).
EigenvectorsType eigenvectors() const
Returns the eigenvectors of given matrix.
Scalar coeff(Index row, Index col) const
Map< Matrix< T, Dynamic, Dynamic, ColMajor >, 0, OuterStride<> > matrix(T *data, int rows, int cols, int stride)
void sort_epairs(Vector &evalues, Matrix &evectors, int SelectionRule)
Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > Vector
@ GEIGS_CHOLESKY
Using Cholesky decomposition to solve generalized eigenvalues.
Eigen::SparseMatrix< Complex > SparseComplexMatrix
Matrix< Scalar, Dynamic, Dynamic > C
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
JacobiRotation< float > G
void setFromTriplets(const InputIterators &begin, const InputIterators &end)
Eigen::SparseMatrix< Scalar > SparseMatrix
Computes eigenvalues and eigenvectors of general matrices.
int checkConvergence_getBlocksize(SparseMatrix &m_residuals, Scalar tolerance_L2, std::vector< int > &columnsToDelete)
Matrix stack_4_matricies(Matrix A, Matrix B, Matrix C, Matrix D)
@ SUCCESSFUL
Computation was successful.
@ SMALLEST_ALGE
Select eigenvalues with smallest algebraic value. Only for symmetric eigen solvers.
ComputationInfo info() const
SparseMatrix m_preconditioner
void applyConstraintsInPlace(SparseMatrix &X, SparseMatrix &Y, SparseMatrix &B)
Jet< T, N > sqrt(const Jet< T, N > &f)
void setConstraints(const SparseMatrix &Y)
Matrix< RealScalar, Dynamic, Dynamic > M
const EigenvalueType & eigenvalues() const
Returns the eigenvalues of given matrix.
gtsam
Author(s):
autogenerated on Sun Dec 22 2024 04:11:56