lin_svd_gen
Computes the SVD, A = USVT, of a real rectangular matrix A. An approximate generalized inverse and rank of A also can be computed.
Synopsis
#include <imsl.h>
float *imsl_f_lin_svd_gen (int m, int n, float a[], …, 0)
The type double function is imsl_d_lin_svd_gen.
Required Arguments
int m (Input)
Number of rows in the matrix.
Number of rows in the matrix.
int n (Input)
Number of columns in the matrix.
Number of columns in the matrix.
float a[] (Input)
Array of size m × n containing the matrix.
Array of size m × n containing the matrix.
Return Value
If no optional arguments are used, imsl_f_lin_svd_gen returns a pointer to an array of size min (m, n) containing the ordered singular values of the matrix. To release this space, use imsl_free. If no value can be computed, then NULL is returned.
Synopsis with Optional Arguments
#include <imsl.h>
float *imsl_f_lin_svd_gen (int m, int n, float a[],
IMSL_METHOD, int imeth
IMSL_A_COL_DIM, int a_col_dim,
IMSL_RETURN_USER, float s[],
IMSL_RANK, float tol, int *rank,
IMSL_U, float **p_u,
IMSL_U_USER, float u[],
IMSL_U_COL_DIM, int u_col_dim,
IMSL_V, float **p_v,
IMSL_V_USER, float v[],
IMSL_V_COL_DIM, int v_col_dim,
IMSL_INVERSE, float **p_gen_inva,
IMSL_INVERSE_USER, float gen_inva[],
IMSL_INV_COL_DIM, int gen_inva_col_dim,
0)
Optional Arguments
IMSL_METHOD, int imeth (Input)
The method used in the computation of the singular values and vectors.
The method used in the computation of the singular values and vectors.
Computational Method | ||
|---|---|---|
imeth | IMSL | LAPACK |
0 | Uses QR method to determine singular values and singular vectors. | Uses QR method if singular vectors are de-sired and the dqds algorithm (Fernando and Parlett, 1994) otherwise. |
1 | Uses the same algorithm as for imeth = 0. | Uses the dqds algorithm (Fernando and Parlett, 1994) if singular values only are desired and a divide-and-conquer algorithm if singular vectors are desired. |
Default: imeth = 0
NOTE: The LAPACK algorithms can be used if a vendor supplied library that supports LAPACK is available. Examples are Intel’s® Math Kernel Library (MKL) or Sun’s™ High Performance Library. Otherwise, only the native IMSL algorithm is available.
IMSL_A_COL_DIM, int a_col_dim (Input)
The column dimension of the array a.
Default: a_col_dim = n
The column dimension of the array a.
Default: a_col_dim = n
IMSL_RETURN_USER, float s[] (Output)
A user-allocated array of size min (m, n) containing the singular values of A in nonincreasing order. If IMSL_RETURN_USER is used, the return value of imsl_f_lin_svd_gen is s.
A user-allocated array of size min (m, n) containing the singular values of A in nonincreasing order. If IMSL_RETURN_USER is used, the return value of imsl_f_lin_svd_gen is s.
IMSL_RANK, float tol, int *rank (Input/Output)
float tol (Input)
Scalar containing the tolerance used to determine when a singular value is negligible and replaced by the value zero. If tol > 0, then a singular value si i is considered negligible if si.i ≤ tol. If tol < 0, then a singular value si.i is considered negligible if si.i ≤ ∣tol∣*∥A∥∞. In this case, ∣tol∣ should be an estimate of relative error or uncertainty in the data.
Default: tol = 100.0 * imsl_f_machine(4)
Scalar containing the tolerance used to determine when a singular value is negligible and replaced by the value zero. If tol > 0, then a singular value si i is considered negligible if si.i ≤ tol. If tol < 0, then a singular value si.i is considered negligible if si.i ≤ ∣tol∣*∥A∥∞. In this case, ∣tol∣ should be an estimate of relative error or uncertainty in the data.
Default: tol = 100.0 * imsl_f_machine(4)
int *rank (Output)
Integer containing an estimate of the rank of A.
Integer containing an estimate of the rank of A.
IMSL_U, float **p_u (Output)
The address of a pointer to an array of size m × min (m, n) containing the min (m, n) left singular vectors of A. On return, the necessary space is allocated by imsl_f_lin_svd_gen. Typically, float *p_u is declared, and &p_u is used as an argument.
The address of a pointer to an array of size m × min (m, n) containing the min (m, n) left singular vectors of A. On return, the necessary space is allocated by imsl_f_lin_svd_gen. Typically, float *p_u is declared, and &p_u is used as an argument.
IMSL_U_USER, float u[] (Output)
The address of a pointer to an array of size m × min (m, n) containing the min (m, n) left singular vectors of A. The left singular vectors can be returned using the storage locations of the array a. Note that the return of the left and right singular vectors cannot use the storage locations of a simultaneously.
The address of a pointer to an array of size m × min (m, n) containing the min (m, n) left singular vectors of A. The left singular vectors can be returned using the storage locations of the array a. Note that the return of the left and right singular vectors cannot use the storage locations of a simultaneously.
IMSL_U_COL_DIM, int u_col_dim (Input)
The column dimension of the array containing the left singular vectors.
Default: u_col_dim = min (m, n)
The column dimension of the array containing the left singular vectors.
Default: u_col_dim = min (m, n)
IMSL_V, float **p_v (Output)
The address of a pointer to array of size n × min (m, n) containing the right singular vectors of A. On return, the necessary space is allocated by imsl_f_lin_svd_gen. Typically, float *p_v is declared, and &p_v is used as an argument.
The address of a pointer to array of size n × min (m, n) containing the right singular vectors of A. On return, the necessary space is allocated by imsl_f_lin_svd_gen. Typically, float *p_v is declared, and &p_v is used as an argument.
IMSL_V_USER, float v[] (Output)
The address of a pointer to array of size n × min (m, n) containing the right singular vectors of A. The right singular vectors can be returned using the storage locations of the array a. Note that the return of the left and right singular vectors cannot use the storage locations of a simultaneously.
The address of a pointer to array of size n × min (m, n) containing the right singular vectors of A. The right singular vectors can be returned using the storage locations of the array a. Note that the return of the left and right singular vectors cannot use the storage locations of a simultaneously.
IMSL_V_COL_DIM, int v_col_dim (Input)
The column dimension of the array containing the right singular vectors.
Default: v_col_dim = min (m, n)
The column dimension of the array containing the right singular vectors.
Default: v_col_dim = min (m, n)
IMSL_INVERSE, float **p_gen_inva (Output)
The address of a pointer to an array of size n × m containing the generalized inverse of the matrix A. On return, the necessary space is allocated by imsl_f_lin_svd_gen. Typically, float *p_gen_inva is declared, and &p_gen_inva is used as an argument.
The address of a pointer to an array of size n × m containing the generalized inverse of the matrix A. On return, the necessary space is allocated by imsl_f_lin_svd_gen. Typically, float *p_gen_inva is declared, and &p_gen_inva is used as an argument.
IMSL_INVERSE_USER, float gen_inva[] (Output)
A user-allocated array of size n × m containing the general inverse of the matrix A.
A user-allocated array of size n × m containing the general inverse of the matrix A.
IMSL_INV_COL_DIM, int gen_inva_col_dim (Input)
The column dimension of the array containing the general inverse of the matrix A.
Default: gen_inva_col_dim = m
The column dimension of the array containing the general inverse of the matrix A.
Default: gen_inva_col_dim = m
Description
The function imsl_f_lin_svd_gen computes the singular value decomposition of a real matrix A. It first reduces the matrix A to a bidiagonal matrix B by pre- and post-multiplying Householder transformations. Then, the singular value decomposition of B is computed using the implicit-shifted QR algorithm. An estimate of the rank of the matrix A is obtained by finding the smallest integer k such that sk,k ≤ tol or sk,k ≤ ∣tol∣*∥A∥∞. Since si+1,i+1 ≤ si,i, it follows that all the si,i satisfy the same inequality for i = k, …, min (m, n) − 1. The rank is set to the value k − 1. If A = USVT, its generalized inverse is A+ = VS+ UT. Here,
Only singular values that are not negligible are reciprocated. If IMSL_INVERSE or IMSL_INVERSE_USER is specified, the function first computes the singular value decomposition of the matrix A. The generalized inverse is then computed. The function imsl_f_lin_svd_gen fails if the QR algorithm does not converge after 30 iterations isolating an individual singular value.
Examples
Example 1
This example computes the singular values of a real 6 × 4 matrix.
#include <imsl.h>
float a[] = {1.0, 2.0, 1.0, 4.0,
3.0, 2.0, 1.0, 3.0,
4.0, 3.0, 1.0, 4.0,
2.0, 1.0, 3.0, 1.0,
1.0, 5.0, 2.0, 2.0,
1.0, 2.0, 2.0, 3.0};
int main()
{
int m = 6, n = 4;
float *s;
/* Compute singular values */
s = imsl_f_lin_svd_gen (m, n, a, 0);
/* Print singular values */
imsl_f_write_matrix ("Singular values", 1, n, s, 0);
}
Output
Singular values
1 2 3 4
11.49 3.27 2.65 2.09
Example 2
This example computes the singular value decomposition of the 6 × 4 real matrix A. The singular values are returned in the user-provided array. The matrices U and V are returned in the space provided by the function imsl_f_lin_svd_gen.
#include <imsl.h>
float a[] = {1.0, 2.0, 1.0, 4.0,
3.0, 2.0, 1.0, 3.0,
4.0, 3.0, 1.0, 4.0,
2.0, 1.0, 3.0, 1.0,
1.0, 5.0, 2.0, 2.0,
1.0, 2.0, 2.0, 3.0};
int main()
{
int m = 6, n = 4;
float s[4], *p_u, *p_v;
/* Compute SVD */
imsl_f_lin_svd_gen (m, n, a,
IMSL_RETURN_USER, s,
IMSL_U, &p_u,
IMSL_V, &p_v,
0);
/* Print decomposition*/
imsl_f_write_matrix ("Singular values, S", 1, n, s, 0);
imsl_f_write_matrix ("Left singular vectors, U", m, n, p_u, 0);
imsl_f_write_matrix ("Right singular vectors, V", n, n, p_v, 0);
}
Output
Singular values, S
1 2 3 4
11.49 3.27 2.65 2.09
Left singular vectors, U
1 2 3 4
1 -0.3805 -0.1197 -0.4391 0.5654
2 -0.4038 -0.3451 0.0566 -0.2148
3 -0.5451 -0.4293 -0.0514 -0.4321
4 -0.2648 0.0683 0.8839 0.2153
5 -0.4463 0.8168 -0.1419 -0.3213
6 -0.3546 0.1021 0.0043 0.5458
Right singular vectors, V
1 2 3 4
1 -0.4443 -0.5555 0.4354 -0.5518
2 -0.5581 0.6543 -0.2775 -0.4283
3 -0.3244 0.3514 0.7321 0.4851
4 -0.6212 -0.3739 -0.4444 0.5261
Example 3
This example computes the rank and generalized inverse of a 3 × 2 matrix A. The rank and the 2 × 3 generalized inverse matrix A+ are printed.
#include <imsl.h>
#include <stdio.h>
float a[] =
{1.0, 0.0,
1.0, 1.0,
100.0, -50.0};
int main()
{
int m = 3, n = 2;
float tol;
float gen_inva[6];
float *s;
int rank;
/* Compute generalized inverse */
tol = 1.e-4;
s = imsl_f_lin_svd_gen (m, n, a,
IMSL_RANK, tol, &rank,
IMSL_INVERSE_USER, gen_inva,
IMSL_INV_COL_DIM, m,
0);
/* Print rank, singular values and */
/* generalized inverse. */
printf ("Rank of matrix = %2d", rank);
imsl_f_write_matrix ("Singular values", 1, n, s, 0);
imsl_f_write_matrix ("Generalized inverse", n, m, gen_inva,
IMSL_A_COL_DIM, m,
0);
}
Output
Rank of matrix = 2
Singular values
1 2
111.8 1.4
Generalized inverse
1 2 3
1 0.100 0.300 0.006
2 0.200 0.600 -0.008
Fatal Errors
IMSL_SLOWCONVERGENT_MATRIX | Convergence cannot be reached after 30 iterations. |
IMSL_UPDATE_PROCESS_FAILED | The algorithm failed to compute a singular value. The update process of divide-and-conquer failed. |
