LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
zlaunhr_col_getrfnp.f
Go to the documentation of this file.
1 *> \brief \b ZLAUNHR_COL_GETRFNP
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZLAUNHR_COL_GETRFNP + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlaunhr_col_getrfnp.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlaunhr_col_getrfnp.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlaunhr_col_getrfnp.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZLAUNHR_COL_GETRFNP( M, N, A, LDA, D, INFO )
22 *
23 * .. Scalar Arguments ..
24 * INTEGER INFO, LDA, M, N
25 * ..
26 * .. Array Arguments ..
27 * COMPLEX*16 A( LDA, * ), D( * )
28 * ..
29 *
30 *
31 *> \par Purpose:
32 * =============
33 *>
34 *> \verbatim
35 *>
36 *> ZLAUNHR_COL_GETRFNP computes the modified LU factorization without
37 *> pivoting of a complex general M-by-N matrix A. The factorization has
38 *> the form:
39 *>
40 *> A - S = L * U,
41 *>
42 *> where:
43 *> S is a m-by-n diagonal sign matrix with the diagonal D, so that
44 *> D(i) = S(i,i), 1 <= i <= min(M,N). The diagonal D is constructed
45 *> as D(i)=-SIGN(A(i,i)), where A(i,i) is the value after performing
46 *> i-1 steps of Gaussian elimination. This means that the diagonal
47 *> element at each step of "modified" Gaussian elimination is
48 *> at least one in absolute value (so that division-by-zero not
49 *> not possible during the division by the diagonal element);
50 *>
51 *> L is a M-by-N lower triangular matrix with unit diagonal elements
52 *> (lower trapezoidal if M > N);
53 *>
54 *> and U is a M-by-N upper triangular matrix
55 *> (upper trapezoidal if M < N).
56 *>
57 *> This routine is an auxiliary routine used in the Householder
58 *> reconstruction routine ZUNHR_COL. In ZUNHR_COL, this routine is
59 *> applied to an M-by-N matrix A with orthonormal columns, where each
60 *> element is bounded by one in absolute value. With the choice of
61 *> the matrix S above, one can show that the diagonal element at each
62 *> step of Gaussian elimination is the largest (in absolute value) in
63 *> the column on or below the diagonal, so that no pivoting is required
64 *> for numerical stability [1].
65 *>
66 *> For more details on the Householder reconstruction algorithm,
67 *> including the modified LU factorization, see [1].
68 *>
69 *> This is the blocked right-looking version of the algorithm,
70 *> calling Level 3 BLAS to update the submatrix. To factorize a block,
71 *> this routine calls the recursive routine ZLAUNHR_COL_GETRFNP2.
72 *>
73 *> [1] "Reconstructing Householder vectors from tall-skinny QR",
74 *> G. Ballard, J. Demmel, L. Grigori, M. Jacquelin, H.D. Nguyen,
75 *> E. Solomonik, J. Parallel Distrib. Comput.,
76 *> vol. 85, pp. 3-31, 2015.
77 *> \endverbatim
78 *
79 * Arguments:
80 * ==========
81 *
82 *> \param[in] M
83 *> \verbatim
84 *> M is INTEGER
85 *> The number of rows of the matrix A. M >= 0.
86 *> \endverbatim
87 *>
88 *> \param[in] N
89 *> \verbatim
90 *> N is INTEGER
91 *> The number of columns of the matrix A. N >= 0.
92 *> \endverbatim
93 *>
94 *> \param[in,out] A
95 *> \verbatim
96 *> A is COMPLEX*16 array, dimension (LDA,N)
97 *> On entry, the M-by-N matrix to be factored.
98 *> On exit, the factors L and U from the factorization
99 *> A-S=L*U; the unit diagonal elements of L are not stored.
100 *> \endverbatim
101 *>
102 *> \param[in] LDA
103 *> \verbatim
104 *> LDA is INTEGER
105 *> The leading dimension of the array A. LDA >= max(1,M).
106 *> \endverbatim
107 *>
108 *> \param[out] D
109 *> \verbatim
110 *> D is COMPLEX*16 array, dimension min(M,N)
111 *> The diagonal elements of the diagonal M-by-N sign matrix S,
112 *> D(i) = S(i,i), where 1 <= i <= min(M,N). The elements can be
113 *> only ( +1.0, 0.0 ) or (-1.0, 0.0 ).
114 *> \endverbatim
115 *>
116 *> \param[out] INFO
117 *> \verbatim
118 *> INFO is INTEGER
119 *> = 0: successful exit
120 *> < 0: if INFO = -i, the i-th argument had an illegal value
121 *> \endverbatim
122 *>
123 * Authors:
124 * ========
125 *
126 *> \author Univ. of Tennessee
127 *> \author Univ. of California Berkeley
128 *> \author Univ. of Colorado Denver
129 *> \author NAG Ltd.
130 *
131 *> \ingroup complex16GEcomputational
132 *
133 *> \par Contributors:
134 * ==================
135 *>
136 *> \verbatim
137 *>
138 *> November 2019, Igor Kozachenko,
139 *> Computer Science Division,
140 *> University of California, Berkeley
141 *>
142 *> \endverbatim
143 *
144 * =====================================================================
145  SUBROUTINE zlaunhr_col_getrfnp( M, N, A, LDA, D, INFO )
146  IMPLICIT NONE
147 *
148 * -- LAPACK computational routine --
149 * -- LAPACK is a software package provided by Univ. of Tennessee, --
150 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
151 *
152 * .. Scalar Arguments ..
153  INTEGER INFO, LDA, M, N
154 * ..
155 * .. Array Arguments ..
156  COMPLEX*16 A( LDA, * ), D( * )
157 * ..
158 *
159 * =====================================================================
160 *
161 * .. Parameters ..
162  COMPLEX*16 CONE
163  parameter( cone = ( 1.0d+0, 0.0d+0 ) )
164 * ..
165 * .. Local Scalars ..
166  INTEGER IINFO, J, JB, NB
167 * ..
168 * .. External Subroutines ..
170 * ..
171 * .. External Functions ..
172  INTEGER ILAENV
173  EXTERNAL ilaenv
174 * ..
175 * .. Intrinsic Functions ..
176  INTRINSIC max, min
177 * ..
178 * .. Executable Statements ..
179 *
180 * Test the input parameters.
181 *
182  info = 0
183  IF( m.LT.0 ) THEN
184  info = -1
185  ELSE IF( n.LT.0 ) THEN
186  info = -2
187  ELSE IF( lda.LT.max( 1, m ) ) THEN
188  info = -4
189  END IF
190  IF( info.NE.0 ) THEN
191  CALL xerbla( 'ZLAUNHR_COL_GETRFNP', -info )
192  RETURN
193  END IF
194 *
195 * Quick return if possible
196 *
197  IF( min( m, n ).EQ.0 )
198  $ RETURN
199 *
200 * Determine the block size for this environment.
201 *
202 
203  nb = ilaenv( 1, 'ZLAUNHR_COL_GETRFNP', ' ', m, n, -1, -1 )
204 
205  IF( nb.LE.1 .OR. nb.GE.min( m, n ) ) THEN
206 *
207 * Use unblocked code.
208 *
209  CALL zlaunhr_col_getrfnp2( m, n, a, lda, d, info )
210  ELSE
211 *
212 * Use blocked code.
213 *
214  DO j = 1, min( m, n ), nb
215  jb = min( min( m, n )-j+1, nb )
216 *
217 * Factor diagonal and subdiagonal blocks.
218 *
219  CALL zlaunhr_col_getrfnp2( m-j+1, jb, a( j, j ), lda,
220  $ d( j ), iinfo )
221 *
222  IF( j+jb.LE.n ) THEN
223 *
224 * Compute block row of U.
225 *
226  CALL ztrsm( 'Left', 'Lower', 'No transpose', 'Unit', jb,
227  $ n-j-jb+1, cone, a( j, j ), lda, a( j, j+jb ),
228  $ lda )
229  IF( j+jb.LE.m ) THEN
230 *
231 * Update trailing submatrix.
232 *
233  CALL zgemm( 'No transpose', 'No transpose', m-j-jb+1,
234  $ n-j-jb+1, jb, -cone, a( j+jb, j ), lda,
235  $ a( j, j+jb ), lda, cone, a( j+jb, j+jb ),
236  $ lda )
237  END IF
238  END IF
239  END DO
240  END IF
241  RETURN
242 *
243 * End of ZLAUNHR_COL_GETRFNP
244 *
245  END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:187
subroutine ztrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
ZTRSM
Definition: ztrsm.f:180
recursive subroutine zlaunhr_col_getrfnp2(M, N, A, LDA, D, INFO)
ZLAUNHR_COL_GETRFNP2
subroutine zlaunhr_col_getrfnp(M, N, A, LDA, D, INFO)
ZLAUNHR_COL_GETRFNP