Definition of the lineal algebra routines used in the CuikSuite.
Simple 2D/3D operations are defined in geom.c
Basic vector/matrix operations are defined in basic_algebra.c
Here we define high level linear algebra routines. These are the only routines if we ever change the support linear algebra library (right now GSL is used).
Definition in file algebra.c.
|
void | dgpadm_ (int *ideg, int *m, double *t, double *H, int *ldh, double *wsp, int *lwsp, int *ipiv, int *iexph, int *ns, int *iflag) |
| Interface to the dgpadm function in EXPOKIT. More...
|
|
void | dgchbv_ (int *m, double *t, double *H, int *ldh, double *y, double *wsp, int *iwsp, int *iflag) |
| Interface to the dgchbv function in EXPOKIT. More...
|
|
unsigned int | AnalyzeKernel (unsigned int nr, unsigned int nc, double *mT, unsigned int dof, double epsilon, boolean computeRank, boolean checkRank, boolean getT, boolean getBase, boolean *singular, unsigned int *rank, boolean **IR, double **T) |
| Analyzes the kernel of a matrix. More...
|
|
unsigned int | BasisColumnSpace (double epsilon, unsigned int nr, unsigned int nc, double *A, boolean **IC, unsigned int *nrT, unsigned int *ncT, double **T) |
| Generates an orthonormal basis of the column space of a matrix. More...
|
|
unsigned int | MatrixExponential (unsigned int m, double *A, double t, double *eA) |
| Exponential of a matrix. More...
|
|
unsigned int | MatrixExponentialVector (unsigned int m, double *A, double t, double *v) |
| Exponential of a matrix and product by a vector. More...
|
|
unsigned int | FindRank (double epsilon, unsigned int nr, unsigned int nc, double *mT) |
| Determines the row-rank of a matrix. More...
|
|
unsigned int | FindKernel (unsigned int nr, unsigned int nc, double *mT, unsigned int dof, boolean check, double epsilon, double **T) |
| Computes the kernel of a matrix. More...
|
|
unsigned int | FindKernelAndRank (unsigned int nr, unsigned int nc, double *mT, unsigned int *rank, double epsilon, double **T) |
| Computes the kernel and the rank of a matrix. More...
|
|
unsigned int | FindKernelAndIndependentRows (unsigned int nr, unsigned int nc, double *mT, unsigned int dof, double epsilon, boolean *singular, boolean **IR, double **T) |
| Computes the kernel of a matrix and determines the independent rows of this matrix. More...
|
|
unsigned int | FindIndependentRows (unsigned int nr, unsigned int nc, double *mT, unsigned int dof, double epsilon, boolean *singular, boolean **IR) |
| Computes the independent rows of this matrix. More...
|
|
double * | GetNewtonMatrixBuffer (TNewton *n) |
| Buffer to store the Newton matrix. More...
|
|
double * | GetNewtonRHBuffer (TNewton *n) |
| Buffer to store the Newton right hand. More...
|
|
void | NewtonSetMatrix (unsigned int i, unsigned int j, double v, TNewton *n) |
| Defines the matrix being used in a Newton step. More...
|
|
void | NewtonSetRH (unsigned int i, double v, TNewton *n) |
| Defines the vector being used in a Newton step. More...
|
|
double * | GetLSMatrixBuffer (TLinearSystem *ls) |
| Buffer to store the A matrix. More...
|
|
double * | GetLSRHBuffer (TLinearSystem *ls) |
| Buffer to store the linear system right hand (RH). More...
|
|
double * | GetLSSolutionBuffer (TLinearSystem *ls) |
| Buffer to store the linear system solution. More...
|
|
void | LSSetMatrix (unsigned int i, unsigned int j, double v, TLinearSystem *ls) |
| Defines the matrix being used in a linear system. More...
|
|
void | LSSetRH (unsigned int i, unsigned int j, double v, TLinearSystem *ls) |
| Defines the vector being used in a linear system. More...
|
|
void | InitBroyden (unsigned int nr, unsigned int nc, TBroyden *b) |
| Defines a Broyden structure. More...
|
|
void | ResetBroyden (TBroyden *b) |
| Resets a Broyden structure. More...
|
|
double * | GetBroydenMatrixBuffer (TBroyden *b) |
| Buffer to store the Broyden matrix. More...
|
|
double * | GetBroydenRHBuffer (TBroyden *b) |
| Buffer to store the Broyden right hand. More...
|
|
int | BroydenStep (double epsilon, unsigned int *tp, double *x, TBroyden *b) |
| One step in a Broyden iteration. More...
|
|
void | BroydenUpdateJacobian (double epsilon, TBroyden *b) |
| Updates the approximated Jacobian. More...
|
|
void | DeleteBroyden (TBroyden *b) |
| Releases a Broyden structure. More...
|
|
◆ dgpadm_()
void dgpadm_ |
( |
int * |
ideg, |
|
|
int * |
m, |
|
|
double * |
t, |
|
|
double * |
H, |
|
|
int * |
ldh, |
|
|
double * |
wsp, |
|
|
int * |
lwsp, |
|
|
int * |
ipiv, |
|
|
int * |
iexph, |
|
|
int * |
ns, |
|
|
int * |
iflag |
|
) |
| |
Description taken directly from EXPOKIT:
Computes exp(t*H), the matrix exponential of a general matrix in full, using the irreducible rational Pade approximation to the exponential function exp(x) = r(x) = (+/-)( I + 2*(q(x)/p(x)) ), combined with scaling-and-squaring.
- Parameters
-
ideg | The degre of the diagonal Pade to be used. A value of 6 is generally satisfactory. |
m | The size of the matrix H. |
t | The time. |
H | The matrix (in column major order!) |
ldh | The leading size of H (distance between elements in the same row). |
wsp | Workspace. This must be of size at least 4*m*m+ideg+1. |
lwsp | The size of wsp. |
ipiv | Workspace of size m. (index permutation?). |
iexph | Index in wsp where the result is stored, i.e., exp(tH) is located at wsp(iexph ... iexph+m*m-1) Note that fortran buffers are indexed from 1 (have to substract one). |
ns | Number of scaling-squaring used. |
iflag | Information flag (0 = no problem). |
Referenced by MatrixExponential().
◆ dgchbv_()
void dgchbv_ |
( |
int * |
m, |
|
|
double * |
t, |
|
|
double * |
H, |
|
|
int * |
ldh, |
|
|
double * |
y, |
|
|
double * |
wsp, |
|
|
int * |
iwsp, |
|
|
int * |
iflag |
|
) |
| |
Description taken directly from EXPOKIT:
dgchbv computes y = exp(t*H)*y using the partial fraction expansion of the uniform rational Chebyshev approximation to exp(-x) of type (14,14). H is a General matrix. About 14-digit accuracy is expected if the matrix H is negative definite. The algorithm may behave poorly otherwise.
- Parameters
-
m | The size of the matrix H. |
t | The time. |
H | The matrix (in column major order!) |
ldh | The leading size of H (distance between elements in the same row). |
y | The vector. Overwritten in the output. |
wsp | Workspace. This is of size 2*m*(m+2). The '2*' is because it strores complex numbers and not just doubles. |
iwsp | Workspace of size m. (index permutation?) |
iflag | Information flag (0 = no problem). |
Referenced by MatrixExponentialVector().
◆ AnalyzeKernel()
unsigned int AnalyzeKernel |
( |
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
mT, |
|
|
unsigned int |
dof, |
|
|
double |
epsilon, |
|
|
boolean |
computeRank, |
|
|
boolean |
checkRank, |
|
|
boolean |
getT, |
|
|
boolean |
getBase, |
|
|
boolean * |
singular, |
|
|
unsigned int * |
rank, |
|
|
boolean ** |
IR, |
|
|
double ** |
T |
|
) |
| |
Analyzes the kernel of a matrix and extract different information, as requested by the caller. This function for many purposes
- Parameters
-
nr | Number of rows of the matrix (no transposed). |
nc | Number of columns of the matrix (no transposed). |
mT | The TRANSPOSED matrix stored as a vector. |
dof | Expected dimension of the kernel. Can be zero if the getRank is TRUE. |
epsilon | Values below epsilon are taken as zero. |
computeRank | If the rank has to be computed from the kernel analysis. Otherwise the information provided by 'dof' is taken as good. |
checkRank | TRUE if an error has to be triggered if the rank is different from the expected one. This only has effect if computeRank is FALSE. |
getT | TRUE if we have to return a basis of the kernel. |
getBase | TRUE if we have to return a basis of the input matrix (selected rows). |
singular | TRUE if the matrix is singular (has more null eigen values than the expected ones). Output. |
rank | Rank of the input matrix. Computed in this function if computeRank is TRUE. Otherwise it is just deduced relying on the 'dof' parameter. |
IR | The set of independent rows as a boolean vector with as many entries as rows in the input matrix and TRUE for the independent rows. The space for this vector is allocated here but must be deallocated externally. If the matrix is singular this contains the most likely basis of the matrix (up to the numerical accuracy). Caution must be taken to use this output in this case. This is only allocated if getBase is TRUE. |
T | The output kernel. This is a (nc x dof) matrix (stored as a vector). The space for this matrix is allocated in this function but must be de-allocated externally. Only allocated if getT is TRUE. |
- Returns
- 0 if all the operations are correct, 1 if the kernel is larger than expected, 2 if it is smaller than expected, 3 if there is an error in the QR decomposition.
Referenced by FindIndependentRows(), FindKernel(), FindKernelAndIndependentRows(), FindKernelAndRank(), and FindRank().
◆ BasisColumnSpace()
unsigned int BasisColumnSpace |
( |
double |
epsilon, |
|
|
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
A, |
|
|
boolean ** |
IC, |
|
|
unsigned int * |
nrT, |
|
|
unsigned int * |
ncT, |
|
|
double ** |
T |
|
) |
| |
Orthonormalizes the columns of a matrix via QR decomposition.
This is an alternative to OrthonormalizeColumns which is probably faster and numerically more stable.
- Parameters
-
epsilon | Numerical accuracy. Values below epsilon are taken as zero. |
nr | Number of rows of the matrix to analyze. |
nc | Number of columns of the matrxi to analyze. |
A | The matrix whose columns space is analyzed. Please, note that this matrix is modified inside the function and that at the output it can not be easily interpreted any more. |
IC | Boolean array identifying the independent columns in A. Allocated internally. |
nrT | Number of rows of the orthonormal basis. It is nr in all the cases. |
ncT | Number of columsn of the orthonormal basis. It can be lower than nc. This is the number of TRUE values in IC. |
T | The computed orthonormal basis allocated internally and stored as a vector. |
◆ MatrixExponential()
unsigned int MatrixExponential |
( |
unsigned int |
m, |
|
|
double * |
A, |
|
|
double |
t, |
|
|
double * |
eA |
|
) |
| |
Computes the exponential of a matrix (at a given time step).
For a constant matrix A, the error of this function grows with the time 't'. Use with caution for large t's (i.e., 5, 10, is already large).
- Parameters
-
m | Number of rows/columsn of the matrix. |
A | The matrix. |
t | The time. |
eA | The exponential of the matrix: exp(A*t). |
- Returns
- 1 if the process was succesful.
Definition at line 1583 of file algebra.c.
References dgpadm_(), IdentityMatrix(), and NEW.
Referenced by dGt(), and LQRPolicy().
◆ MatrixExponentialVector()
unsigned int MatrixExponentialVector |
( |
unsigned int |
m, |
|
|
double * |
A, |
|
|
double |
t, |
|
|
double * |
v |
|
) |
| |
Computes the product of the exponential of a matrix with a vector at a given time step.
For a constant matrix A, the error of this function grows with the time 't'. Use with caution for large t's because the error in this case scales much faster than when using MatrixExponential. In this case the error can be significative just for t=1. If accuracy is required (over speed) just use MatrixExponential and multiply the result by vector v.
- Parameters
-
m | Number of rows/columsn of the matrix. |
A | The matrix. |
t | The time. |
v | The input vector. Re-used to store the output: exp(A*t)*v. |
- Returns
- 1 if the process was succesful.
Definition at line 1660 of file algebra.c.
References dgchbv_(), and NEW.
◆ FindRank()
unsigned int FindRank |
( |
double |
epsilon, |
|
|
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
mT |
|
) |
| |
Determines the rank of a matrix, i.e. the dimension of the space spanned by the rows/column of the matrix.
For a given problem, the number of variables minus the rank of the Jacobian gives the dimensionality of the solution space, assuming that the Jacobian is evaluated in a regular point. The dimensionality of the solution space is the same as that of its tangent space.
IMPORTANT: The intput matrix may be modified inside this function. The caller CAN NOT assume that it remains unchanged.
- Parameters
-
epsilon | Numerical accuracy. |
nr | Number of rows of the matrix (no transposed). |
nc | Number of columns of the matrix (no transposed). |
mT | The TRANSPOSED matrix. This matrix is overwriten inside this function!! |
- Returns
- The rank of the matrix. NO_UINT if the rank can not be computed.
Definition at line 1724 of file algebra.c.
References AnalyzeKernel(), FALSE, and TRUE.
Referenced by ManifoldDimension().
◆ FindKernel()
unsigned int FindKernel |
( |
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
mT, |
|
|
unsigned int |
dof, |
|
|
boolean |
check, |
|
|
double |
epsilon, |
|
|
double ** |
T |
|
) |
| |
Defines a basis of the null space of a matrix.
IMPORTANT: The intput matrix may be modified inside this function. The caller CAN NOT assume that it remains unchanged.
- Parameters
-
nr | Number of rows of the matrix (no transposed). |
nc | Numbe of columns of the matrix (no transposed). |
mT | The TRANSPOSED matrix stored as a vector. This matrix is overwriten inside this function!! |
dof | Expected dimension of the kernel. |
check | If TRUE the function introduces some consistancy checks (whether the kernel dimensionality is larger or smaller than the expected one). |
epsilon | Values below epsilon are taken as zero. |
T | The output kernel. This is a (nc x dof) matrix (stored as a vector). The space for this matrix is allocated in this function but must be de-allocated externally. |
- Returns
- 0 if all the operations are correct, 1 if there the kernel is too large (i.e., the given point is singular), 2 if the chart could not be defined since the kernel is too small at the given point, and 3 if there is an error while performing QR decomposition. These outputs come directly from AnalyzeKernel.
Definition at line 1738 of file algebra.c.
References AnalyzeKernel(), FALSE, and TRUE.
Referenced by FindRightNullVector(), and RefineSingularPoint().
◆ FindKernelAndRank()
unsigned int FindKernelAndRank |
( |
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
mT, |
|
|
unsigned int * |
rank, |
|
|
double |
epsilon, |
|
|
double ** |
T |
|
) |
| |
Defines a basis of the null space of a matrix, without any clue about the dimension of this null space. The returne rank is actually the number of variables minus the dimension of such space.
IMPORTANT: The intput matrix may be modified inside this function. The caller CAN NOT assume that it remains unchanged.
- Parameters
-
nr | Number of rows of the matrix (no transposed). |
nc | Numbe of columns of the matrix (no transposed). |
mT | The TRANSPOSED matrix stored as a vector. This matrix is overwriten inside this function!! |
rank | The rank of the input matrix: number of variables (columns) minus dimension of the output null space (num columns of the matrix representing the matrix). |
epsilon | Values below epsilon are taken as zero. |
T | The output kernel. This is a (nc x dof) matrix (stored as a vector). The space for this matrix is allocated in this function but must be de-allocated externally. |
- Returns
- 0 if all the operations are correct and 3 if there is an error while performing QR decomposition. These outputs come directly from AnalyzeKernel.
Definition at line 1752 of file algebra.c.
References AnalyzeKernel(), FALSE, and TRUE.
◆ FindKernelAndIndependentRows()
unsigned int FindKernelAndIndependentRows |
( |
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
mT, |
|
|
unsigned int |
dof, |
|
|
double |
epsilon, |
|
|
boolean * |
singular, |
|
|
boolean ** |
IR, |
|
|
double ** |
T |
|
) |
| |
Defines a basis of the null space of a matrix and determines a subset of the rows of the matrix that are independent.
This is useful because in our case most (all?) the matrices have redundancy (i.e., rows that are linearly dependent on other rows). However, for some purposes we need to determine a subset of the rows that are linearly independent.
IMPORTANT: The intput matrix may be modified inside this function. The caller CAN NOT assume that it remains unchanged.
- Parameters
-
nr | Number of rows of the matrix (no transposed). |
nc | Numbe of columns of the matrix (no transposed). |
mT | The TRANSPOSED matrix stored as a vector. This matrix is overwriten inside this function!! |
dof | Expected dimension of the kernel. If zero, the function tries to determine the rank automatically. |
epsilon | Values below epsilon are taken as zero. |
singular | TRUE if the matrix is singular (has more null eigen values than the expected ones). Output. |
IR | The set of independent rows as a boolean vector with as many entriees as rows in the input matrix and TRUE for the independent rows. The space for this vector is allocated here but must be deallocated externally. If the matrix is singular this contains the most likely basis of the matrix (up to the numerical accuracy). Caution must be taken to use this output in this case. |
T | The output kernel. This is a (nc x dof) matrix (stored as a vector). The space for this matrix is allocated in this function but must be de-allocated externally. |
- Returns
- 0 if all the operations are correct, 1 if there the kernel is too large (i.e., the given point is singular), 2 if the chart could not be defined since the kernel is too small at the given point, and 3 if there is an error while performing QR decomposition. These outputs come directly from AnalyzeKernel.
Definition at line 1765 of file algebra.c.
References AnalyzeKernel(), FALSE, and TRUE.
Referenced by ComputeJacobianKernelBasis().
◆ FindIndependentRows()
unsigned int FindIndependentRows |
( |
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
double * |
mT, |
|
|
unsigned int |
dof, |
|
|
double |
epsilon, |
|
|
boolean * |
singular, |
|
|
boolean ** |
IR |
|
) |
| |
Determines a subset of the rows of a matrix that are independent.
IMPORTANT: The intput matrix may be modified inside this function. The caller CAN NOT assume that it remains unchanged.
- Parameters
-
nr | Number of rows of the matrix (no transposed). |
nc | Numbe of columns of the matrix (no transposed). |
mT | The TRANSPOSED matrix stored as a vector. This matrix is overwriten inside this function!! |
dof | Expected dimension of the kernel. If zero, the function tries to determine the rank automatically. |
epsilon | Values below epsilon are taken as zero. |
singular | TRUE if the matrix is singular (has more null eigen values than the expected ones). Output. |
IR | The set of independent rows as a boolean vector with as many entriees as rows in the input matrix and TRUE for the independent rows. The space for this vector is allocated here but must be deallocated externally. If the matrix is singular this contains the most likely basis of the matrix (up to the numerical accuracy). Caution must be taken to use this output in this case. |
- Returns
- 0 if all the operations are correct, 1 if there the kernel is too large (i.e., the given point is singular), 2 if the chart could not be defined since the kernel is too small at the given point, and 3 if there is an error while performing QR decomposition. These outputs come directly from AnalyzeKernel.
Definition at line 1779 of file algebra.c.
References AnalyzeKernel(), FALSE, and TRUE.
Referenced by GetPositionJacobian().
◆ GetNewtonMatrixBuffer()
double* GetNewtonMatrixBuffer |
( |
TNewton * |
n | ) |
|
|
inline |
◆ GetNewtonRHBuffer()
double* GetNewtonRHBuffer |
( |
TNewton * |
n | ) |
|
|
inline |
◆ NewtonSetMatrix()
void NewtonSetMatrix |
( |
unsigned int |
i, |
|
|
unsigned int |
j, |
|
|
double |
v, |
|
|
TNewton * |
n |
|
) |
| |
|
inline |
Sets one element of the matrix to be used in one Newton step. This matrix is typically initilized externally, but here we provide a mehtod to set it.
- Parameters
-
i | The row. |
j | The column. |
v | The new value. |
n | The Newton structure to set. |
Definition at line 1804 of file algebra.c.
References RC2INDEX.
Referenced by CuikNewtonInBox().
◆ NewtonSetRH()
void NewtonSetRH |
( |
unsigned int |
i, |
|
|
double |
v, |
|
|
TNewton * |
n |
|
) |
| |
|
inline |
Sets one element of the vector to be used in one Newton step. This vector is typically initilized externally, but here we provide a mehtod to set it.
- Parameters
-
i | The index in the vector. |
v | The new value. |
n | The Newton structure to set. |
Definition at line 1809 of file algebra.c.
Referenced by CuikNewtonInBox().
◆ GetLSMatrixBuffer()
double* GetLSMatrixBuffer |
( |
TLinearSystem * |
ls | ) |
|
|
inline |
◆ GetLSRHBuffer()
double* GetLSRHBuffer |
( |
TLinearSystem * |
ls | ) |
|
|
inline |
◆ GetLSSolutionBuffer()
double* GetLSSolutionBuffer |
( |
TLinearSystem * |
ls | ) |
|
|
inline |
◆ LSSetMatrix()
void LSSetMatrix |
( |
unsigned int |
i, |
|
|
unsigned int |
j, |
|
|
double |
v, |
|
|
TLinearSystem * |
ls |
|
) |
| |
|
inline |
Sets one element of the matrix to be used in a linear system.
- Parameters
-
i | The row. |
j | The column. |
v | The new value. |
ls | The linear system structure. |
Definition at line 1829 of file algebra.c.
References RC2INDEX.
◆ LSSetRH()
void LSSetRH |
( |
unsigned int |
i, |
|
|
unsigned int |
j, |
|
|
double |
v, |
|
|
TLinearSystem * |
ls |
|
) |
| |
|
inline |
Sets one element of the right hand vector of a linear system.
- Parameters
-
i | The row index. |
j | The column index. |
v | The new value. |
ls | The linear system structure. |
Definition at line 1834 of file algebra.c.
References RC2INDEX.
◆ InitBroyden()
void InitBroyden |
( |
unsigned int |
nr, |
|
|
unsigned int |
nc, |
|
|
TBroyden * |
b |
|
) |
| |
◆ ResetBroyden()
Resets the Broyden structure (the internal iteration) without releasing memory.
- Parameters
-
b | The Broyden structure to initialize. |
Definition at line 1862 of file algebra.c.
Referenced by NextDynamicState().
◆ GetBroydenMatrixBuffer()
double* GetBroydenMatrixBuffer |
( |
TBroyden * |
b | ) |
|
|
inline |
Buffer to store the Broyden matrix.
This buffer must be accessed using the RC2INDEX macro since the matrix can be stored row major or column major depending on the underlying lineal algebra library being used.
- Parameters
-
- Returns
- A pointer to the buffer where to store the matrix.
Definition at line 1867 of file algebra.c.
References Error().
Referenced by InitDynamicSpace().
◆ GetBroydenRHBuffer()
double* GetBroydenRHBuffer |
( |
TBroyden * |
b | ) |
|
|
inline |
Buffer to store the Broyden RH.
- Parameters
-
- Returns
- A pointer to the buffer where to store the RH.
Definition at line 1874 of file algebra.c.
Referenced by InitDynamicSpace().
◆ BroydenStep()
int BroydenStep |
( |
double |
epsilon, |
|
|
unsigned int * |
tp, |
|
|
double * |
x, |
|
|
TBroyden * |
b |
|
) |
| |
Computes and applies the correction in one step of a Broyden iteration. Calls the Newton iteration with the current Jacobian matrix and updates this matrix.
The result is stored in the right-hand side buffer.
- Parameters
-
epsilon | Numerical accuracy. |
tp | Topology of the variables. |
x | The current approximation to the system solution. This is updated internally. |
b | The Broyden information. |
- Returns
- 0 if no error is encountered in the Broyden step.
Definition at line 1879 of file algebra.c.
References LSSolve(), NEW, Norm(), PI2PI, PrintMatrix(), PrintVector(), and TOPOLOGY_S.
Referenced by NextDynamicState().
◆ BroydenUpdateJacobian()
void BroydenUpdateJacobian |
( |
double |
epsilon, |
|
|
TBroyden * |
b |
|
) |
| |
◆ DeleteBroyden()
|
Follow us!