LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ slagge()

subroutine slagge ( integer  M,
integer  N,
integer  KL,
integer  KU,
real, dimension( * )  D,
real, dimension( lda, * )  A,
integer  LDA,
integer, dimension( 4 )  ISEED,
real, dimension( * )  WORK,
integer  INFO 
)

SLAGGE

Purpose:
 SLAGGE generates a real general m by n matrix A, by pre- and post-
 multiplying a real diagonal matrix D with random orthogonal matrices:
 A = U*D*V. The lower and upper bandwidths may then be reduced to
 kl and ku by additional orthogonal transformations.
Parameters
[in]M
          M is INTEGER
          The number of rows of the matrix A.  M >= 0.
[in]N
          N is INTEGER
          The number of columns of the matrix A.  N >= 0.
[in]KL
          KL is INTEGER
          The number of nonzero subdiagonals within the band of A.
          0 <= KL <= M-1.
[in]KU
          KU is INTEGER
          The number of nonzero superdiagonals within the band of A.
          0 <= KU <= N-1.
[in]D
          D is REAL array, dimension (min(M,N))
          The diagonal elements of the diagonal matrix D.
[out]A
          A is REAL array, dimension (LDA,N)
          The generated m by n matrix A.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= M.
[in,out]ISEED
          ISEED is INTEGER array, dimension (4)
          On entry, the seed of the random number generator; the array
          elements must be between 0 and 4095, and ISEED(4) must be
          odd.
          On exit, the seed is updated.
[out]WORK
          WORK is REAL array, dimension (M+N)
[out]INFO
          INFO is INTEGER
          = 0: successful exit
          < 0: if INFO = -i, the i-th argument had an illegal value
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 112 of file slagge.f.

