--- rpl/lapack/lapack/zungbr.f	2011/07/22 07:38:21	1.8
+++ rpl/lapack/lapack/zungbr.f	2011/11/21 20:43:23	1.9
@@ -1,9 +1,166 @@
+*> \brief \b ZUNGBR
+*
+*  =========== DOCUMENTATION ===========
+*
+* Online html documentation available at 
+*            http://www.netlib.org/lapack/explore-html/ 
+*
+*> \htmlonly
+*> Download ZUNGBR + dependencies 
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zungbr.f"> 
+*> [TGZ]</a> 
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zungbr.f"> 
+*> [ZIP]</a> 
+*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zungbr.f"> 
+*> [TXT]</a>
+*> \endhtmlonly 
+*
+*  Definition:
+*  ===========
+*
+*       SUBROUTINE ZUNGBR( VECT, M, N, K, A, LDA, TAU, WORK, LWORK, INFO )
+* 
+*       .. Scalar Arguments ..
+*       CHARACTER          VECT
+*       INTEGER            INFO, K, LDA, LWORK, M, N
+*       ..
+*       .. Array Arguments ..
+*       COMPLEX*16         A( LDA, * ), TAU( * ), WORK( * )
+*       ..
+*  
+*
+*> \par Purpose:
+*  =============
+*>
+*> \verbatim
+*>
+*> ZUNGBR generates one of the complex unitary matrices Q or P**H
+*> determined by ZGEBRD when reducing a complex matrix A to bidiagonal
+*> form: A = Q * B * P**H.  Q and P**H are defined as products of
+*> elementary reflectors H(i) or G(i) respectively.
+*>
+*> If VECT = 'Q', A is assumed to have been an M-by-K matrix, and Q
+*> is of order M:
+*> if m >= k, Q = H(1) H(2) . . . H(k) and ZUNGBR returns the first n
+*> columns of Q, where m >= n >= k;
+*> if m < k, Q = H(1) H(2) . . . H(m-1) and ZUNGBR returns Q as an
+*> M-by-M matrix.
+*>
+*> If VECT = 'P', A is assumed to have been a K-by-N matrix, and P**H
+*> is of order N:
+*> if k < n, P**H = G(k) . . . G(2) G(1) and ZUNGBR returns the first m
+*> rows of P**H, where n >= m >= k;
+*> if k >= n, P**H = G(n-1) . . . G(2) G(1) and ZUNGBR returns P**H as
+*> an N-by-N matrix.
+*> \endverbatim
+*
+*  Arguments:
+*  ==========
+*
+*> \param[in] VECT
+*> \verbatim
+*>          VECT is CHARACTER*1
+*>          Specifies whether the matrix Q or the matrix P**H is
+*>          required, as defined in the transformation applied by ZGEBRD:
+*>          = 'Q':  generate Q;
+*>          = 'P':  generate P**H.
+*> \endverbatim
+*>
+*> \param[in] M
+*> \verbatim
+*>          M is INTEGER
+*>          The number of rows of the matrix Q or P**H to be returned.
+*>          M >= 0.
+*> \endverbatim
+*>
+*> \param[in] N
+*> \verbatim
+*>          N is INTEGER
+*>          The number of columns of the matrix Q or P**H to be returned.
+*>          N >= 0.
+*>          If VECT = 'Q', M >= N >= min(M,K);
+*>          if VECT = 'P', N >= M >= min(N,K).
+*> \endverbatim
+*>
+*> \param[in] K
+*> \verbatim
+*>          K is INTEGER
+*>          If VECT = 'Q', the number of columns in the original M-by-K
+*>          matrix reduced by ZGEBRD.
+*>          If VECT = 'P', the number of rows in the original K-by-N
+*>          matrix reduced by ZGEBRD.
+*>          K >= 0.
+*> \endverbatim
+*>
+*> \param[in,out] A
+*> \verbatim
+*>          A is COMPLEX*16 array, dimension (LDA,N)
+*>          On entry, the vectors which define the elementary reflectors,
+*>          as returned by ZGEBRD.
+*>          On exit, the M-by-N matrix Q or P**H.
+*> \endverbatim
+*>
+*> \param[in] LDA
+*> \verbatim
+*>          LDA is INTEGER
+*>          The leading dimension of the array A. LDA >= M.
+*> \endverbatim
+*>
+*> \param[in] TAU
+*> \verbatim
+*>          TAU is COMPLEX*16 array, dimension
+*>                                (min(M,K)) if VECT = 'Q'
+*>                                (min(N,K)) if VECT = 'P'
+*>          TAU(i) must contain the scalar factor of the elementary
+*>          reflector H(i) or G(i), which determines Q or P**H, as
+*>          returned by ZGEBRD in its array argument TAUQ or TAUP.
+*> \endverbatim
+*>
+*> \param[out] WORK
+*> \verbatim
+*>          WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
+*>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
+*> \endverbatim
+*>
+*> \param[in] LWORK
+*> \verbatim
+*>          LWORK is INTEGER
+*>          The dimension of the array WORK. LWORK >= max(1,min(M,N)).
+*>          For optimum performance LWORK >= min(M,N)*NB, where NB
+*>          is the optimal blocksize.
+*>
+*>          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.
+*> \endverbatim
+*>
+*> \param[out] INFO
+*> \verbatim
+*>          INFO is INTEGER
+*>          = 0:  successful exit
+*>          < 0:  if INFO = -i, the i-th argument had an illegal value
+*> \endverbatim
+*
+*  Authors:
+*  ========
+*
+*> \author Univ. of Tennessee 
+*> \author Univ. of California Berkeley 
+*> \author Univ. of Colorado Denver 
+*> \author NAG Ltd. 
+*
+*> \date November 2011
+*
+*> \ingroup complex16GBcomputational
+*
+*  =====================================================================
       SUBROUTINE ZUNGBR( VECT, M, N, K, A, LDA, TAU, WORK, LWORK, INFO )
 *
