FFPACK Class Reference
[linbox/ffpack]

Set of elimination based routines for dense linear algebra with matrices over finite prime field of characteristic less than 2^26. More...

#include <ffpack.h>

Inherits FFLAS.

List of all members.

Static Public Member Functions

template<class Field >
static size_t Rank (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda)
 using LQUP factorization.
template<class Field >
static bool IsSingular (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda)
 using LQUP factorization with early termination.
template<class Field >
static Field::Element Det (const Field &F, const size_t M, const size_t N, typename Field::Element *A, const size_t lda)
 using LQUP factorization with early termination.
template<class Field >
static void fgetrs (const Field &F, const 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)
template<class Field >
static Field::Element * fgetrs (const Field &F, const 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)
template<class Field >
static size_t fgesv (const Field &F, const 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.
template<class Field >
static size_t fgesv (const Field &F, const 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.
template<class Field >
static 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 linear system using LQUP factorization.
template<class Field >
static Field::Element * Invert (const Field &F, const size_t M, typename Field::Element *A, const size_t lda, int &nullity)
 Invert a matrix or return its nullity.
template<class Field >
static 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 a matrix or return its nullity.
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static size_t LUdivine (const Field &F, const FFLAS_DIAG Diag, const 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)
 LQUP factorization.
template<class Field >
static void ftrtri (const Field &F, const FFLAS_UPLO Uplo, const FFLAS_DIAG Diag, const size_t N, typename Field::Element *A, const size_t lda)
template<class Field >
static void ftrtrm (const Field &F, const FFLAS_DIAG diag, const size_t N, typename Field::Element *A, const size_t lda)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static 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)
template<class Field >
static void applyP (const Field &F, const FFLAS_SIDE Side, const 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)
template<class Field , class Polynomial >
static 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=FfpackLUK)
template<class Field , class Polynomial >
static 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_MINPOLY_TAG MinTag, const size_t kg_mc, const size_t kg_mb, const size_t kg_j)

Detailed Description

Set of elimination based routines for dense linear algebra with matrices over finite prime field of characteristic less than 2^26.

This class only provides a set of static member functions. No instantiation is allowed.

It enlarges the set of BLAS routines of the class FFLAS, with higher level routines based on elimination.


Member Function Documentation

static size_t Rank ( const Field F,
const size_t  M,
const size_t  N,
typename Field::Element *  A,
const size_t  lda 
) [inline, static]

using LQUP factorization.

Computes the rank of the given matrix using a LQUP factorization. The input matrix is modified.

Parameters:
M row dimension of the matrix
N column dimension of the matrix
A input matrix
lda leading dimension of A
static bool IsSingular ( const Field F,
const size_t  M,
const size_t  N,
typename Field::Element *  A,
const size_t  lda 
) [inline, static]

using LQUP factorization with early termination.

Returns true if the given matrix is singular. The method is a block elimination with early termination The input matrix is modified.

Parameters:
M row dimension of the matrix
N column dimension of the matrix
A input matrix
lda leading dimension of A
static Field::Element Det ( const Field F,
const size_t  M,
const size_t  N,
typename Field::Element *  A,
const size_t  lda 
) [inline, static]

using LQUP factorization with early termination.

Returns the determinant of the given matrix. The method is a block elimination with early termination The input matrix is modified.

Parameters:
M row dimension of the matrix
N column dimension of the matrix
A input matrix
lda leading dimension of A
static void fgetrs ( const Field F,
const 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 
) [inline, static]

Solve the system A X = B or X A = B, using the LQUP decomposition of A already computed inplace with LUdivine(FflasNoTrans, FflasNonUnit). Version for A square. If A is rank deficient, a solution is returned if the system is consistent, Otherwise an info is 1

Parameters:
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 Succes of the computation: 0 if successfull, >0 if system is inconsistent
static Field::Element* fgetrs ( const Field F,
const 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 
) [inline, static]

Solve the system A X = B or X A = B, using the LQUP decomposition of A already computed inplace with LUdivine(FflasNoTrans, FflasNonUnit). Version for A rectangular. If A is rank deficient, a solution is returned if the system is consistent, Otherwise an info is 1

Parameters:
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 = FflasLeft) or row (if Side = 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
static size_t fgesv ( const Field F,
const 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 
) [inline, static]

Square system solver.

Parameters:
Field 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
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 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
Returns:
the rank of the system

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

static size_t fgesv ( const Field F,
const 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 
) [inline, static]

Rectangular system solver.

Parameters:
Field 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 = FflasLeft) or row (if Side = FflasRight) of the matrices X and B
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 contains B, finally contains the solution X.
ldb leading dimension of B Success of the computation: 0 if successfull, >0 if system is inconsistent
Returns:
the rank of the system

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

static 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 
) [inline, static]

Solve linear system using LQUP factorization.

Solve the system Ax=b, using LQUP factorization and two triangular system resolutions. The input matrix is modified.

Parameters:
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 b
static Field::Element* Invert ( const Field F,
const size_t  M,
typename Field::Element *  A,
const size_t  lda,
int &  nullity 
) [inline, static]

Invert a matrix or return its nullity.

