LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dtpqrt.f
Go to the documentation of this file.
00001 *> \brief \b DTPQRT
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download DTPQRT + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtpqrt.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtpqrt.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtpqrt.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE DTPQRT( M, N, L, NB, A, LDA, B, LDB, T, LDT, WORK,
00022 *                          INFO )
00023 * 
00024 *       .. Scalar Arguments ..
00025 *       INTEGER INFO, LDA, LDB, LDT, N, M, L, NB
00026 *       ..
00027 *       .. Array Arguments ..
00028 *       DOUBLE PRECISION A( LDA, * ), B( LDB, * ), T( LDT, * ), WORK( * )
00029 *       ..
00030 *  
00031 *
00032 *> \par Purpose:
00033 *  =============
00034 *>
00035 *> \verbatim
00036 *>
00037 *> DTPQRT computes a blocked QR factorization of a real 
00038 *> "triangular-pentagonal" matrix C, which is composed of a 
00039 *> triangular block A and pentagonal block B, using the compact 
00040 *> WY representation for Q.
00041 *> \endverbatim
00042 *
00043 *  Arguments:
00044 *  ==========
00045 *
00046 *> \param[in] M
00047 *> \verbatim
00048 *>          M is INTEGER
00049 *>          The number of rows of the matrix B.  
00050 *>          M >= 0.
00051 *> \endverbatim
00052 *>
00053 *> \param[in] N
00054 *> \verbatim
00055 *>          N is INTEGER
00056 *>          The number of columns of the matrix B, and the order of the
00057 *>          triangular matrix A.
00058 *>          N >= 0.
00059 *> \endverbatim
00060 *>
00061 *> \param[in] L
00062 *> \verbatim
00063 *>          L is INTEGER
00064 *>          The number of rows of the upper trapezoidal part of B.
00065 *>          MIN(M,N) >= L >= 0.  See Further Details.
00066 *> \endverbatim
00067 *>
00068 *> \param[in] NB
00069 *> \verbatim
00070 *>          NB is INTEGER
00071 *>          The block size to be used in the blocked QR.  N >= NB >= 1.
00072 *> \endverbatim
00073 *>
00074 *> \param[in,out] A
00075 *> \verbatim
00076 *>          A is DOUBLE PRECISION array, dimension (LDA,N)
00077 *>          On entry, the upper triangular N-by-N matrix A.
00078 *>          On exit, the elements on and above the diagonal of the array
00079 *>          contain the upper triangular matrix R.
00080 *> \endverbatim
00081 *>
00082 *> \param[in] LDA
00083 *> \verbatim
00084 *>          LDA is INTEGER
00085 *>          The leading dimension of the array A.  LDA >= max(1,N).
00086 *> \endverbatim
00087 *>
00088 *> \param[in,out] B
00089 *> \verbatim
00090 *>          B is DOUBLE PRECISION array, dimension (LDB,N)
00091 *>          On entry, the pentagonal M-by-N matrix B.  The first M-L rows 
00092 *>          are rectangular, and the last L rows are upper trapezoidal.
00093 *>          On exit, B contains the pentagonal matrix V.  See Further Details.
00094 *> \endverbatim
00095 *>
00096 *> \param[in] LDB
00097 *> \verbatim
00098 *>          LDB is INTEGER
00099 *>          The leading dimension of the array B.  LDB >= max(1,M).
00100 *> \endverbatim
00101 *>
00102 *> \param[out] T
00103 *> \verbatim
00104 *>          T is DOUBLE PRECISION array, dimension (LDT,N)
00105 *>          The upper triangular block reflectors stored in compact form
00106 *>          as a sequence of upper triangular blocks.  See Further Details.
00107 *> \endverbatim
00108 *>          
00109 *> \param[in] LDT
00110 *> \verbatim
00111 *>          LDT is INTEGER
00112 *>          The leading dimension of the array T.  LDT >= NB.
00113 *> \endverbatim
00114 *>
00115 *> \param[out] WORK
00116 *> \verbatim
00117 *>          WORK is DOUBLE PRECISION array, dimension (NB*N)
00118 *> \endverbatim
00119 *>
00120 *> \param[out] INFO
00121 *> \verbatim
00122 *>          INFO is INTEGER
00123 *>          = 0:  successful exit
00124 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
00125 *> \endverbatim
00126 *
00127 *  Authors:
00128 *  ========
00129 *
00130 *> \author Univ. of Tennessee 
00131 *> \author Univ. of California Berkeley 
00132 *> \author Univ. of Colorado Denver 
00133 *> \author NAG Ltd. 
00134 *
00135 *> \date April 2012
00136 *
00137 *> \ingroup doubleOTHERcomputational
00138 *
00139 *> \par Further Details:
00140 *  =====================
00141 *>
00142 *> \verbatim
00143 *>
00144 *>  The input matrix C is a (N+M)-by-N matrix  
00145 *>
00146 *>               C = [ A ]
00147 *>                   [ B ]        
00148 *>
00149 *>  where A is an upper triangular N-by-N matrix, and B is M-by-N pentagonal
00150 *>  matrix consisting of a (M-L)-by-N rectangular matrix B1 on top of a L-by-N
00151 *>  upper trapezoidal matrix B2:
00152 *>
00153 *>               B = [ B1 ]  <- (M-L)-by-N rectangular
00154 *>                   [ B2 ]  <-     L-by-N upper trapezoidal.
00155 *>
00156 *>  The upper trapezoidal matrix B2 consists of the first L rows of a
00157 *>  N-by-N upper triangular matrix, where 0 <= L <= MIN(M,N).  If L=0, 
00158 *>  B is rectangular M-by-N; if M=L=N, B is upper triangular.  
00159 *>
00160 *>  The matrix W stores the elementary reflectors H(i) in the i-th column
00161 *>  below the diagonal (of A) in the (N+M)-by-N input matrix C
00162 *>
00163 *>               C = [ A ]  <- upper triangular N-by-N
00164 *>                   [ B ]  <- M-by-N pentagonal
00165 *>
00166 *>  so that W can be represented as
00167 *>
00168 *>               W = [ I ]  <- identity, N-by-N
00169 *>                   [ V ]  <- M-by-N, same form as B.
00170 *>
00171 *>  Thus, all of information needed for W is contained on exit in B, which
00172 *>  we call V above.  Note that V has the same form as B; that is, 
00173 *>
00174 *>               V = [ V1 ] <- (M-L)-by-N rectangular
00175 *>                   [ V2 ] <-     L-by-N upper trapezoidal.
00176 *>
00177 *>  The columns of V represent the vectors which define the H(i)'s.  
00178 *>
00179 *>  The number of blocks is B = ceiling(N/NB), where each
00180 *>  block is of order NB except for the last block, which is of order 
00181 *>  IB = N - (B-1)*NB.  For each of the B blocks, a upper triangular block
00182 *>  reflector factor is computed: T1, T2, ..., TB.  The NB-by-NB (and IB-by-IB 
00183 *>  for the last block) T's are stored in the NB-by-N matrix T as
00184 *>
00185 *>               T = [T1 T2 ... TB].
00186 *> \endverbatim
00187 *>
00188 *  =====================================================================
00189       SUBROUTINE DTPQRT( M, N, L, NB, A, LDA, B, LDB, T, LDT, WORK,
00190      $                   INFO )
00191 *
00192 *  -- LAPACK computational routine (version 3.4.1) --
00193 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00194 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00195 *     April 2012
00196 *
00197 *     .. Scalar Arguments ..
00198       INTEGER INFO, LDA, LDB, LDT, N, M, L, NB
00199 *     ..
00200 *     .. Array Arguments ..
00201       DOUBLE PRECISION A( LDA, * ), B( LDB, * ), T( LDT, * ), WORK( * )
00202 *     ..
00203 *
00204 * =====================================================================
00205 *
00206 *     ..
00207 *     .. Local Scalars ..
00208       INTEGER    I, IB, LB, MB, IINFO
00209 *     ..
00210 *     .. External Subroutines ..
00211       EXTERNAL   DTPQRT2, DTPRFB, XERBLA
00212 *     ..
00213 *     .. Executable Statements ..
00214 *
00215 *     Test the input arguments
00216 *
00217       INFO = 0
00218       IF( M.LT.0 ) THEN
00219          INFO = -1
00220       ELSE IF( N.LT.0 ) THEN
00221          INFO = -2
00222       ELSE IF( L.LT.0 .OR. L.GT.MIN(M,N) ) THEN
00223          INFO = -3
00224       ELSE IF( NB.LT.1 .OR. NB.GT.N ) THEN
00225          INFO = -4
00226       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
00227          INFO = -6
00228       ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
00229          INFO = -8
00230       ELSE IF( LDT.LT.NB ) THEN
00231          INFO = -10
00232       END IF
00233       IF( INFO.NE.0 ) THEN
00234          CALL XERBLA( 'DTPQRT', -INFO )
00235          RETURN
00236       END IF
00237 *
00238 *     Quick return if possible
00239 *
00240       IF( M.EQ.0 .OR. N.EQ.0 ) RETURN
00241 *
00242       DO I = 1, N, NB
00243 *     
00244 *     Compute the QR factorization of the current block
00245 *
00246          IB = MIN( N-I+1, NB )
00247          MB = MIN( M-L+I+IB-1, M )
00248          IF( I.GE.L ) THEN
00249             LB = 0
00250          ELSE
00251             LB = MB-M+L-I+1
00252          END IF
00253 *
00254          CALL DTPQRT2( MB, IB, LB, A(I,I), LDA, B( 1, I ), LDB, 
00255      $                 T(1, I ), LDT, IINFO )
00256 *
00257 *     Update by applying H**T to B(:,I+IB:N) from the left
00258 *
00259          IF( I+IB.LE.N ) THEN
00260             CALL DTPRFB( 'L', 'T', 'F', 'C', MB, N-I-IB+1, IB, LB,
00261      $                    B( 1, I ), LDB, T( 1, I ), LDT, 
00262      $                    A( I, I+IB ), LDA, B( 1, I+IB ), LDB, 
00263      $                    WORK, IB )
00264          END IF
00265       END DO
00266       RETURN
00267 *     
00268 *     End of DTPQRT
00269 *
00270       END
 All Files Functions