00001 /* ctzrqf.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 complex c_b1 = {1.f,0.f}; 00019 static integer c__1 = 1; 00020 00021 /* Subroutine */ int ctzrqf_(integer *m, integer *n, complex *a, integer *lda, 00022 complex *tau, integer *info) 00023 { 00024 /* System generated locals */ 00025 integer a_dim1, a_offset, i__1, i__2; 00026 complex q__1, q__2; 00027 00028 /* Builtin functions */ 00029 void r_cnjg(complex *, complex *); 00030 00031 /* Local variables */ 00032 integer i__, k, m1; 00033 extern /* Subroutine */ int cgerc_(integer *, integer *, complex *, 00034 complex *, integer *, complex *, integer *, complex *, integer *); 00035 complex alpha; 00036 extern /* Subroutine */ int cgemv_(char *, integer *, integer *, complex * 00037 , complex *, integer *, complex *, integer *, complex *, complex * 00038 , integer *), ccopy_(integer *, complex *, integer *, 00039 complex *, integer *), caxpy_(integer *, complex *, complex *, 00040 integer *, complex *, integer *), clacgv_(integer *, complex *, 00041 integer *), clarfp_(integer *, complex *, complex *, integer *, 00042 complex *), xerbla_(char *, integer *); 00043 00044 00045 /* -- LAPACK routine (version 3.2) -- */ 00046 /* -- LAPACK is a software package provided by Univ. of Tennessee, -- */ 00047 /* -- Univ. of California Berkeley, Univ. of Colorado Denver 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 CTZRZF. */ 00059 00060 /* CTZRQF 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 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 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_("CTZRQF", &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.f, tau[i__2].i = 0.f; 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 r_cnjg(&q__1, &a[k + k * a_dim1]); 00184 a[i__1].r = q__1.r, a[i__1].i = q__1.i; 00185 i__1 = *n - *m; 00186 clacgv_(&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 clarfp_(&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 r_cnjg(&q__1, &tau[k]); 00195 tau[i__1].r = q__1.r, tau[i__1].i = q__1.i; 00196 00197 i__1 = k; 00198 if ((tau[i__1].r != 0.f || tau[i__1].i != 0.f) && 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 ccopy_(&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 cgemv_("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 r_cnjg(&q__2, &tau[k]); 00223 q__1.r = -q__2.r, q__1.i = -q__2.i; 00224 caxpy_(&i__1, &q__1, &tau[1], &c__1, &a[k * a_dim1 + 1], & 00225 c__1); 00226 i__1 = k - 1; 00227 i__2 = *n - *m; 00228 r_cnjg(&q__2, &tau[k]); 00229 q__1.r = -q__2.r, q__1.i = -q__2.i; 00230 cgerc_(&i__1, &i__2, &q__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 CTZRQF */ 00240 00241 } /* ctzrqf_ */