GeneralBrokenLines V03-01-03
using EIGEN
|
(Symmetric) Bordered Band Matrix. More...
#include <BorderedBandMatrix.h>
Public Member Functions | |
BorderedBandMatrix () | |
Create bordered band matrix. More... | |
virtual | ~BorderedBandMatrix () |
void | resize (unsigned int nSize, unsigned int nBorder=1, unsigned int nBand=7) |
Resize bordered band matrix. More... | |
void | setZero () |
Set content to 0. More... | |
void | solveAndInvertBorderedBand (const VVector &aRightHandSide, VVector &aSolution) |
Solve linear equation system, partially calculate inverse. More... | |
void | addBlockMatrix (double aWeight, const std::vector< unsigned int > *anIndex, const std::vector< double > *aVector) |
Add symmetric block matrix. More... | |
void | addBlockMatrix (double aWeight, unsigned int nSimple, unsigned int *anIndex, double *aVector) |
Add symmetric block matrix. More... | |
Eigen::MatrixXd | getBlockMatrix (const std::vector< unsigned int > anIndex) const |
Retrieve symmetric block matrix. More... | |
Eigen::MatrixXd | getBlockMatrix (unsigned int aSize, unsigned int *anIndex) const |
Retrieve symmetric block matrix. More... | |
void | printMatrix () const |
Print bordered band matrix. More... | |
double | getBandCondition () const |
Get condition from band (decomposition) More... | |
Private Member Functions | |
void | decomposeBand () |
(root free) Cholesky decomposition of band part: C=LDL^T More... | |
VVector | solveBand (const VVector &aRightHandSide) const |
Solve for band part. More... | |
VMatrix | solveBand (const VMatrix &aRightHandSide) const |
solve band part for mixed part (border rows). More... | |
VMatrix | invertBand () |
Invert band part. More... | |
VMatrix | bandOfAVAT (const VMatrix &anArray, const VSymMatrix &aSymArray) const |
Calculate band part of: 'anArray * aSymArray * anArray.T'. More... | |
Private Attributes | |
unsigned int | numSize |
Matrix size. More... | |
unsigned int | numBorder |
Border size. More... | |
unsigned int | numBand |
Band width. More... | |
unsigned int | numCol |
Band matrix size. More... | |
VSymMatrix | theBorder |
Border part. More... | |
VMatrix | theMixed |
Mixed part. More... | |
VMatrix | theBand |
Band part. More... | |
(Symmetric) Bordered Band Matrix.
Separate storage of border, mixed and band parts (as vector<double>).
* Example for matrix size=8 with border size and band width of two * * +- -+ * | B11 B12 M13 M14 M15 M16 M17 M18 | * | B12 B22 M23 M24 M25 M26 M27 M28 | * | M13 M23 C33 C34 C35 0. 0. 0. | * | M14 M24 C34 C44 C45 C46 0. 0. | * | M15 M25 C35 C45 C55 C56 C57 0. | * | M16 M26 0. C46 C56 C66 C67 C68 | * | M17 M27 0. 0. C57 C67 C77 C78 | * | M18 M28 0. 0. 0. C68 C78 C88 | * +- -+ * * Is stored as:: * * +- -+ +- -+ * | B11 B12 | | M13 M14 M15 M16 M17 M18 | * | B12 B22 | | M23 M24 M25 M26 M27 M28 | * +- -+ +- -+ * * +- -+ * | C33 C44 C55 C66 C77 C88 | * | C34 C45 C56 C67 C78 0. | * | C35 C46 C57 C68 0. 0. | * +- -+ *
Definition at line 76 of file BorderedBandMatrix.h.
gbl::BorderedBandMatrix::BorderedBandMatrix | ( | ) |
Create bordered band matrix.
Definition at line 37 of file BorderedBandMatrix.cpp.
|
virtual |
Definition at line 41 of file BorderedBandMatrix.cpp.
void gbl::BorderedBandMatrix::addBlockMatrix | ( | double | aWeight, |
const std::vector< unsigned int > * | anIndex, | ||
const std::vector< double > * | aVector | ||
) |
Add symmetric block matrix.
Add (extended) block matrix defined by 'aVector * aWeight * aVector.T' to bordered band matrix: BBmatrix(anIndex(i),anIndex(j)) += aVector(i) * aWeight * aVector(j).
aWeight | [in] Weight |
anIndex | [in] List of rows/colums to be used |
aVector | [in] Vector |
Definition at line 77 of file BorderedBandMatrix.cpp.
References numBand, numBorder, theBand, theBorder, and theMixed.
Referenced by gbl::GblTrajectory::buildLinearEquationSystem().
void gbl::BorderedBandMatrix::addBlockMatrix | ( | double | aWeight, |
unsigned int | aSize, | ||
unsigned int * | anIndex, | ||
double * | aVector | ||
) |
Add symmetric block matrix.
Add (extended) block matrix defined by 'aVector * aWeight * aVector.T' to bordered band matrix: BBmatrix(anIndex(i),anIndex(j)) += aVector(i) * aWeight * aVector(j).
aWeight | [in] Weight |
aSize | [in] Size of block matrix |
anIndex | [in] List of rows/colums to be used |
aVector | [in] Vector |
Definition at line 111 of file BorderedBandMatrix.cpp.
References numBand, numBorder, theBand, theBorder, and theMixed.
|
private |
Calculate band part of: 'anArray * aSymArray * anArray.T'.
Definition at line 398 of file BorderedBandMatrix.cpp.
References numBand, numBorder, and numCol.
Referenced by solveAndInvertBorderedBand().
|
private |
(root free) Cholesky decomposition of band part: C=LDL^T
Decompose band matrix into diagonal matrix D and lower triangular band matrix L (diagonal=1). Overwrite band matrix with D and off-diagonal part of L.
2 | : matrix is singular. |
3 | : matrix is not positive definite. |
Definition at line 281 of file BorderedBandMatrix.cpp.
References numBand, numCol, and theBand.
Referenced by solveAndInvertBorderedBand().
double gbl::BorderedBandMatrix::getBandCondition | ( | ) | const |
Get condition from band (decomposition)
Definition at line 260 of file BorderedBandMatrix.cpp.
References numCol, and theBand.
Referenced by gbl::GblTrajectory::getBandCondition().
MatrixXd gbl::BorderedBandMatrix::getBlockMatrix | ( | const std::vector< unsigned int > | anIndex | ) | const |
Retrieve symmetric block matrix.
Get (compressed) block from bordered band matrix: aMatrix(i,j) = BBmatrix(anIndex(i),anIndex(j)).
anIndex | [in] List of rows/colums to be used |
Definition at line 138 of file BorderedBandMatrix.cpp.
References numBorder, theBand, theBorder, and theMixed.
Referenced by gbl::GblTrajectory::getExtResults(), gbl::GblTrajectory::getResAndErr(), and gbl::GblTrajectory::getResults().
MatrixXd gbl::BorderedBandMatrix::getBlockMatrix | ( | unsigned int | aSize, |
unsigned int * | anIndex | ||
) | const |
Retrieve symmetric block matrix.
Get (compressed) block from bordered band matrix: aMatrix(i,j) = BBmatrix(anIndex(i),anIndex(j)).
aSize | [in] Matrix size |
anIndex | [in] Array of rows/colums to be used |
Definition at line 167 of file BorderedBandMatrix.cpp.
|
private |
Invert band part.
Definition at line 374 of file BorderedBandMatrix.cpp.
References numBand, numCol, and theBand.
Referenced by solveAndInvertBorderedBand().
void gbl::BorderedBandMatrix::printMatrix | ( | ) | const |
Print bordered band matrix.
Definition at line 250 of file BorderedBandMatrix.cpp.
References gbl::VMatrix::print(), gbl::VSymMatrix::print(), theBand, theBorder, and theMixed.
Referenced by gbl::GblTrajectory::printTrajectory().
void gbl::BorderedBandMatrix::resize | ( | unsigned int | nSize, |
unsigned int | nBorder = 1 , |
||
unsigned int | nBand = 7 |
||
) |
Resize bordered band matrix.
nSize | [in] Size of matrix |
nBorder | [in] Size of border (=1 for q/p + additional local parameters) |
nBand | [in] Band width (usually = 5, for simplified jacobians = 4, +2 for steps) |
Definition at line 50 of file BorderedBandMatrix.cpp.
References numBand, numBorder, numCol, numSize, gbl::VSymMatrix::resize(), gbl::VMatrix::resize(), theBand, theBorder, and theMixed.
Referenced by gbl::GblTrajectory::buildLinearEquationSystem().
void gbl::BorderedBandMatrix::setZero | ( | ) |
Set content to 0.
Definition at line 62 of file BorderedBandMatrix.cpp.
References gbl::VMatrix::setZero(), gbl::VSymMatrix::setZero(), theBand, theBorder, and theMixed.
Referenced by gbl::GblTrajectory::buildLinearEquationSystem().
void gbl::BorderedBandMatrix::solveAndInvertBorderedBand | ( | const VVector & | aRightHandSide, |
VVector & | aSolution | ||
) |
Solve linear equation system, partially calculate inverse.
Solve linear equation A*x=b system with bordered band matrix A, calculate bordered band part of inverse of A. Use decomposition in border and band part for block matrix algebra:
| A Ct | | x1 | | b1 | , A is the border part | | * | | = | | , Ct is the mixed part | C D | | x2 | | b2 | , D is the band part
Explicit inversion of D is avoided by using solution X of D*X=C (X=D^-1*C, obtained from Cholesky decomposition and forward/backward substitution)
| x1 | | E*b1 - E*Xt*b2 | , E^-1 = A-Ct*D^-1*C = A-Ct*X | | = | | | x2 | | x - X*x1 | , x is solution of D*x=b2 (x=D^-1*b2)
Inverse matrix is:
| E -E*Xt | | | , only band part of (D^-1 + X*E*Xt) | -X*E D^-1 + X*E*Xt | is calculated
[in] | aRightHandSide | Right hand side (vector) 'b' of A*x=b |
[out] | aSolution | Solution (vector) x of A*x=b |
Definition at line 217 of file BorderedBandMatrix.cpp.
References bandOfAVAT(), decomposeBand(), gbl::VVector::getVec(), gbl::VSymMatrix::invert(), invertBand(), numBorder, numCol, gbl::VVector::putVec(), solveBand(), theBand, theBorder, theMixed, and gbl::VMatrix::transpose().
Referenced by gbl::GblTrajectory::fit().
solve band part for mixed part (border rows).
Solve C*X=B for mixed part using decomposition C=LDL^T and forward and backward substitution.
[in] | aRightHandSide | Right hand side (matrix) 'B' of C*X=B |
Definition at line 345 of file BorderedBandMatrix.cpp.
References gbl::VMatrix::getNumCols(), gbl::VMatrix::getNumRows(), numBorder, and theBand.
Solve for band part.
Solve C*x=b for band part using decomposition C=LDL^T and forward (L*z=b) and backward substitution (L^T*x=D^-1*z).
[in] | aRightHandSide | Right hand side (vector) 'b' of C*x=b |
Definition at line 316 of file BorderedBandMatrix.cpp.
References gbl::VMatrix::getNumCols(), gbl::VMatrix::getNumRows(), and theBand.
Referenced by solveAndInvertBorderedBand().
|
private |
Band width.
Definition at line 100 of file BorderedBandMatrix.h.
Referenced by addBlockMatrix(), bandOfAVAT(), decomposeBand(), invertBand(), and resize().
|
private |
Border size.
Definition at line 99 of file BorderedBandMatrix.h.
Referenced by addBlockMatrix(), bandOfAVAT(), getBlockMatrix(), resize(), solveAndInvertBorderedBand(), and solveBand().
|
private |
Band matrix size.
Definition at line 101 of file BorderedBandMatrix.h.
Referenced by bandOfAVAT(), decomposeBand(), getBandCondition(), invertBand(), resize(), and solveAndInvertBorderedBand().
|
private |
|
private |
Band part.
Definition at line 104 of file BorderedBandMatrix.h.
Referenced by addBlockMatrix(), decomposeBand(), getBandCondition(), getBlockMatrix(), invertBand(), printMatrix(), resize(), setZero(), solveAndInvertBorderedBand(), and solveBand().
|
private |
Border part.
Definition at line 102 of file BorderedBandMatrix.h.
Referenced by addBlockMatrix(), getBlockMatrix(), printMatrix(), resize(), setZero(), and solveAndInvertBorderedBand().
|
private |
Mixed part.
Definition at line 103 of file BorderedBandMatrix.h.
Referenced by addBlockMatrix(), getBlockMatrix(), printMatrix(), resize(), setZero(), and solveAndInvertBorderedBand().