LAPACK++  2022.05.00
LAPACK C++ API
Constrained

Functions

int64_t lapack::ggglm (int64_t n, int64_t m, int64_t p, double *A, int64_t lda, double *B, int64_t ldb, double *D, double *X, double *Y)
 
int64_t lapack::ggglm (int64_t n, int64_t m, int64_t p, float *A, int64_t lda, float *B, int64_t ldb, float *D, float *X, float *Y)
 
int64_t lapack::ggglm (int64_t n, int64_t m, int64_t p, std::complex< double > *A, int64_t lda, std::complex< double > *B, int64_t ldb, std::complex< double > *D, std::complex< double > *X, std::complex< double > *Y)
 Solves a general Gauss-Markov linear model (GLM) problem: More...
 
int64_t lapack::ggglm (int64_t n, int64_t m, int64_t p, std::complex< float > *A, int64_t lda, std::complex< float > *B, int64_t ldb, std::complex< float > *D, std::complex< float > *X, std::complex< float > *Y)
 
int64_t lapack::gglse (int64_t m, int64_t n, int64_t p, double *A, int64_t lda, double *B, int64_t ldb, double *C, double *D, double *X)
 
int64_t lapack::gglse (int64_t m, int64_t n, int64_t p, float *A, int64_t lda, float *B, int64_t ldb, float *C, float *D, float *X)
 
int64_t lapack::gglse (int64_t m, int64_t n, int64_t p, std::complex< double > *A, int64_t lda, std::complex< double > *B, int64_t ldb, std::complex< double > *C, std::complex< double > *D, std::complex< double > *X)
 Solves the linear equality-constrained least squares (LSE) problem: More...
 
int64_t lapack::gglse (int64_t m, int64_t n, int64_t p, std::complex< float > *A, int64_t lda, std::complex< float > *B, int64_t ldb, std::complex< float > *C, std::complex< float > *D, std::complex< float > *X)
 

Detailed Description

Function Documentation

◆ ggglm()

int64_t lapack::ggglm ( int64_t  n,
int64_t  m,
int64_t  p,
std::complex< double > *  A,
int64_t  lda,
std::complex< double > *  B,
int64_t  ldb,
std::complex< double > *  D,
std::complex< double > *  X,
std::complex< double > *  Y 
)

Solves a general Gauss-Markov linear model (GLM) problem:

\[ \min_x || y ||_2 \text{ subject to } d = A x + B y, \]

where A is an n-by-m matrix, B is an n-by-p matrix, and d is a given n-vector. It is assumed that m <= n <= m+p, and

\[ rank(A) = m \text{ and } rank( [A, B] ) = n. \]

Under these assumptions, the constrained equation is always consistent, and there is a unique solution x and a minimal 2-norm solution y, which is obtained using a generalized QR factorization of the matrices (A, B) given by

\[ A = Q \begin{bmatrix} R \\ 0 \end{bmatrix}, \quad B = Q T Z. \]

In particular, if matrix B is square nonsingular, then the problem GLM is equivalent to the following weighted linear least squares problem

\[ \min_x || B^{-1} (d - A x) ||_2. \]

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>.

Parameters
[in]nThe number of rows of the matrices A and B. n >= 0.
[in]mThe number of columns of the matrix A. 0 <= m <= n.
[in]pThe number of columns of the matrix B. p >= n-m.
[in,out]AThe n-by-m matrix A, stored in an lda-by-m array. On entry, the n-by-m matrix A. On exit, the upper triangular part of the array A contains the m-by-m upper triangular matrix R.
[in]ldaThe leading dimension of the array A. lda >= max(1,n).
[in,out]BThe n-by-p matrix B, stored in an ldb-by-p array. On entry, the n-by-p matrix B. On exit, if n <= p, the upper triangle of the subarray B(1:n,p-n+1:p) contains the n-by-n upper triangular matrix T; if n > p, the elements on and above the (n-p)th subdiagonal contain the n-by-p upper trapezoidal matrix T.
[in]ldbThe leading dimension of the array B. ldb >= max(1,n).
[in,out]DThe vector D of length n. On entry, D is the left hand side of the GLM equation. On exit, D is destroyed.
[out]XThe vector X of length m.
[out]YThe vector Y of length p. On exit, X and Y are the solutions of the GLM problem.
Returns
= 0: successful exit.
= 1: the upper triangular factor R associated with A in the generalized QR factorization of the pair (A, B) is singular, so that rank(A) < m; the least squares solution could not be computed.
= 2: the bottom (n-m) by (n-m) part of the upper trapezoidal factor T associated with B in the generalized QR factorization of the pair (A, B) is singular, so that \(rank( [A, B] ) < n\); the least squares solution could not be computed.

◆ gglse()

int64_t lapack::gglse ( int64_t  m,
int64_t  n,
int64_t  p,
std::complex< double > *  A,
int64_t  lda,
std::complex< double > *  B,
int64_t  ldb,
std::complex< double > *  C,
std::complex< double > *  D,
std::complex< double > *  X 
)

Solves the linear equality-constrained least squares (LSE) problem:

\[ \min_x || c - A x ||_2 \text{ subject to } B x = d \]

where A is an m-by-n matrix, B is a p-by-n matrix, c is a given m-vector, and d is a given p-vector. It is assumed that p <= n <= m+p, and

\[ rank(B) = p \]

and

\[ rank\left( \begin{bmatrix} A \\ B \end{bmatrix} \right) = n. \]

These conditions ensure that the LSE problem has a unique solution, which is obtained using a generalized RQ factorization of the matrices (B, A) given by

\[ B = \begin{bmatrix} 0 & R \end{bmatrix} Q, \quad A = Z T Q. \]

Overloaded versions are available for float, double, std::complex<float>, and std::complex<double>.

Parameters
[in]mThe number of rows of the matrix A. m >= 0.
[in]nThe number of columns of the matrices A and B. n >= 0.
[in]pThe number of rows of the matrix B. 0 <= p <= n <= m+p.
[in,out]AThe m-by-n matrix A, stored in an lda-by-n array. On entry, the m-by-n matrix A. On exit, the elements on and above the diagonal of the array contain the min(m,n)-by-n upper trapezoidal matrix T.
[in]ldaThe leading dimension of the array A. lda >= max(1,m).
[in,out]BThe p-by-n matrix B, stored in an ldb-by-n array. On entry, the p-by-n matrix B. On exit, the upper triangle of the subarray B(1:p,n-p+1:n) contains the p-by-p upper triangular matrix R.
[in]ldbThe leading dimension of the array B. ldb >= max(1,p).
[in,out]CThe vector C of length m. On entry, C contains the right hand side vector for the least squares part of the LSE problem. On exit, the residual sum of squares for the solution is given by the sum of squares of elements n-p+1 to m of vector C.
[in,out]DThe vector D of length p. On entry, D contains the right hand side vector for the constrained equation. On exit, D is destroyed.
[out]XThe vector X of length n. On exit, X is the solution of the LSE problem.
Returns
= 0: successful exit.
= 1: the upper triangular factor R associated with B in the generalized RQ factorization of the pair (B, A) is singular, so that rank(B) < p; the least squares solution could not be computed.
= 2: the (n-p) by (n-p) part of the upper trapezoidal factor T associated with A in the generalized RQ factorization of the pair (B, A) is singular, so that

\[ rank\left( \begin{bmatrix} A \\ B \end{bmatrix} \right) < n; \]

the least squares solution could not be computed.