zgetsqrhrt()

subroutine zgetsqrhrt	(	integer	m,
		integer	n,
		integer	mb1,
		integer	nb1,
		integer	nb2,
		complex*16, dimension( lda, * )	a,
		integer	lda,
		complex*16, dimension( ldt, * )	t,
		integer	ldt,
		complex*16, dimension( * )	work,
		integer	lwork,
		integer	info
	)

ZGETSQRHRT

Download ZGETSQRHRT + dependencies [TGZ] [ZIP] [TXT]

Purpose:

 ZGETSQRHRT computes a NB2-sized column blocked QR-factorization
 of a complex M-by-N matrix A with M >= N,

    A = Q * R.

 The routine uses internally a NB1-sized column blocked and MB1-sized
 row blocked TSQR-factorization and perfors the reconstruction
 of the Householder vectors from the TSQR output. The routine also
 converts the R_tsqr factor from the TSQR-factorization output into
 the R factor that corresponds to the Householder QR-factorization,

    A = Q_tsqr * R_tsqr = Q * R.

 The output Q and R factors are stored in the same format as in ZGEQRT
 (Q is in blocked compact WY-representation). See the documentation
 of ZGEQRT for more details on the format.

Parameters

[in]	M	M is INTEGER The number of rows of the matrix A. M >= 0.
[in]	N	N is INTEGER The number of columns of the matrix A. M >= N >= 0.
[in]	MB1	MB1 is INTEGER The row block size to be used in the blocked TSQR. MB1 > N.
[in]	NB1	NB1 is INTEGER The column block size to be used in the blocked TSQR. N >= NB1 >= 1.
[in]	NB2	NB2 is INTEGER The block size to be used in the blocked QR that is output. NB2 >= 1.
[in,out]	A	A is COMPLEX*16 array, dimension (LDA,N) On entry: an M-by-N matrix A. On exit: a) the elements on and above the diagonal of the array contain the N-by-N upper-triangular matrix R corresponding to the Householder QR; b) the elements below the diagonal represent Q by the columns of blocked V (compact WY-representation).
[in]	LDA	LDA is INTEGER The leading dimension of the array A. LDA >= max(1,M).
[out]	T	T is COMPLEX*16 array, dimension (LDT,N)) The upper triangular block reflectors stored in compact form as a sequence of upper triangular blocks.
[in]	LDT	LDT is INTEGER The leading dimension of the array T. LDT >= NB2.
[out]	WORK	(workspace) COMPLEX*16 array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]	LWORK	The dimension of the array WORK. LWORK >= MAX( LWT + LW1, MAX( LWT+N*N+LW2, LWT+N*N+N ) ), where NUM_ALL_ROW_BLOCKS = CEIL((M-N)/(MB1-N)), NB1LOCAL = MIN(NB1,N). LWT = NUM_ALL_ROW_BLOCKS * N * NB1LOCAL, LW1 = NB1LOCAL * N, LW2 = NB1LOCAL * MAX( NB1LOCAL, ( N - NB1LOCAL ) ), If LWORK = -1, then 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 by XERBLA.
[out]	INFO	INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Contributors:

 November 2020, Igor Kozachenko,
                Computer Science Division,
                University of California, Berkeley

Definition at line 177 of file zgetsqrhrt.f.

      IMPLICIT NONE
*
*  -- LAPACK computational routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      INTEGER           INFO, LDA, LDT, LWORK, M, N, NB1, NB2, MB1
*     ..
*     .. Array Arguments ..
      COMPLEX*16        A( LDA, * ), T( LDT, * ), WORK( * )
*     ..
*
*  =====================================================================
*
*     .. Parameters ..
      COMPLEX*16         CONE
      parameter( cone = ( 1.0d+0, 0.0d+0 ) )
*     ..
*     .. Local Scalars ..
      LOGICAL            LQUERY
      INTEGER            I, IINFO, J, LW1, LW2, LWT, LDWT, LWORKOPT,
     $                   NB1LOCAL, NB2LOCAL, NUM_ALL_ROW_BLOCKS
*     ..
*     .. External Subroutines ..
      EXTERNAL           zcopy, zlatsqr, zungtsqr_row, zunhr_col,
     $                   xerbla
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          ceiling, dble, dcmplx, max, min
*     ..
*     .. Executable Statements ..
*
*     Test the input arguments
*
      info = 0
      lquery  = lwork.EQ.-1
      IF( m.LT.0 ) THEN
         info = -1
      ELSE IF( n.LT.0 .OR. m.LT.n ) THEN
         info = -2
      ELSE IF( mb1.LE.n ) THEN
         info = -3
      ELSE IF( nb1.LT.1 ) THEN
         info = -4
      ELSE IF( nb2.LT.1 ) THEN
         info = -5
      ELSE IF( lda.LT.max( 1, m ) ) THEN
         info = -7
      ELSE IF( ldt.LT.max( 1,  min( nb2, n ) ) ) THEN
         info = -9
      ELSE