-*  -- LAPACK routine (version 3.2) --
+*  -- LAPACK computational routine (version 3.4.0) --
 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
-*     November 2006
+*     November 2011
 *
 *     .. Scalar Arguments ..
       CHARACTER          VECT
@@ -13,86 +170,6 @@
       COMPLEX*16         A( LDA, * ), TAU( * ), WORK( * )
 *     ..
 *
-*  Purpose
-*  =======
-*
-*  ZUNGBR generates one of the complex unitary matrices Q or P**H
-*  determined by ZGEBRD when reducing a complex matrix A to bidiagonal
-*  form: A = Q * B * P**H.  Q and P**H are defined as products of
-*  elementary reflectors H(i) or G(i) respectively.
-*
-*  If VECT = 'Q', A is assumed to have been an M-by-K matrix, and Q
-*  is of order M:
-*  if m >= k, Q = H(1) H(2) . . . H(k) and ZUNGBR returns the first n
-*  columns of Q, where m >= n >= k;
-*  if m < k, Q = H(1) H(2) . . . H(m-1) and ZUNGBR returns Q as an
-*  M-by-M matrix.
-*
-*  If VECT = 'P', A is assumed to have been a K-by-N matrix, and P**H
-*  is of order N:
-*  if k < n, P**H = G(k) . . . G(2) G(1) and ZUNGBR returns the first m
-*  rows of P**H, where n >= m >= k;
-*  if k >= n, P**H = G(n-1) . . . G(2) G(1) and ZUNGBR returns P**H as
-*  an N-by-N matrix.
-*
-*  Arguments
-*  =========
-*
-*  VECT    (input) CHARACTER*1
-*          Specifies whether the matrix Q or the matrix P**H is
-*          required, as defined in the transformation applied by ZGEBRD:
-*          = 'Q':  generate Q;
-*          = 'P':  generate P**H.
-*
-*  M       (input) INTEGER
-*          The number of rows of the matrix Q or P**H to be returned.
-*          M >= 0.
-*
-*  N       (input) INTEGER
-*          The number of columns of the matrix Q or P**H to be returned.
-*          N >= 0.
-*          If VECT = 'Q', M >= N >= min(M,K);
-*          if VECT = 'P', N >= M >= min(N,K).
-*
-*  K       (input) INTEGER
-*          If VECT = 'Q', the number of columns in the original M-by-K
-*          matrix reduced by ZGEBRD.
-*          If VECT = 'P', the number of rows in the original K-by-N
-*          matrix reduced by ZGEBRD.
-*          K >= 0.
-*
-*  A       (input/output) COMPLEX*16 array, dimension (LDA,N)
-*          On entry, the vectors which define the elementary reflectors,
-*          as returned by ZGEBRD.
-*          On exit, the M-by-N matrix Q or P**H.
-*
-*  LDA     (input) INTEGER
-*          The leading dimension of the array A. LDA >= M.
-*
-*  TAU     (input) COMPLEX*16 array, dimension
-*                                (min(M,K)) if VECT = 'Q'
-*                                (min(N,K)) if VECT = 'P'
-*          TAU(i) must contain the scalar factor of the elementary
-*          reflector H(i) or G(i), which determines Q or P**H, as
-*          returned by ZGEBRD in its array argument TAUQ or TAUP.
-*
-*  WORK    (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
-*          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
-*
-*  LWORK   (input) INTEGER
-*          The dimension of the array WORK. LWORK >= max(1,min(M,N)).
-*          For optimum performance LWORK >= min(M,N)*NB, where NB
-*          is the optimal blocksize.
-*
-*          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.
-*
-*  INFO    (output) INTEGER
-*          = 0:  successful exit
-*          < 0:  if INFO = -i, the i-th argument had an illegal value
-*
 *  =====================================================================
 *
 *     .. Parameters ..
@@ -140,13 +217,27 @@
       END IF
 *
       IF( INFO.EQ.0 ) THEN
+         WORK( 1 ) = 1
          IF( WANTQ ) THEN
-            NB = ILAENV( 1, 'ZUNGQR', ' ', M, N, K, -1 )
+            IF( M.GE.K ) THEN
+               CALL ZUNGQR( M, N, K, A, LDA, TAU, WORK, -1, IINFO )
+            ELSE
+               IF( M.GT.1 ) THEN
+                  CALL ZUNGQR( M-1, M-1, M-1, A( 2, 2 ), LDA, TAU, WORK,
+     $                         -1, IINFO )
+               END IF
+            END IF
          ELSE
-            NB = ILAENV( 1, 'ZUNGLQ', ' ', M, N, K, -1 )
+            IF( K.LT.N ) THEN
+               CALL ZUNGLQ( M, N, K, A, LDA, TAU, WORK, -1, IINFO )
+            ELSE
+               IF( N.GT.1 ) THEN
+                  CALL ZUNGLQ( N-1, N-1, N-1, A( 2, 2 ), LDA, TAU, WORK,
+     $                         -1, IINFO )
+               END IF
+            END IF
          END IF
-         LWKOPT = MAX( 1, MN )*NB
-         WORK( 1 ) = LWKOPT
+         LWKOPT = WORK( 1 )
       END IF
 *
       IF( INFO.NE.0 ) THEN