00001 /* ztzrqf.f -- translated by f2c (version 20061008). 00002 You must link the resulting object file with libf2c: 00003 on Microsoft Windows system, link with libf2c.lib; 00004 on Linux or Unix systems, link with .../path/to/libf2c.a -lm 00005 or, if you install libf2c.a in a standard place, with -lf2c -lm 00006 -- in that order, at the end of the command line, as in 00007 cc *.o -lf2c -lm 00008 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g., 00009 00010 http://www.netlib.org/f2c/libf2c.zip 00011 */ 00012 00013 #include "f2c.h" 00014 #include "blaswrap.h" 00015 00016 /* Table of constant values */ 00017 00018 static doublecomplex c_b1 = {1.,0.}; 00019 static integer c__1 = 1; 00020 00021 /* Subroutine */ int ztzrqf_(integer *m, integer *n, doublecomplex *a, 00022 integer *lda, doublecomplex *tau, integer *info) 00023 { 00024 /* System generated locals */ 00025 integer a_dim1, a_offset, i__1, i__2; 00026 doublecomplex z__1, z__2; 00027 00028 /* Builtin functions */ 00029 void d_cnjg(doublecomplex *, doublecomplex *); 00030 00031 /* Local variables */ 00032 integer i__, k, m1; 00033 doublecomplex alpha; 00034 extern /* Subroutine */ int zgerc_(integer *, integer *, doublecomplex *, 00035 doublecomplex *, integer *, doublecomplex *, integer *, 00036 doublecomplex *, integer *), zgemv_(char *, integer *, integer *, 00037 doublecomplex *, doublecomplex *, integer *, doublecomplex *, 00038 integer *, doublecomplex *, doublecomplex *, integer *), 00039 zcopy_(integer *, doublecomplex *, integer *, doublecomplex *, 00040 integer *), zaxpy_(integer *, doublecomplex *, doublecomplex *, 00041 integer *, doublecomplex *, integer *), xerbla_(char *, integer *), zlacgv_(integer *, doublecomplex *, integer *), zlarfp_( 00042 integer *, doublecomplex *, doublecomplex *, integer *, 00043 doublecomplex *); 00044 00045 00046 /* -- LAPACK routine (version 3.2) -- */ 00047 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ 00048 /* November 2006 */ 00049 00050 /* .. Scalar Arguments .. */ 00051 /* .. */ 00052 /* .. Array Arguments .. */ 00053 /* .. */ 00054 00055 /* Purpose */ 00056 /* ======= */ 00057 00058 /* This routine is deprecated and has been replaced by routine ZTZRZF. */ 00059 00060 /* ZTZRQF reduces the M-by-N ( M<=N ) complex upper trapezoidal matrix A */ 00061 /* to upper triangular form by means of unitary transformations. */ 00062 00063 /* The upper trapezoidal matrix A is factored as */ 00064 00065 /* A = ( R 0 ) * Z, */ 00066 00067 /* where Z is an N-by-N unitary matrix and R is an M-by-M upper */ 00068 /* triangular matrix. */ 00069 00070 /* Arguments */ 00071 /* ========= */ 00072 00073 /* M (input) INTEGER */ 00074 /* The number of rows of the matrix A. M >= 0. */ 00075 00076 /* N (input) INTEGER */ 00077 /* The number of columns of the matrix A. N >= M. */ 00078 00079 /* A (input/output) COMPLEX*16 array, dimension (LDA,N) */ 00080 /* On entry, the leading M-by-N upper trapezoidal part of the */ 00081 /* array A must contain the matrix to be factorized. */ 00082 /* On exit, the leading M-by-M upper triangular part of A */ 00083 /* contains the upper triangular matrix R, and elements M+1 to */ 00084 /* N of the first M rows of A, with the array TAU, represent the */ 00085 /* unitary matrix Z as a product of M elementary reflectors. */ 00086 00087 /* LDA (input) INTEGER */ 00088 /* The leading dimension of the array A. LDA >= max(1,M). */ 00089 00090 /* TAU (output) COMPLEX*16 array, dimension (M) */ 00091 /* The scalar factors of the elementary reflectors. */ 00092 00093 /* INFO (output) INTEGER */ 00094 /* = 0: successful exit */ 00095 /* < 0: if INFO = -i, the i-th argument had an illegal value */ 00096 00097 /* Further Details */ 00098 /* =============== */ 00099 00100 /* The factorization is obtained by Householder's method. The kth */ 00101 /* transformation matrix, Z( k ), whose conjugate transpose is used to */ 00102 /* introduce zeros into the (m - k + 1)th row of A, is given in the form */ 00103 00104 /* Z( k ) = ( I 0 ), */ 00105 /* ( 0 T( k ) ) */ 00106 00107 /* where */ 00108 00109 /* T( k ) = I - tau*u( k )*u( k )', u( k ) = ( 1 ), */ 00110 /* ( 0 ) */ 00111 /* ( z( k ) ) */ 00112 00113 /* tau is a scalar and z( k ) is an ( n - m ) element vector. */ 00114 /* tau and z( k ) are chosen to annihilate the elements of the kth row */ 00115 /* of X. */ 00116 00117 /* The scalar tau is returned in the kth element of TAU and the vector */ 00118 /* u( k ) in the kth row of A, such that the elements of z( k ) are */ 00119 /* in a( k, m + 1 ), ..., a( k, n ). The elements of R are returned in */ 00120 /* the upper triangular part of A. */ 00121 00122 /* Z is given by */ 00123 00124 /* Z = Z( 1 ) * Z( 2 ) * ... * Z( m ). */ 00125 00126 /* ===================================================================== */ 00127 00128 /* .. Parameters .. */ 00129 /* .. */ 00130 /* .. Local Scalars .. */ 00131 /* .. */ 00132 /* .. Intrinsic Functions .. */ 00133 /* .. */ 00134 /* .. External Subroutines .. */ 00135 /* .. */ 00136 /* .. Executable Statements .. */ 00137 00138 /* Test the input parameters. */ 00139 00140 /* Parameter adjustments */ 00141 a_dim1 = *lda; 00142 a_offset = 1 + a_dim1; 00143 a -= a_offset; 00144 --tau; 00145 00146 /* Function Body */ 00147 *info = 0; 00148 if (*m < 0) { 00149 *info = -1; 00150 } else if (*n < *m) { 00151 *info = -2; 00152 } else if (*lda < max(1,*m)) { 00153 *info = -4; 00154 } 00155 if (*info != 0) { 00156 i__1 = -(*info); 00157 xerbla_("ZTZRQF", &i__1); 00158 return 0; 00159 } 00160 00161 /* Perform the factorization. */ 00162 00163 if (*m == 0) { 00164 return 0; 00165 } 00166 if (*m == *n) { 00167 i__1 = *n; 00168 for (i__ = 1; i__ <= i__1; ++i__) { 00169 i__2 = i__; 00170 tau[i__2].r = 0., tau[i__2].i = 0.; 00171 /* L10: */ 00172 } 00173 } else { 00174 /* Computing MIN */ 00175 i__1 = *m + 1; 00176 m1 = min(i__1,*n); 00177 for (k = *m; k >= 1; --k) { 00178 00179 /* Use a Householder reflection to zero the kth row of A. */ 00180 /* First set up the reflection. */ 00181 00182 i__1 = k + k * a_dim1; 00183 d_cnjg(&z__1, &a[k + k * a_dim1]); 00184 a[i__1].r = z__1.r, a[i__1].i = z__1.i; 00185 i__1 = *n - *m; 00186 zlacgv_(&i__1, &a[k + m1 * a_dim1], lda); 00187 i__1 = k + k * a_dim1; 00188 alpha.r = a[i__1].r, alpha.i = a[i__1].i; 00189 i__1 = *n - *m + 1; 00190 zlarfp_(&i__1, &alpha, &a[k + m1 * a_dim1], lda, &tau[k]); 00191 i__1 = k + k * a_dim1; 00192 a[i__1].r = alpha.r, a[i__1].i = alpha.i; 00193 i__1 = k; 00194 d_cnjg(&z__1, &tau[k]); 00195 tau[i__1].r = z__1.r, tau[i__1].i = z__1.i; 00196 00197 i__1 = k; 00198 if ((tau[i__1].r != 0. || tau[i__1].i != 0.) && k > 1) { 00199 00200 /* We now perform the operation A := A*P( k )'. */ 00201 00202 /* Use the first ( k - 1 ) elements of TAU to store a( k ), */ 00203 /* where a( k ) consists of the first ( k - 1 ) elements of */ 00204 /* the kth column of A. Also let B denote the first */ 00205 /* ( k - 1 ) rows of the last ( n - m ) columns of A. */ 00206 00207 i__1 = k - 1; 00208 zcopy_(&i__1, &a[k * a_dim1 + 1], &c__1, &tau[1], &c__1); 00209 00210 /* Form w = a( k ) + B*z( k ) in TAU. */ 00211 00212 i__1 = k - 1; 00213 i__2 = *n - *m; 00214 zgemv_("No transpose", &i__1, &i__2, &c_b1, &a[m1 * a_dim1 + 00215 1], lda, &a[k + m1 * a_dim1], lda, &c_b1, &tau[1], & 00216 c__1); 00217 00218 /* Now form a( k ) := a( k ) - conjg(tau)*w */ 00219 /* and B := B - conjg(tau)*w*z( k )'. */ 00220 00221 i__1 = k - 1; 00222 d_cnjg(&z__2, &tau[k]); 00223 z__1.r = -z__2.r, z__1.i = -z__2.i; 00224 zaxpy_(&i__1, &z__1, &tau[1], &c__1, &a[k * a_dim1 + 1], & 00225 c__1); 00226 i__1 = k - 1; 00227 i__2 = *n - *m; 00228 d_cnjg(&z__2, &tau[k]); 00229 z__1.r = -z__2.r, z__1.i = -z__2.i; 00230 zgerc_(&i__1, &i__2, &z__1, &tau[1], &c__1, &a[k + m1 * 00231 a_dim1], lda, &a[m1 * a_dim1 + 1], lda); 00232 } 00233 /* L20: */ 00234 } 00235 } 00236 00237 return 0; 00238 00239 /* End of ZTZRQF */ 00240 00241 } /* ztzrqf_ */