LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlaorhr_col_getrfnp.f
Go to the documentation of this file.
1*> \brief \b DLAORHR_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 DLAORHR_COL_GETRFNP + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaorhr_col_getrfnp.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaorhr_col_getrfnp.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaorhr_col_getrfnp.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DLAORHR_COL_GETRFNP( M, N, A, LDA, D, INFO )
22*
23* .. Scalar Arguments ..
24* INTEGER INFO, LDA, M, N
25* ..
26* .. Array Arguments ..
27* DOUBLE PRECISION A( LDA, * ), D( * )
28* ..
29*
30*
31*> \par Purpose:
32* =============
33*>
34*> \verbatim
35*>
36*> DLAORHR_COL_GETRFNP computes the modified LU factorization without
37*> pivoting of a real 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 DORHR_COL. In DORHR_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 DLAORHR_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 DOUBLE PRECISION 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 DOUBLE PRECISION 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
113*> be only plus or minus one.
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 launhr_col_getrfnp
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 dlaorhr_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 DOUBLE PRECISION A( LDA, * ), D( * )
157* ..
158*
159* =====================================================================
160*
161* .. Parameters ..
162 DOUBLE PRECISION ONE
163 parameter( one = 1.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( 'DLAORHR_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, 'DLAORHR_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 dlaorhr_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 dlaorhr_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 dtrsm( 'Left', 'Lower', 'No transpose', 'Unit', jb,
227 $ n-j-jb+1, one, 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 dgemm( 'No transpose', 'No transpose', m-j-jb+1,
234 $ n-j-jb+1, jb, -one, a( j+jb, j ), lda,
235 $ a( j, j+jb ), lda, one, 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 DLAORHR_COL_GETRFNP
244*
245 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
DGEMM
Definition dgemm.f:188
recursive subroutine dlaorhr_col_getrfnp2(m, n, a, lda, d, info)
DLAORHR_COL_GETRFNP2
subroutine dlaorhr_col_getrfnp(m, n, a, lda, d, info)
DLAORHR_COL_GETRFNP
subroutine dtrsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
DTRSM
Definition dtrsm.f:181