*
*        Test the input LWORK for the dimension of the array WORK.
*        This workspace is used to store array:
*        a) Matrix T and WORK for ZLATSQR;
*        b) N-by-N upper-triangular factor R_tsqr;
*        c) Matrix T and array WORK for ZUNGTSQR_ROW;
*        d) Diagonal D for ZUNHR_COL.
*
         IF( lwork.LT.n*n+1 .AND. .NOT.lquery ) THEN
            info = -11
         ELSE
*
*           Set block size for column blocks
*
            nb1local = min( nb1, n )
*
            num_all_row_blocks = max( 1,
     $                   ceiling( dble( m - n ) / dble( mb1 - n ) ) )
*
*           Length and leading dimension of WORK array to place
*           T array in TSQR.
*
            lwt = num_all_row_blocks * n * nb1local
 
            ldwt = nb1local
*
*           Length of TSQR work array
*
            lw1 = nb1local * n
*
*           Length of ZUNGTSQR_ROW work array.
*
            lw2 = nb1local * max( nb1local, ( n - nb1local ) )
*
            lworkopt = max( lwt + lw1, max( lwt+n*n+lw2, lwt+n*n+n ) )
*
            IF( ( lwork.LT.max( 1, lworkopt ) ).AND.(.NOT.lquery) ) THEN
               info = -11
            END IF
*
         END IF
      END IF
*
*     Handle error in the input parameters and return workspace query.
*
      IF( info.NE.0 ) THEN
         CALL xerbla( 'ZGETSQRHRT', -info )
         RETURN
      ELSE IF ( lquery ) THEN
         work( 1 ) = dcmplx( lworkopt )
         RETURN
      END IF
*
*     Quick return if possible
*
      IF( min( m, n ).EQ.0 ) THEN
         work( 1 ) = dcmplx( lworkopt )
         RETURN
      END IF
*
      nb2local = min( nb2, n )
*
*
*     (1) Perform TSQR-factorization of the M-by-N matrix A.
*
      CALL zlatsqr( m, n, mb1, nb1local, a, lda, work, ldwt,
     $              work(lwt+1), lw1, iinfo )
*
*     (2) Copy the factor R_tsqr stored in the upper-triangular part
*         of A into the square matrix in the work array
*         WORK(LWT+1:LWT+N*N) column-by-column.
*
      DO j = 1, n
         CALL zcopy( j, a( 1, j ), 1, work( lwt + n*(j-1)+1 ), 1 )
      END DO
*
*     (3) Generate a M-by-N matrix Q with orthonormal columns from
*     the result stored below the diagonal in the array A in place.
*
 
      CALL zungtsqr_row( m, n, mb1, nb1local, a, lda, work, ldwt,
     $                   work( lwt+n*n+1 ), lw2, iinfo )
*
*     (4) Perform the reconstruction of Householder vectors from
*     the matrix Q (stored in A) in place.
*
      CALL zunhr_col( m, n, nb2local, a, lda, t, ldt,
     $                work( lwt+n*n+1 ), iinfo )
*
*     (5) Copy the factor R_tsqr stored in the square matrix in the
*     work array WORK(LWT+1:LWT+N*N) into the upper-triangular
*     part of A.
*
*     (6) Compute from R_tsqr the factor R_hr corresponding to
*     the reconstructed Householder vectors, i.e. R_hr = S * R_tsqr.
*     This multiplication by the sign matrix S on the left means
*     changing the sign of I-th row of the matrix R_tsqr according
*     to sign of the I-th diagonal element DIAG(I) of the matrix S.
*     DIAG is stored in WORK( LWT+N*N+1 ) from the ZUNHR_COL output.
*
*     (5) and (6) can be combined in a single loop, so the rows in A
*     are accessed only once.
*
      DO i = 1, n
         IF( work( lwt+n*n+i ).EQ.-cone ) THEN
            DO j = i, n
               a( i, j ) = -cone * work( lwt+n*(j-1)+i )
            END DO
         ELSE
            CALL zcopy( n-i+1, work(lwt+n*(i-1)+i), n, a( i, i ), lda )
         END IF
      END DO
*
      work( 1 ) = dcmplx( lworkopt )
      RETURN
*
*     End of ZGETSQRHRT
*

Here is the call graph for this function:

Here is the caller graph for this function: