FflasFfpack
|
Finite Field PACK Set of elimination based routines for dense linear algebra. More...
Namespaces | |
namespace | Protected |
Data Structures | |
class | CharpolyFailed |
class | callLUdivine_small |
class | callLUdivine_small< double > |
class | callLUdivine_small< float > |
class | UnparametricField |
class | Modular |
class | ModularBalanced |
class | ModularBalanced< double > |
it is forbiden to have char 2 More... | |
class | ModularBalanced< float > |
class | ModularBalanced< int32_t > |
class | ModularBalanced< int64_t > |
class | Modular< double > |
class | Modular< float > |
class | Modular< int32_t > |
Specialization of Modular to int32_t element type with efficient dot product. More... | |
class | Modular< int64_t > |
Specialization of Modular to int64_t element type with efficient dot product. More... | |
class | ModularRandIter |
class | ModularBalancedRandIter |
class | NonzeroRandIter |
Random iterator for nonzero random numbers. More... | |
class | UnparametricOperations |
class | Failure |
A precondtion failed. More... | |
Enumerations | |
enum | FFPACK_LUDIVINE_TAG { FfpackLQUP =1, FfpackSingular =2 } |
enum | FFPACK_CHARPOLY_TAG { FfpackLUK =1, FfpackKG =2, FfpackHybrid =3, FfpackKGFast =4, FfpackDanilevski =5, FfpackArithProg =6, FfpackKGFastG =7 } |
enum | FFPACK_MINPOLY_TAG { FfpackDense =1, FfpackKGF =2 } |
Functions | |
template<class Field > | |
void | applyP (const Field &F, const FFLAS::FFLAS_SIDE Side, const FFLAS::FFLAS_TRANSPOSE Trans, const size_t M, const int ibeg, const int iend, typename Field::Element *A, const size_t lda, const size_t *P) |
Apply a permutation submatrix of P (between ibeg and iend) to a matrix to (iend-ibeg) vectors of size M stored in A (as column for NoTrans and rows for Trans). More... | |
template<class Field > | |
size_t | Rank (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda) |
Computes the rank of the given matrix using a LQUP factorization. More... | |
template<class Field > | |
bool | IsSingular (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda) |
Returns true if the given matrix is singular. More... | |
template<class Field > | |
Field::Element | Det (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda) |
Returns the determinant of the given matrix. More... | |
template<class Field > | |
void | solveLB2 (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, const size_t R, typename Field::Element *L, const size_t ldl, const size_t *Q, typename Field::Element *B, const size_t ldb) |
Solve L X = B in place. More... | |
template<class Field > | |
void | fgetrs (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, const size_t R, typename Field::Element *A, const size_t lda, const size_t *P, const size_t *Q, typename Field::Element *B, const size_t ldb, int *info) |
Solve the system ![]() ![]() | |
template<class Field > | |
Field::Element * | fgetrs (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, const size_t NRHS, const size_t R, typename Field::Element *A, const size_t lda, const size_t *P, const size_t *Q, typename Field::Element *X, const size_t ldx, const typename Field::Element *B, const size_t ldb, int *info) |
Solve the system A X = B or X A = B. More... | |
template<class Field > | |
size_t | fgesv (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *B, const size_t ldb, int *info) |
Square system solver. More... | |
template<class Field > | |
size_t | fgesv (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, const size_t NRHS, typename Field::Element *A, const size_t lda, typename Field::Element *X, const size_t ldx, const typename Field::Element *B, const size_t ldb, int *info) |
Rectangular system solver. More... | |
template<class Field > | |
Field::Element * | Solve (const Field &F, const size_t M, typename Field::Element *A, const size_t lda, typename Field::Element *x, const int incx, const typename Field::Element *b, const int incb) |
Solve the system Ax=b. More... | |
template<class Field > | |
size_t | NullSpaceBasis (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *&NS, size_t &ldn, size_t &NSdim) |
Computes a basis of the Left/Right nullspace of the matrix A. More... | |
template<class Field > | |
size_t | RowRankProfile (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *&rkprofile) |
Computes the row rank profile of A. More... | |
template<class Field > | |
size_t | ColumnRankProfile (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *&rkprofile) |
Computes the column rank profile of A. More... | |
template<class Field > | |
size_t | RowRankProfileSubmatrixIndices (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *&rowindices, size_t *&colindices, size_t &R) |
RowRankProfileSubmatrixIndices. More... | |
template<class Field > | |
size_t | ColRankProfileSubmatrixIndices (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *&rowindices, size_t *&colindices, size_t &R) |
Computes the indices of the submatrix r*r X of A whose columns correspond to the column rank profile of A. More... | |
template<class Field > | |
size_t | RowRankProfileSubmatrix (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *&X, size_t &R) |
Compute the r*r submatrix X of A, by picking the row rank profile rows of A. More... | |
template<class Field > | |
size_t | ColRankProfileSubmatrix (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *&X, size_t &R) |
Compute the ![]() | |
template<class Field > | |
Field::Element * | LQUPtoInverseOfFullRankMinor (const Field &F, const size_t rank, typename Field::Element *A_factors, const size_t lda, const size_t *QtPointer, typename Field::Element *X, const size_t ldx) |
LQUPtoInverseOfFullRankMinor. More... | |
template<class Field > | |
size_t | TURBO (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Q, const size_t cutoff) |
template<class Field > | |
size_t | LUdivine (const Field &F, const FFLAS::FFLAS_DIAG Diag, const FFLAS::FFLAS_TRANSPOSE trans, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Qt, const FFPACK_LUDIVINE_TAG LuTag=FfpackLQUP, const size_t cutoff=__FFPACK_LUDIVINE_CUTOFF) |
Compute the LQUP factorization of the given matrix. More... | |
template<class Field > | |
size_t | LUpdate (const Field &F, const FFLAS::FFLAS_DIAG Diag, const FFLAS::FFLAS_TRANSPOSE trans, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, const size_t R, const size_t K, typename Field::Element *B, const size_t ldb, size_t *P, size_t *Q, const FFPACK::FFPACK_LUDIVINE_TAG LuTag=FFPACK::FfpackLQUP, const size_t cutoff=__FFPACK_LUDIVINE_CUTOFF) |
LUpdate. More... | |
template<class Field > | |
size_t | LUdivine_small (const Field &F, const FFLAS::FFLAS_DIAG Diag, const FFLAS::FFLAS_TRANSPOSE trans, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Q, const FFPACK_LUDIVINE_TAG LuTag=FfpackLQUP) |
LUdivine small case. More... | |
template<class Field > | |
size_t | LUdivine_gauss (const Field &F, const FFLAS::FFLAS_DIAG Diag, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Q, const FFPACK_LUDIVINE_TAG LuTag=FfpackLQUP) |
LUdivine gauss. More... | |
template<class Field > | |
void | ftrtri (const Field &F, const FFLAS::FFLAS_UPLO Uplo, const FFLAS::FFLAS_DIAG Diag, const size_t N, typename Field::Element *A, const size_t lda) |
Compute the inverse of a triangular matrix. More... | |
template<class Field > | |
void | ftrtrm (const Field &F, const FFLAS::FFLAS_DIAG diag, const size_t N, typename Field::Element *A, const size_t lda) |
Compute the product UL. More... | |
template<class Field > | |
size_t | ColumnEchelonForm (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Qt, bool transform=true) |
Compute the Column Echelon form of the input matrix in-place. More... | |
template<class Field > | |
size_t | RowEchelonForm (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Qt, const bool transform=false) |
Compute the Row Echelon form of the input matrix in-place. More... | |
template<class Field > | |
size_t | ReducedColumnEchelonForm (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Qt, const bool transform=true) |
Compute the Reduced Column Echelon form of the input matrix in-place. More... | |
template<class Field > | |
size_t | ReducedRowEchelonForm (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Qt, const bool transform=true) |
Compute the Reduced Row Echelon form of the input matrix in-place. More... | |
template<class Field > | |
size_t | ReducedRowEchelonForm2 (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Qt, const bool transform=true) |
Variant by the block recursive algorithm. More... | |
template<class Field > | |
size_t | REF (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, const size_t colbeg, const size_t rowbeg, const size_t colsize, size_t *Qt, size_t *P) |
REF. More... | |
template<class Field > | |
Field::Element * | Invert (const Field &F, const size_t M, typename Field::Element *A, const size_t lda, int &nullity) |
Invert the given matrix in place or computes its nullity if it is singular. More... | |
template<class Field > | |
Field::Element * | Invert (const Field &F, const size_t M, const typename Field::Element *A, const size_t lda, typename Field::Element *X, const size_t ldx, int &nullity) |
Invert the given matrix in place or computes its nullity if it is singular. More... | |
template<class Field > | |
Field::Element * | Invert2 (const Field &F, const size_t M, typename Field::Element *A, const size_t lda, typename Field::Element *X, const size_t ldx, int &nullity) |
Invert the given matrix or computes its nullity if it is singular. More... | |
template<class Field , class Polynomial > | |
std::list< Polynomial > & | CharPoly (const Field &F, std::list< Polynomial > &charp, const size_t N, typename Field::Element *A, const size_t lda, const FFPACK_CHARPOLY_TAG CharpTag=FfpackArithProg) |
Compute the characteristic polynomial of A using Krylov Method, and LUP factorization of the Krylov matrix. More... | |
template<class Polynomial , class Field > | |
Polynomial & | mulpoly (const Field &F, Polynomial &res, const Polynomial &P1, const Polynomial &P2) |
template<class Field , class Polynomial > | |
std::list< Polynomial > & | CharPoly (const Field &F, Polynomial &charp, const size_t N, typename Field::Element *A, const size_t lda, const FFPACK_CHARPOLY_TAG CharpTag=FfpackArithProg) |
template<class Field , class Polynomial > | |
Polynomial & | MinPoly (const Field &F, Polynomial &minP, const size_t N, const typename Field::Element *A, const size_t lda, typename Field::Element *X, const size_t ldx, size_t *P, const FFPACK::FFPACK_MINPOLY_TAG MinTag=FFPACK::FfpackDense, const size_t kg_mc=0, const size_t kg_mb=0, const size_t kg_j=0) |
Compute the minimal polynomial. More... | |
template<class Field > | |
void | solveLB (const Field &F, const FFLAS::FFLAS_SIDE Side, const size_t M, const size_t N, const size_t R, typename Field::Element *L, const size_t ldl, const size_t *Q, typename Field::Element *B, const size_t ldb) |
Solve L X = B or X L = B in place. More... | |
template<class Field > | |
void | trinv_left (const Field &F, const size_t N, const typename Field::Element *L, const size_t ldl, typename Field::Element *X, const size_t ldx) |
template<class Field > | |
size_t | KrylovElim (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, size_t *P, size_t *Q, const size_t deg, size_t *iterates, size_t *inviterates, const size_t maxit, size_t virt) |
template<class Field > | |
size_t | SpecRankProfile (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, const size_t deg, size_t *rankProfile) |
template<class Field , class Polynomial > | |
std::list< Polynomial > & | CharpolyArithProg (const Field &F, std::list< Polynomial > &frobeniusForm, const size_t N, typename Field::Element *A, const size_t lda, const size_t c) |
template<class Field > | |
void | CompressRows (Field &F, const size_t M, typename Field::Element *A, const size_t lda, typename Field::Element *tmp, const size_t ldtmp, const size_t *d, const size_t nb_blocs) |
template<class Field > | |
void | CompressRowsQK (Field &F, const size_t M, typename Field::Element *A, const size_t lda, typename Field::Element *tmp, const size_t ldtmp, const size_t *d, const size_t deg, const size_t nb_blocs) |
template<class Field > | |
void | DeCompressRows (Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *tmp, const size_t ldtmp, const size_t *d, const size_t nb_blocs) |
template<class Field > | |
void | DeCompressRowsQK (Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *tmp, const size_t ldtmp, const size_t *d, const size_t deg, const size_t nb_blocs) |
template<class Field > | |
void | CompressRowsQA (Field &F, const size_t M, typename Field::Element *A, const size_t lda, typename Field::Element *tmp, const size_t ldtmp, const size_t *d, const size_t nb_blocs) |
template<class Field > | |
void | DeCompressRowsQA (Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda, typename Field::Element *tmp, const size_t ldtmp, const size_t *d, const size_t nb_blocs) |
template<class Field , class Polynomial > | |
std::list< Polynomial > & | Danilevski (const Field &F, std::list< Polynomial > &charp, const size_t N, typename Field::Element *A, const size_t lda) |
template<class Field > | |
Field::Element * | buildMatrix (const Field &F, const typename Field::Element *E, const typename Field::Element *C, const size_t lda, const size_t *B, const size_t *T, const size_t me, const size_t mc, const size_t lambda, const size_t mu) |
template<class Field , class Polynomial > | |
Polynomial & | MinPoly (const Field &F, Polynomial &minP, const size_t N, const typename Field::Element *A, const size_t lda, typename Field::Element *U, size_t ldu, typename Field::Element *X, size_t ldx, size_t *P) |
template<class T > | |
std::ostream & | operator<< (std::ostream &o, const Modular< T > &F) |
template<class T > | |
std::ostream & | operator<< (std::ostream &o, const ModularBalanced< T > &F) |
template<class T > | |
std::ostream & | operator<< (std::ostream &o, const UnparametricField< T > &F) |
template<typename Target , typename Source > | |
Target & | Caster (Target &t, const Source &s) |
template<class Field > | |
Field::Element * | RandomMatrix (const Field &F, typename Field::Element *A, size_t m, size_t n, size_t lda) |
Random Matrix. More... | |
int | RandInt (int a, int b) |
Random integer in range. More... | |
size_t | RandInt (size_t a, size_t b) |
Random integer in range. More... | |
template<class Field > | |
Field::Element * | RandomMatrixWithRank (const Field &F, typename Field::Element *A, size_t r, size_t m, size_t n, size_t lda) |
Random Matrix with prescribed rank. More... | |
Finite Field PACK Set of elimination based routines for dense linear algebra.
This namespace enlarges the set of BLAS routines of the class FFLAS, with higher level routines based on elimination.
enum FFPACK_LUDIVINE_TAG |
enum FFPACK_CHARPOLY_TAG |
enum FFPACK_MINPOLY_TAG |
void FFPACK::applyP | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const FFLAS::FFLAS_TRANSPOSE | Trans, | ||
const size_t | M, | ||
const int | ibeg, | ||
const int | iend, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t * | P | ||
) |
Apply a permutation submatrix of P (between ibeg and iend) to a matrix to (iend-ibeg) vectors of size M stored in A (as column for NoTrans and rows for Trans).
Side==FFLAS::FflasLeft for row permutation Side==FFLAS::FflasRight for a column permutation Trans==FFLAS::FflasTrans for the inverse permutation of P
F | |
Side | |
Trans | |
M | |
ibeg | |
iend | |
A | |
lda | |
P |
size_t FFPACK::Rank | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda | ||
) |
Computes the rank of the given matrix using a LQUP factorization.
The input matrix is modified.
F | field |
M | row dimension of the matrix |
N | column dimension of the matrix |
A | input matrix |
lda | leading dimension of A |
bool FFPACK::IsSingular | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda | ||
) |
Returns true if the given matrix is singular.
The method is a block elimination with early termination
using LQUP factorization with early termination.
F | field |
M | row dimension of the matrix |
N | column dimension of the matrix |
A | input matrix |
lda | leading dimension of A |
Field::Element FFPACK::Det | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda | ||
) |
Returns the determinant of the given matrix.
The method is a block elimination with early termination
F | field |
M | row dimension of the matrix |
N | column dimension of the matrix |
A | input matrix |
lda | leading dimension of Ausing LQUP factorization with early termination. |
void solveLB2 | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
const size_t | R, | ||
typename Field::Element * | L, | ||
const size_t | ldl, | ||
const size_t * | Q, | ||
typename Field::Element * | B, | ||
const size_t | ldb | ||
) |
Solve L X = B in place.
L is M*M or N*N, B is M*N. Only the R non trivial column of L are stored in the M*R matrix L
void FFPACK::fgetrs | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
const size_t | R, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t * | P, | ||
const size_t * | Q, | ||
typename Field::Element * | B, | ||
const size_t | ldb, | ||
int * | info | ||
) |
Solve the system or
.
Solving using the LQUP
decomposition of A
already computed inplace with LUdivine(FFLAS::FflasNoTrans, FFLAS::FflasNonUnit)
. Version for A square. If A is rank deficient, a solution is returned if the system is consistent, Otherwise an info is 1
F | field |
Side | Determine wheter the resolution is left or right looking. |
M | row dimension of B |
N | col dimension of B |
R | rank of A |
A | input matrix |
lda | leading dimension of A |
P | column permutation of the LQUP decomposition of A |
Q | column permutation of the LQUP decomposition of A |
B | Right/Left hand side matrix. Initially stores B , finally stores the solution X . |
ldb | leading dimension of B |
info | Success of the computation: 0 if successfull, >0 if system is inconsistent |
Field::Element* FFPACK::fgetrs | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
const size_t | NRHS, | ||
const size_t | R, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t * | P, | ||
const size_t * | Q, | ||
typename Field::Element * | X, | ||
const size_t | ldx, | ||
const typename Field::Element * | B, | ||
const size_t | ldb, | ||
int * | info | ||
) |
Solve the system A X = B or X A = B.
Solving using the LQUP decomposition of A already computed inplace with LUdivine(FFLAS::FflasNoTrans, FFLAS::FflasNonUnit). Version for A rectangular. If A is rank deficient, a solution is returned if the system is consistent, Otherwise an info is 1
F | field |
Side | Determine wheter the resolution is left or right looking. |
M | row dimension of A |
N | col dimension of A |
NRHS | number of columns (if Side = FFLAS::FflasLeft) or row (if Side = FFLAS::FflasRight) of the matrices X and B |
R | rank of A |
A | input matrix |
lda | leading dimension of A |
P | column permutation of the LQUP decomposition of A |
Q | column permutation of the LQUP decomposition of A |
X | solution matrix |
ldx | leading dimension of X |
B | Right/Left hand side matrix. |
ldb | leading dimension of B |
info | Succes of the computation: 0 if successfull, >0 if system is inconsistent |
size_t FFPACK::fgesv | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | B, | ||
const size_t | ldb, | ||
int * | info | ||
) |
Square system solver.
F | The computation domain |
Side | Determine wheter the resolution is left or right looking |
M | row dimension of B |
N | col dimension of B |
A | input matrix |
lda | leading dimension of A |
B | Right/Left hand side matrix. Initially contains B, finally contains the solution X. |
ldb | leading dimension of B |
info | Success of the computation: 0 if successfull, >0 if system is inconsistent |
Solve the system A X = B or X A = B. Version for A square. If A is rank deficient, a solution is returned if the system is consistent, Otherwise an info is 1
size_t FFPACK::fgesv | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
const size_t | NRHS, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | X, | ||
const size_t | ldx, | ||
const typename Field::Element * | B, | ||
const size_t | ldb, | ||
int * | info | ||
) |
Rectangular system solver.
F | The computation domain |
Side | Determine wheter the resolution is left or right looking |
M | row dimension of A |
N | col dimension of A |
NRHS | number of columns (if Side = FFLAS::FflasLeft) or row (if Side = FFLAS::FflasRight) of the matrices X and B |
A | input matrix |
lda | leading dimension of A |
B | Right/Left hand side matrix. Initially contains B, finally contains the solution X. |
ldb | leading dimension of B |
X | |
ldx | |
info | Success of the computation: 0 if successfull, >0 if system is inconsistent |
Solve the system A X = B or X A = B. Version for A square. If A is rank deficient, a solution is returned if the system is consistent, Otherwise an info is 1
Field::Element* FFPACK::Solve | ( | const Field & | F, |
const size_t | M, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | x, | ||
const int | incx, | ||
const typename Field::Element * | b, | ||
const int | incb | ||
) |
Solve the system Ax=b.
Solving using LQUP factorization and two triangular system resolutions. The input matrix is modified.
F | The computation domain |
M | row dimension of the matrix |
A | input matrix |
lda | leading dimension of A |
x | solution vector |
incx | increment of x |
b | right hand side vector |
incb | increment of bSolve linear system using LQUP factorization. |
size_t FFPACK::NullSpaceBasis | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element *& | NS, | ||
size_t & | ldn, | ||
size_t & | NSdim | ||
) |
Computes a basis of the Left/Right nullspace of the matrix A.
return the dimension of the nullspace.
F | The computation domain |
Side | |
M | |
N | |
A | input matrix of dimension M x N, A is modified |
lda | |
NS | output matrix of dimension N x NSdim (allocated here) |
ldn | |
NSdim | the dimension of the Nullspace (N-rank(A)) |
size_t FFPACK::RowRankProfile | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t *& | rkprofile | ||
) |
Computes the row rank profile of A.
F | |
M | |
N | |
A | input matrix of dimension M x N |
lda | |
rkprofile | return the rank profile as an array of row indexes, of dimension r=rank(A) |
rkprofile is allocated during the computation.
size_t FFPACK::ColumnRankProfile | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t *& | rkprofile | ||
) |
Computes the column rank profile of A.
F | |
M | |
N | |
A | input matrix of dimension |
lda | |
rkprofile | return the rank profile as an array of row indexes, of dimension r=rank(A) |
A is modified rkprofile is allocated during the computation.
size_t FFPACK::RowRankProfileSubmatrixIndices | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t *& | rowindices, | ||
size_t *& | colindices, | ||
size_t & | R | ||
) |
RowRankProfileSubmatrixIndices.
Computes the indices of the submatrix r*r X of A whose rows correspond to the row rank profile of A.
F | ||
M | ||
N | ||
A | input matrix of dimension | |
rowindices | array of the row indices of X in A | |
colindices | array of the col indices of X in A | |
lda | ||
[out] | R | rowindices and colindices are allocated during the computation. A is modified |
size_t FFPACK::ColRankProfileSubmatrixIndices | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t *& | rowindices, | ||
size_t *& | colindices, | ||
size_t & | R | ||
) |
Computes the indices of the submatrix r*r X of A whose columns correspond to the column rank profile of A.
F | ||
M | ||
N | ||
A | input matrix of dimension | |
rowindices | array of the row indices of X in A | |
colindices | array of the col indices of X in A | |
lda | ||
[out] | R | rowindices and colindices are allocated during the computation. |
size_t FFPACK::RowRankProfileSubmatrix | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element *& | X, | ||
size_t & | R | ||
) |
Compute the r*r submatrix X of A, by picking the row rank profile rows of A.
F | ||
M | ||
N | ||
A | input matrix of dimension M x N | |
lda | ||
X | the output matrix | |
[out] | R | A is not modified X is allocated during the computation. |
size_t FFPACK::ColRankProfileSubmatrix | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element *& | X, | ||
size_t & | R | ||
) |
Compute the submatrix X of A, by picking the row rank profile rows of A.
F | ||
M | ||
N | ||
A | input matrix of dimension M x N | |
lda | ||
X | the output matrix | |
[out] | R | A is not modified X is allocated during the computation. |
Field::Element* FFPACK::LQUPtoInverseOfFullRankMinor | ( | const Field & | F, |
const size_t | rank, | ||
typename Field::Element * | A_factors, | ||
const size_t | lda, | ||
const size_t * | QtPointer, | ||
typename Field::Element * | X, | ||
const size_t | ldx | ||
) |
LQUPtoInverseOfFullRankMinor.
Suppose A has been factorized as L.Q.U.P, with rank r. Then Qt.A.Pt has an invertible leading principal r x r submatrix This procedure efficiently computes the inverse of this minor and puts it into X.
F | |
rank | rank of the matrix. |
A_factors | matrix containing the L and U entries of the factorization |
lda | |
QtPointer | theLQUP->getQ()->getPointer() (note: getQ returns Qt!) |
X | desired location for output |
ldx |
|
inline |
!!!!! Attention a appliquer Q4, Q2, Q3, Q3b a gauche !!!!!!!
!!!!! Attention a appliquer Q4, Q2, Q3, Q3b a gauche !!!!!!!
|
inline |
Compute the LQUP factorization of the given matrix.
Using a block algorithm and return its rank. The permutations P and Q are represented using LAPACK's convention.
F | field |
Diag | precise whether U should have a unit diagonal or not |
trans | UNKOWN TAG, probably the LU of ![]() |
M | matrix row dimension |
N | matrix column dimension |
A | input matrix |
lda | leading dimension of A |
P | the column permutation |
Qt | the transpose of the row permutation Q |
LuTag | flag for setting the earling termination if the matrix is singular |
cutoff | UNKOWN TAG, probably a switch to a faster algo below cutoff |
A
LSP
Matrix Decomposition Revisited, 2006LUdivine
, une divine factorisation LU
, 2002size_t LUpdate | ( | const Field & | F, |
const FFLAS::FFLAS_DIAG | Diag, | ||
const FFLAS::FFLAS_TRANSPOSE | trans, | ||
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t | R, | ||
const size_t | K, | ||
typename Field::Element * | B, | ||
const size_t | ldb, | ||
size_t * | P, | ||
size_t * | Q, | ||
const FFPACK::FFPACK_LUDIVINE_TAG | LuTag, | ||
const size_t | cutoff | ||
) |
LUpdate.
Updates an existing LU factorisation with more rows.
F | Field on which arithmetic is done |
Diag | Is L unit ? (if so, FFLAS::FflasUnit ) |
trans | Not used yet, should be FFLAS::FflasNoTrans |
M | rows in A |
N | cols in A |
A | A is already in LU factorisation |
lda | leading dimension of A |
R | rank of A |
K | rows in B |
B | more rows to append to A |
ldb | leading dimension of B (not tested if != lda) |
P | permutation for LU in A . Should be big enough so it can store the permutation for LU of A and B |
Q | same as P |
LuTag | see LUdivine . |
cutoff | see LUdivine . |
A.append(B)
|
inline |
LUdivine small case.
|
inline |
LUdivine gauss.
void FFPACK::ftrtri | ( | const Field & | F, |
const FFLAS::FFLAS_UPLO | Uplo, | ||
const FFLAS::FFLAS_DIAG | Diag, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda | ||
) |
Compute the inverse of a triangular matrix.
F | |
Uplo | whether the matrix is upper of lower triangular |
Diag | whether the matrix if unit diagonal |
N | |
A | |
lda |
void FFPACK::ftrtrm | ( | const Field & | F, |
const FFLAS::FFLAS_DIAG | diag, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda | ||
) |
Compute the product UL.
Product UL of the upper, resp lower triangular matrices U and L stored one above the other in the square matrix A. Diag == Unit if the matrix U is unit diagonal
F | |
diag | |
N | |
A | |
lda |
size_t ColumnEchelonForm | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t * | P, | ||
size_t * | Qt, | ||
bool | transform = true |
||
) |
Compute the Column Echelon form of the input matrix in-place.
After the computation A = [ M \ V ] such that AU = C is a column echelon decomposition of A, with U = P^T [ V ] and C = M + Q [ Ir ] [ 0 In-r ] [ 0 ] Qt = Q^T If transform=false, the matrix U is not computed. See also test-colechelon for an example of use
F | |
M | |
N | |
A | |
lda | |
P | |
Qt | |
transform |
size_t RowEchelonForm | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t * | P, | ||
size_t * | Qt, | ||
const bool | transform = false |
||
) |
Compute the Row Echelon form of the input matrix in-place.
After the computation A = [ L \ M ] such that L A = R is a row echelon decomposition of A, with L = [ L 0 ] P and R = M + [Ir 0] Q^T [ In-r] Qt = Q^T If transform=false, the matrix L is not computed. See also test-rowechelon for an example of use
F | |
M | |
N | |
A | |
lda | |
P | |
Qt | |
transform |
size_t ReducedColumnEchelonForm | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t * | P, | ||
size_t * | Qt, | ||
const bool | transform = true |
||
) |
Compute the Reduced Column Echelon form of the input matrix in-place.
After the computation A = [ V ] such that AU = R is a reduced col echelon [ M 0 ] decomposition of A, where U = P^T [ V ] and R = Q [ Ir ] [ 0 In-r ] [ M 0 ] Qt = Q^T If transform=false, the matrix U is not computed and the matrix A = R
F | |
M | |
N | |
A | |
lda | |
P | |
Qt | |
transform |
size_t ReducedRowEchelonForm | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t * | P, | ||
size_t * | Qt, | ||
const bool | transform = true |
||
) |
Compute the Reduced Row Echelon form of the input matrix in-place.
After the computation A = [ V1 M ] such that L A = R is a reduced row echelon [ V2 0 ] decomposition of A, where L = [ V1 0 ] P and R = [ Ir M ] Q^T [ V2 In-r ] [ 0 ] Qt = Q^T If transform=false, the matrix U is not computed and the matrix A = R
F | |
M | |
N | |
A | |
lda | |
P | |
Qt | |
transform |
size_t FFPACK::ReducedRowEchelonForm2 | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
size_t * | P, | ||
size_t * | Qt, | ||
const bool | transform = true |
||
) |
Variant by the block recursive algorithm.
(See A. Storjohann Thesis 2000) !!!!!! Warning !!!!!! This code is NOT WORKING properly for some echelon structures. This is due to a limitation of the way we represent permutation matrices (LAPACK's standard):
size_t REF | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t | colbeg, | ||
const size_t | rowbeg, | ||
const size_t | colsize, | ||
size_t * | Qt, | ||
size_t * | P | ||
) |
REF.
| I | A11 | A12 | | |-—|--—|--—|-—| | |I | *| A22 | | | |0 | 0| A22 | | |-—|--—|--—|-—| | | 0 | A32 | | |-—|--—|--—|-—|
where the transformation matrix is stored at the pivot column position
Field::Element* FFPACK::Invert | ( | const Field & | F, |
const size_t | M, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
int & | nullity | ||
) |
Invert the given matrix in place or computes its nullity if it is singular.
An inplace 2n^3 algorithm is used.
F | The computation domain | |
M | order of the matrix | |
[in,out] | A | input matrix ( ![]() |
lda | leading dimension of A | |
nullity | dimension of the kernel of A |
Field::Element* FFPACK::Invert | ( | const Field & | F, |
const size_t | M, | ||
const typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | X, | ||
const size_t | ldx, | ||
int & | nullity | ||
) |
Invert the given matrix in place or computes its nullity if it is singular.
X is preallocated.
F | The computation domain | |
M | order of the matrix | |
[in] | A | input matrix ( ![]() |
lda | leading dimension of A | |
[out] | X | output matrix |
ldx | leading dimension of X | |
nullity | dimension of the kernel of A |
Field::Element* FFPACK::Invert2 | ( | const Field & | F, |
const size_t | M, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | X, | ||
const size_t | ldx, | ||
int & | nullity | ||
) |
Invert the given matrix or computes its nullity if it is singular.
An 2n^3 algorithm is used. This routine can be % faster than Invert but is not totally inplace. X is preallocated.
F | ||
M | order of the matrix | |
[in,out] | A | input matrix ( ![]() |
lda | leading dimension of A | |
[out] | X | output matrix |
ldx | leading dimension of X | |
nullity | dimension of the kernel of A |
std::list< Polynomial > & CharPoly | ( | const Field & | F, |
std::list< Polynomial > & | charp, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const FFPACK_CHARPOLY_TAG | CharpTag = FfpackArithProg |
||
) |
Compute the characteristic polynomial of A using Krylov Method, and LUP factorization of the Krylov matrix.
Polynomial& FFPACK::mulpoly | ( | const Field & | F, |
Polynomial & | res, | ||
const Polynomial & | P1, | ||
const Polynomial & | P2 | ||
) |
std::list<Polynomial>& FFPACK::CharPoly | ( | const Field & | F, |
Polynomial & | charp, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const FFPACK_CHARPOLY_TAG | CharpTag = FfpackArithProg |
||
) |
Polynomial & MinPoly | ( | const Field & | F, |
Polynomial & | minP, | ||
const size_t | N, | ||
const typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | X, | ||
const size_t | ldx, | ||
size_t * | P, | ||
const FFPACK::FFPACK_MINPOLY_TAG | MinTag = FFPACK::FfpackDense , |
||
const size_t | kg_mc = 0 , |
||
const size_t | kg_mb = 0 , |
||
const size_t | kg_j = 0 |
||
) |
Compute the minimal polynomial.
Minpoly of (A,v) using an LUP factorization of the Krylov Base (v, Av, .., A^kv) U,X must be (n+1)*n U contains the Krylov matrix and X, its LSP factorization
void FFPACK::solveLB | ( | const Field & | F, |
const FFLAS::FFLAS_SIDE | Side, | ||
const size_t | M, | ||
const size_t | N, | ||
const size_t | R, | ||
typename Field::Element * | L, | ||
const size_t | ldl, | ||
const size_t * | Q, | ||
typename Field::Element * | B, | ||
const size_t | ldb | ||
) |
Solve L X = B or X L = B in place.
L is M*M if Side == FFLAS::FflasLeft and N*N if Side == FFLAS::FflasRight, B is M*N. Only the R non trivial column of L are stored in the M*R matrix L Requirement : so that L could be expanded in-place
void FFPACK::trinv_left | ( | const Field & | F, |
const size_t | N, | ||
const typename Field::Element * | L, | ||
const size_t | ldl, | ||
typename Field::Element * | X, | ||
const size_t | ldx | ||
) |
|
inline |
size_t SpecRankProfile | ( | const Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t | deg, | ||
size_t * | rankProfile | ||
) |
std::list< Polynomial > & CharpolyArithProg | ( | const Field & | F, |
std::list< Polynomial > & | frobeniusForm, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
const size_t | c | ||
) |
void CompressRows | ( | Field & | F, |
const size_t | M, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | tmp, | ||
const size_t | ldtmp, | ||
const size_t * | d, | ||
const size_t | nb_blocs | ||
) |
void CompressRowsQK | ( | Field & | F, |
const size_t | M, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | tmp, | ||
const size_t | ldtmp, | ||
const size_t * | d, | ||
const size_t | deg, | ||
const size_t | nb_blocs | ||
) |
void DeCompressRows | ( | Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | tmp, | ||
const size_t | ldtmp, | ||
const size_t * | d, | ||
const size_t | nb_blocs | ||
) |
void DeCompressRowsQK | ( | Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | tmp, | ||
const size_t | ldtmp, | ||
const size_t * | d, | ||
const size_t | deg, | ||
const size_t | nb_blocs | ||
) |
void CompressRowsQA | ( | Field & | F, |
const size_t | M, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | tmp, | ||
const size_t | ldtmp, | ||
const size_t * | d, | ||
const size_t | nb_blocs | ||
) |
void DeCompressRowsQA | ( | Field & | F, |
const size_t | M, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | tmp, | ||
const size_t | ldtmp, | ||
const size_t * | d, | ||
const size_t | nb_blocs | ||
) |
std::list<Polynomial>& FFPACK::Danilevski | ( | const Field & | F, |
std::list< Polynomial > & | charp, | ||
const size_t | N, | ||
typename Field::Element * | A, | ||
const size_t | lda | ||
) |
Field::Element* FFPACK::buildMatrix | ( | const Field & | F, |
const typename Field::Element * | E, | ||
const typename Field::Element * | C, | ||
const size_t | lda, | ||
const size_t * | B, | ||
const size_t * | T, | ||
const size_t | me, | ||
const size_t | mc, | ||
const size_t | lambda, | ||
const size_t | mu | ||
) |
Polynomial& FFPACK::MinPoly | ( | const Field & | F, |
Polynomial & | minP, | ||
const size_t | N, | ||
const typename Field::Element * | A, | ||
const size_t | lda, | ||
typename Field::Element * | U, | ||
size_t | ldu, | ||
typename Field::Element * | X, | ||
size_t | ldx, | ||
size_t * | P | ||
) |
std::ostream& FFPACK::operator<< | ( | std::ostream & | o, |
const Modular< T > & | F | ||
) |
std::ostream& FFPACK::operator<< | ( | std::ostream & | o, |
const ModularBalanced< T > & | F | ||
) |
std::ostream& FFPACK::operator<< | ( | std::ostream & | o, |
const UnparametricField< T > & | F | ||
) |
Target& FFPACK::Caster | ( | Target & | t, |
const Source & | s | ||
) |
Field::Element* FFPACK::RandomMatrix | ( | const Field & | F, |
typename Field::Element * | A, | ||
size_t | m, | ||
size_t | n, | ||
size_t | lda | ||
) |
Random Matrix.
Creates a \c m x \c n matrix with random entries.
F | field |
A | pointer to the matrix (preallocated to at least m x lda field elements) |
m | number of rows in A |
n | number of cols in A |
lda | leading dimension of A |
A
. int FFPACK::RandInt | ( | int | a, |
int | b | ||
) |
Random integer in range.
a | min bound |
b | max bound |
size_t FFPACK::RandInt | ( | size_t | a, |
size_t | b | ||
) |
Random integer in range.
a | min bound |
b | max bound |
Field::Element* FFPACK::RandomMatrixWithRank | ( | const Field & | F, |
typename Field::Element * | A, | ||
size_t | r, | ||
size_t | m, | ||
size_t | n, | ||
size_t | lda | ||
) |
Random Matrix with prescribed rank.
Creates a \c m x \c n matrix with random entries and rank \c r.
F | field |
A | pointer to the matrix (preallocated to at least m x lda field elements) |
r | rank of the matrix to build |
m | number of rows in A |
n | number of cols in A |
lda | leading dimension of A |
A
.