113 *
114 * -- LAPACK auxiliary routine --
115 * -- LAPACK is a software package provided by Univ. of Tennessee, --
116 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
117 *
118 * .. Scalar Arguments ..
119  INTEGER INFO, KL, KU, LDA, M, N
120 * ..
121 * .. Array Arguments ..
122  INTEGER ISEED( 4 )
123  REAL A( LDA, * ), D( * ), WORK( * )
124 * ..
125 *
126 * =====================================================================
127 *
128 * .. Parameters ..
129  REAL ZERO, ONE
130  parameter( zero = 0.0e+0, one = 1.0e+0 )
131 * ..
132 * .. Local Scalars ..
133  INTEGER I, J
134  REAL TAU, WA, WB, WN
135 * ..
136 * .. External Subroutines ..
137  EXTERNAL sgemv, sger, slarnv, sscal, xerbla
138 * ..
139 * .. Intrinsic Functions ..
140  INTRINSIC max, min, sign
141 * ..
142 * .. External Functions ..
143  REAL SNRM2
144  EXTERNAL snrm2
145 * ..
146 * .. Executable Statements ..
147 *
148 * Test the input arguments
149 *
150  info = 0
151  IF( m.LT.0 ) THEN
152  info = -1
153  ELSE IF( n.LT.0 ) THEN
154  info = -2
155  ELSE IF( kl.LT.0 .OR. kl.GT.m-1 ) THEN
156  info = -3
157  ELSE IF( ku.LT.0 .OR. ku.GT.n-1 ) THEN
158  info = -4
159  ELSE IF( lda.LT.max( 1, m ) ) THEN
160  info = -7
161  END IF
162  IF( info.LT.0 ) THEN
163  CALL xerbla( 'SLAGGE', -info )
164  RETURN
165  END IF
166 *
167 * initialize A to diagonal matrix
168 *
169  DO 20 j = 1, n
170  DO 10 i = 1, m
171  a( i, j ) = zero
172  10 CONTINUE
173  20 CONTINUE
174  DO 30 i = 1, min( m, n )
175  a( i, i ) = d( i )
176  30 CONTINUE
177 *
178 * Quick exit if the user wants a diagonal matrix
179 *
180  IF(( kl .EQ. 0 ).AND.( ku .EQ. 0)) RETURN
181 *
182 * pre- and post-multiply A by random orthogonal matrices
183 *
184  DO 40 i = min( m, n ), 1, -1
185  IF( i.LT.m ) THEN
186 *
187 * generate random reflection
188 *
189  CALL slarnv( 3, iseed, m-i+1, work )
190  wn = snrm2( m-i+1, work, 1 )
191  wa = sign( wn, work( 1 ) )
192  IF( wn.EQ.zero ) THEN
193  tau = zero
194  ELSE
195  wb = work( 1 ) + wa
196  CALL sscal( m-i, one / wb, work( 2 ), 1 )
197  work( 1 ) = one
198  tau = wb / wa
199  END IF
200 *
201 * multiply A(i:m,i:n) by random reflection from the left
202 *
203  CALL sgemv( 'Transpose', m-i+1, n-i+1, one, a( i, i ), lda,
204  $ work, 1, zero, work( m+1 ), 1 )
205  CALL sger( m-i+1, n-i+1, -tau, work, 1, work( m+1 ), 1,
206  $ a( i, i ), lda )
207  END IF
208  IF( i.LT.n ) THEN
209 *
210 * generate random reflection
211 *
212  CALL slarnv( 3, iseed, n-i+1, work )
213  wn = snrm2( n-i+1, work, 1 )
214  wa = sign( wn, work( 1 ) )
215  IF( wn.EQ.zero ) THEN
216  tau = zero
217  ELSE
218  wb = work( 1 ) + wa
219  CALL sscal( n-i, one / wb, work( 2 ), 1 )
220  work( 1 ) = one
221  tau = wb / wa
222  END IF
223 *
224 * multiply A(i:m,i:n) by random reflection from the right
225 *
226  CALL sgemv( 'No transpose', m-i+1, n-i+1, one, a( i, i ),
227  $ lda, work, 1, zero, work( n+1 ), 1 )
228  CALL sger( m-i+1, n-i+1, -tau, work( n+1 ), 1, work, 1,
229  $ a( i, i ), lda )
230  END IF
231  40 CONTINUE
232 *
233 * Reduce number of subdiagonals to KL and number of superdiagonals
234 * to KU
235 *
236  DO 70 i = 1, max( m-1-kl, n-1-ku )
237  IF( kl.LE.ku ) THEN
238 *
239 * annihilate subdiagonal elements first (necessary if KL = 0)
240 *
241  IF( i.LE.min( m-1-kl, n ) ) THEN
242 *
243 * generate reflection to annihilate A(kl+i+1:m,i)
244 *
245  wn = snrm2( m-kl-i+1, a( kl+i, i ), 1 )
246  wa = sign( wn, a( kl+i, i ) )
247  IF( wn.EQ.zero ) THEN
248  tau = zero
249  ELSE
250  wb = a( kl+i, i ) + wa
251  CALL sscal( m-kl-i, one / wb, a( kl+i+1, i ), 1 )
252  a( kl+i, i ) = one
253  tau = wb / wa
254  END IF
255 *
256 * apply reflection to A(kl+i:m,i+1:n) from the left
257 *
258  CALL sgemv( 'Transpose', m-kl-i+1, n-i, one,
259  $ a( kl+i, i+1 ), lda, a( kl+i, i ), 1, zero,
260  $ work, 1 )
261  CALL sger( m-kl-i+1, n-i, -tau, a( kl+i, i ), 1, work, 1,
262  $ a( kl+i, i+1 ), lda )
263  a( kl+i, i ) = -wa
264  END IF
265 *
266  IF( i.LE.min( n-1-ku, m ) ) THEN
267 *
268 * generate reflection to annihilate A(i,ku+i+1:n)
269 *
270  wn = snrm2( n-ku-i+1, a( i, ku+i ), lda )
271  wa = sign( wn, a( i, ku+i ) )
272  IF( wn.EQ.zero ) THEN
273  tau = zero
274  ELSE
275  wb = a( i, ku+i ) + wa
276  CALL sscal( n-ku-i, one / wb, a( i, ku+i+1 ), lda )
277  a( i, ku+i ) = one
278  tau = wb / wa
279  END IF
280 *
281 * apply reflection to A(i+1:m,ku+i:n) from the right
282 *
283  CALL sgemv( 'No transpose', m-i, n-ku-i+1, one,
284  $ a( i+1, ku+i ), lda, a( i, ku+i ), lda, zero,
285  $ work, 1 )
286  CALL sger( m-i, n-ku-i+1, -tau, work, 1, a( i, ku+i ),
287  $ lda, a( i+1, ku+i ), lda )
288  a( i, ku+i ) = -wa
289  END IF
290  ELSE
291 *
292 * annihilate superdiagonal elements first (necessary if
293 * KU = 0)
294 *
295  IF( i.LE.min( n-1-ku, m ) ) THEN
296 *
297 * generate reflection to annihilate A(i,ku+i+1:n)
298 *
299  wn = snrm2( n-ku-i+1, a( i, ku+i ), lda )
300  wa = sign( wn, a( i, ku+i ) )
301  IF( wn.EQ.zero ) THEN
302  tau = zero
303  ELSE
304  wb = a( i, ku+i ) + wa
305  CALL sscal( n-ku-i, one / wb, a( i, ku+i+1 ), lda )
306  a( i, ku+i ) = one
307  tau = wb / wa
308  END IF
309 *
310 * apply reflection to A(i+1:m,ku+i:n) from the right
311 *
312  CALL sgemv( 'No transpose', m-i, n-ku-i+1, one,
313  $ a( i+1, ku+i ), lda, a( i, ku+i ), lda, zero,
314  $ work, 1 )
315  CALL sger( m-i, n-ku-i+1, -tau, work, 1, a( i, ku+i ),
316  $ lda, a( i+1, ku+i ), lda )
317  a( i, ku+i ) = -wa
318  END IF
319 *
320  IF( i.LE.min( m-1-kl, n ) ) THEN
321 *
322 * generate reflection to annihilate A(kl+i+1:m,i)
323 *
324  wn = snrm2( m-kl-i+1, a( kl+i, i ), 1 )
325  wa = sign( wn, a( kl+i, i ) )
326  IF( wn.EQ.zero ) THEN
327  tau = zero
328  ELSE
329  wb = a( kl+i, i ) + wa
330  CALL sscal( m-kl-i, one / wb, a( kl+i+1, i ), 1 )
331  a( kl+i, i ) = one
332  tau = wb / wa
333  END IF
334 *
335 * apply reflection to A(kl+i:m,i+1:n) from the left
336 *
337  CALL sgemv( 'Transpose', m-kl-i+1, n-i, one,
338  $ a( kl+i, i+1 ), lda, a( kl+i, i ), 1, zero,
339  $ work, 1 )
340  CALL sger( m-kl-i+1, n-i, -tau, a( kl+i, i ), 1, work, 1,
341  $ a( kl+i, i+1 ), lda )
342  a( kl+i, i ) = -wa
343  END IF
344  END IF
345 *
346  IF (i .LE. n) THEN
347  DO 50 j = kl + i + 1, m
348  a( j, i ) = zero
349  50 CONTINUE
350  END IF
351 *
352  IF (i .LE. m) THEN
353  DO 60 j = ku + i + 1, n
354  a( i, j ) = zero
355  60 CONTINUE
356  END IF
357  70 CONTINUE
358  RETURN
359 *
360 * End of SLAGGE
361 *
subroutine slarnv(IDIST, ISEED, N, X)
SLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition: slarnv.f:97
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79
real(wp) function snrm2(n, x, incx)
SNRM2
Definition: snrm2.f90:89
subroutine sger(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
SGER
Definition: sger.f:130
subroutine sgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
SGEMV
Definition: sgemv.f:156
Here is the call graph for this function:
Here is the caller graph for this function: