?gels

Uses QR or LQ factorization to solve an overdetermined or underdetermined linear system with full rank matrix.

Interface Definition

C interface:

sgels_(const char *trans, const int *m, const int *n, const int *nrhs, float *a, const int *lda, float *b, const int *ldb, float *work, const int *lwork, int *info);

dgels_(const char *trans, const int *m, const int *n, const int *nrhs, double *a, const int *lda, double *b, const int *ldb, double *work, const int *lwork, int *info);

cgels_(const char *trans, const int *m, const int *n, const int *nrhs, float _Complex *a, const int *lda, float _Complex *b, const int *ldb, float _Complex *work, const int *lwork, int *info);

zgels_(const char *trans, const int *m, const int *n, const int *nrhs, double _Complex *a, const int *lda, double _Complex *b, const int *ldb, double _Complex *work, const int *lwork, int *info);

Fortran interface:

SGELS(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO);

DGELS(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO);

CGELS(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO);

ZGELS(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO);

Parameters

Parameters	Type	Description	Input/Output
trans	Character	'N': Matrix A is not transposed. 'C': Matrix A is transposed.	Input
m	Integer	Number of rows in matrix A. m ≥ 0.	Input
n	Integer	Number of columns in matrix A, n ≥ 0.	Input
nrhs	Integer	Number of right-hand side columns, nrhs ≥ 0.	Input
a	A single-precision floating-point array for sgels A double-precision floating-point array for dgels A complex single-precision array for cgels A complex double-precision array for zgels	Array with a size of ldan. On entry: an mn matrix. On exit: If m ≥ n, it stores the QR factorization result returned by GEQRF. If m < n, it stores the LQ factorization result returned by GELQF.	Input/Output
lda	Integer	Leading dimension of matrix A. lda ≥ max(1, m).	Input
b	A single-precision floating-point array for sgels A double-precision floating-point array for dgels A complex single-precision array for cgels A complex double-precision array for zgels	Right-hand side matrix, with a size of ldbnrhs. On entry: If trans* = 'N', it is an mnrhs matrix. If trans* = 'C' or 'T', it is an nnrhs matrix. On exit: When info=0: If trans='N' and m≥n, rows 1 to n of the matrix contain the least squares solution vectors; the residual sum of squares in each column is given by the sum of squares of modulus of elements n+1* to m in that column. If trans='N' and m<n, rows 1 to n of the matrix contain the minimum norm solution vectors. If trans='C' or 'T' and m≥n, rows 1 to n of the matrix contain the minimum norm solution vectors. If trans='C' or 'T' and m<n, rows 1 to n of the matrix contain the least squares solution vectors; the residual sum of squares in each column is given by the sum of squares of modulus of elements n+1 to m in that column.	Input/Output
ldb	Integer	Leading dimension of matrix b.	Input
work	A single-precision floating-point array for sgels A double-precision floating-point array for dgels A complex single-precision array for cgels A complex double-precision array for zgels	Work array. If info = 0, work(0) returns the optimal lwork size.	Output
lwork	Integer	Size of work.	Input
Info	Integer	Function execution status. 0: The execution is successful. Smaller than 0: If info = -i, the i-th parameter has an illegal value. Greater than 0: An algorithm error occurs.	Output

Dependency

#include "klapack.h"

Example

C interface:

const char trans = 'N';
const int m = 4;
const int n = 4;
const int nrhs = 2;
const int lda = 4;
const int ldb = 4;

double a[] = {72.1673, 66.1857, 64.7644, 28.0199, 91.4151,
              6.5180,  62.8483, 72.4323, 46.5760, 8.6928,
              28.9821, 42.1828, 18.6437, 99.8612, 35.6972,
              67.9812, 5.0880,  85.5035, 79.2945, 54.5920,
              28.6869, 49.7512, 7.5186,  28.6929, 84.6041};

double b[] = {9.4532, 1.5204, 2.2127, 0.9891, 7.1778,
              6.8955, 7.2465, 3.5019, 8.2268, 3.5287};

double qwork;
int lwork = -1;
int info = 0;

dgels_(&trans, &m, &n, &nrhs, a, &lda, b, &ldb, &qwork, &lwork, &info);
if (info != 0) {
    printf("Error, info = %d\n", info);
    return info;
}

lwork = (int)qwork;
double *work = (double*)malloc(lwork * sizeof(double));

dgels_(&trans, &m, &n, &nrhs, a, &lda, b, &ldb, work, &lwork, &info);
if (info != 0) {
    printf("Error, info = %d\n", info);
    return info;
}
/* output */
* a
* -120.698833     -108.770377     -57.958881      -100.842591
* 0.343169        75.924679       38.025488       -19.044303
* 0.335800        0.031671        7.685677        67.409117
* 0.145282        -0.313887       -0.556039       33.759113
* b
* 0.153366        0.083995
* -0.538653       0.096517
* 1.080726        -0.174382
* -0.145340       0.022261

Fortran interface:

CHARACTER :: tarns = "N"
PARAMETER (n = 4) 
PARAMETER (m = 4)
PARAMETER (lda = 4) 
PARAMETER (ldb = 4)  
PARAMETER (nrhs = 2) 

INTEGER :: info = 0 
REAL(8) :: a(lda, n) 
REAL(8) :: b(ldb, nrhs)
REAL(8), ALLOCATABLE :: work(:)
REAL(8) :: qwork(1)
INTEGER :: lwork = -1
  
DATA a / 72.1673, 66.1857, 64.7644, 28.0199, 91.4151,
&         6.5180,  62.8483, 72.4323, 46.5760, 8.6928,
&        28.9821, 42.1828, 18.6437, 99.8612, 35.6972,
&        67.9812, 5.0880,  85.5035, 79.2945, 54.5920,
&        28.6869, 49.7512, 7.5186,  28.6929, 84.6041 / 
DATA b / 9.4532, 1.5204, 2.2127, 0.9891, 7.1778,
&        6.8955, 7.2465, 3.5019, 8.2268, 3.5287 /
EXTERNAL DGELS
CALL DGELS(trans, m, n, nrhs, a, lda, b, ldb, qwork, lwork, info)
IF (info.NE.0) THEN 
    CALL EXIT(1) 
END IF
lwork = INT(qwork(1))  
ALLOCATE(work(lwork))
CALL DGELS(trans, m, n, nrhs, a, lda, b, ldb, work, lwork, info)
IF (info.NE.0) THEN 
    CALL EXIT(1) 
END IF
DEALLOCATE(work);
* 
* Output: 
* 
a
* -120.698833     -108.770377     -57.958881      -100.842591
* 0.343169        75.924679       38.025488       -19.044303
* 0.335800        0.031671        7.685677        67.409117
* 0.145282        -0.313887       -0.556039       33.759113
* b
* 0.153366        0.083995
* -0.538653       0.096517
* 1.080726        -0.174382
* -0.145340       0.022261

Parent topic: Least Squares Solution Functions