The routine may be called by the names f08raf, nagf_lapackeig_dorcsd or its LAPACK name dorcsd.
3Description
The orthogonal matrix is partitioned as
where is a submatrix and the dimensions of the other submatrices , and are such that remains .
The CS decomposition of is
where , and are matrices, such that
is an orthogonal matrix containing the orthogonal matrix and the orthogonal matrix ;
is an orthogonal matrix containing the orthogonal matrix and the orthogonal matrix ; and
contains the non-negative diagonal submatrices and satisfying , where and the top left partition is .
The identity matrix is of order and vanishes if .
The identity matrix is of order and vanishes if .
The identity matrix is of order and vanishes if .
The identity matrix is of order and vanishes if .
In each of the four cases at least two of the identity matrices vanish.
The indicated zeros represent augmentations by additional rows or columns (but not both) to the square diagonal matrices formed by and or .
does not need to be stored in full; it is sufficient to return only the values for where and .
The algorithm used to perform the complete decomposition is described fully in Sutton (2009) including discussions of the stability and accuracy of the algorithm.
4References
Anderson E, Bai Z, Bischof C, Blackford S, Demmel J, Dongarra J J, Du Croz J J, Greenbaum A, Hammarling S, McKenney A and Sorensen D (1999) LAPACK Users' Guide (3rd Edition) SIAM, Philadelphia https://www.netlib.org/lapack/lug
Golub G H and Van Loan C F (2012) Matrix Computations (4th Edition) Johns Hopkins University Press, Baltimore
If on exit, contains the values that, together with , define the matrix in intermediate bidiagonal-block form remaining after nonconvergence. info specifies the number of nonzero PHI's.
28: – IntegerInput
On entry: the dimension of the array work as declared in the (sub)program from which f08raf is called.
If , a workspace query is assumed; the routine only calculates the optimal size of the work array, returns this value as the first entry of the work array, and no error message related to lwork is issued.
The minimum workspace required is where . The optimal workspace depends on internal block sizes and the relative dimensions of the problem.
Constraint:
or .
29: – Integer arrayWorkspace
30: – IntegerOutput
On exit: unless the routine detects an error (see Section 6).
6Error Indicators and Warnings
If , argument had an illegal value. An explanatory message is output, and execution of the program is terminated.
The Jacobi-type procedure failed to converge during an internal reduction to bidiagonal-block form. The process requires convergence to values, the value of info gives the number of converged values.
7Accuracy
The computed decomposition is nearly the exact decomposition for the nearby matrix , where
and is the machine precision.
8Parallelism and Performance
Background information to multithreading can be found in the Multithreading documentation.
f08raf is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.
f08raf makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the Users' Note for your implementation for any additional implementation-specific information.
9Further Comments
The total number of floating-point operations required to perform the full decomposition is approximately .
The decomposition is performed both on submatrices of the orthogonal matrix and on separated partition matrices. Code is also provided to perform a recombining check if required.