Invert the given matrix in place or computes its nullity if it is singular. An inplace 2n^3 algorithm is used.

Parameters:
M order of the matrix
A input matrix
lda leading dimension of A
nullity dimension of the kernel of A
static 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 
) [inline, static]

Invert a matrix or return its nullity.

Invert the given matrix or computes its nullity if it is singular. An 2n^3 algorithm is used. The input matrix is modified.

Parameters:
M order of the matrix
A input matrix
lda leading dimension of A
X inverse of A
ldx leading dimension of X
nullity dimension of the kernel of A
static 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 
) [inline, static]

RowRankProfile Computes the row rank profile of A.

Parameters:
A,: input matrix of dimension
rklprofile,: return the rank profile as an array of row indexes, of dimension r=rank(A)

rkprofile is allocated during the computation. Returns R

static 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 
) [inline, static]

ColumnRankProfile Computes the column rank profile of A.

Parameters:
A,: input matrix of dimension
rklprofile,: return the rank profile as an array of row indexes, of dimension r=rank(A)

A is modified rkprofile is allocated during the computation. Returns R

static 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 
) [inline, static]

RowRankProfileSubmatrixIndices Computes the indices of the submatrix r*r X of A whose rows correspond to the row rank profile of A.

Parameters:
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

rowindices and colindices are allocated during the computation. A is modified Returns R

static 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 
) [inline, static]

ColRankProfileSubmatrixIndices Computes the indices of the submatrix r*r X of A whose columns correspond to the column rank profile of A.

Parameters:
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

rowindices and colindices are allocated during the computation. A is modified Returns R

static 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 
) [inline, static]

RowRankProfileSubmatrix Compute the r*r submatrix X of A, by picking the row rank profile rows of A

Parameters:
A,: input matrix of dimension M x N
X,: the output matrix

A is not modified X is allocated during the computation. Returns R

static 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 
) [inline, static]

ColRankProfileSubmatrix Compute the r*r submatrix X of A, by picking the row rank profile rows of A

Parameters:
A,: input matrix of dimension M x N
X,: the output matrix

A is not modified X is allocated during the computation. Returns R

static 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 
) [inline, static]

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. NOTE: It changes the lower entries of A_factors in the process (NB: unless A was nonsingular and square)

Parameters:
rank,: rank of the matrix.
A_factors,: matrix containing the L and U entries of the factorization
QtPointer,: theLQUP->getQ()->getPointer() (note: getQ returns Qt!)
X,: desired location for output
static size_t LUdivine ( const Field F,
const FFLAS_DIAG  Diag,
const 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 
) [inline, static]

LQUP factorization.

Compute the LQUP factorization of the given matrix using a block agorithm and return its rank. The permutations P and Q are represented using LAPACK's convention.

Parameters:
Diag precise whether U should have a unit diagonal or not
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
static void ftrtri ( const Field F,
const FFLAS_UPLO  Uplo,
const FFLAS_DIAG  Diag,
const size_t  N,
typename Field::Element *  A,
const size_t  lda 
) [inline, static]

Compute the inverse of a triangular matrix.

Parameters:
Uplo whether the matrix is upper of lower triangular
Diag whether the matrix if unit diagonal
static void ftrtrm ( const Field F,
const FFLAS_DIAG  diag,
const size_t  N,
typename Field::Element *  A,
const size_t  lda 
) [inline, static]

Compute the 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

static 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 
) [inline, static]

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 + Ir ] and C = M //+ Q [ Ir ] [ 0 In-r ] // [ 0 ] Qt = Q^T

static 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 
) [inline, static]

Compute the Row Echelon form of the input matrix in-place.

After the computation A = [ L \ M ] such that L A = R is a column echelon decomposition of A, with L = [ L+Ir 0 ] P and R = M [ In-r ] Qt = Q^T

static 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 
) [inline, static]

Compute the Reduced Column Echelon form of the input matrix in-place.

After the computation A = [ V ] such that AU = R is a reduced column echelon [ M 0 ] decomposition of A, where U = P^T [ V ] and R = Q [ Ir ] [ 0 In-r ] [ M 0 ] Qt = Q^T

static 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 
) [inline, static]

Compute the Reduced Row Echelon form of the input matrix in-place.

After the computation A = [ V ] such that L A = R is a reduced row echelon [ M 0 ] decomposition of A, where L = [ V ] P^T and R = [ Ir M ] Q [ 0 In-r ] [ 0 ] Qt = Q^T

static void applyP ( const Field F,
const FFLAS_SIDE  Side,
const 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 
) [inline, static]

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==FflasLeft for row permutation Side==FflasRight for a column permutation Trans==FflasTrans for the inverse permutation of P

static 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 = FfpackLUK 
) [inline, static]

Compute the characteristic polynomial of A using Krylov Method, and LUP factorization of the Krylov matrix

static 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_MINPOLY_TAG  MinTag,
const size_t  kg_mc,
const size_t  kg_mb,
const size_t  kg_j 
) [inline, static]

Compute the minimal polynomial 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


The documentation for this class was generated from the following file:

Generated on Thu Oct 1 18:01:06 2009 for linbox by  doxygen 1.6.1