LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ dsfrk()

subroutine dsfrk ( character  transr,
character  uplo,
character  trans,
integer  n,
integer  k,
double precision  alpha,
double precision, dimension( lda, * )  a,
integer  lda,
double precision  beta,
double precision, dimension( * )  c 
)

DSFRK performs a symmetric rank-k operation for matrix in RFP format.

Download DSFRK + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 Level 3 BLAS like routine for C in RFP Format.

 DSFRK performs one of the symmetric rank--k operations

    C := alpha*A*A**T + beta*C,

 or

    C := alpha*A**T*A + beta*C,

 where alpha and beta are real scalars, C is an n--by--n symmetric
 matrix and A is an n--by--k matrix in the first case and a k--by--n
 matrix in the second case.
Parameters
[in]TRANSR
          TRANSR is CHARACTER*1
          = 'N':  The Normal Form of RFP A is stored;
          = 'T':  The Transpose Form of RFP A is stored.
[in]UPLO
          UPLO is CHARACTER*1
           On  entry, UPLO specifies whether the upper or lower
           triangular part of the array C is to be referenced as
           follows:

              UPLO = 'U' or 'u'   Only the upper triangular part of C
                                  is to be referenced.

              UPLO = 'L' or 'l'   Only the lower triangular part of C
                                  is to be referenced.

           Unchanged on exit.
[in]TRANS
          TRANS is CHARACTER*1
           On entry, TRANS specifies the operation to be performed as
           follows:

              TRANS = 'N' or 'n'   C := alpha*A*A**T + beta*C.

              TRANS = 'T' or 't'   C := alpha*A**T*A + beta*C.

           Unchanged on exit.
[in]N
          N is INTEGER
           On entry, N specifies the order of the matrix C. N must be
           at least zero.
           Unchanged on exit.
[in]K
          K is INTEGER
           On entry with TRANS = 'N' or 'n', K specifies the number
           of  columns of the matrix A, and on entry with TRANS = 'T'
           or 't', K specifies the number of rows of the matrix A. K
           must be at least zero.
           Unchanged on exit.
[in]ALPHA
          ALPHA is DOUBLE PRECISION
           On entry, ALPHA specifies the scalar alpha.
           Unchanged on exit.
[in]A
          A is DOUBLE PRECISION array, dimension (LDA,ka)
           where KA
           is K  when TRANS = 'N' or 'n', and is N otherwise. Before
           entry with TRANS = 'N' or 'n', the leading N--by--K part of
           the array A must contain the matrix A, otherwise the leading
           K--by--N part of the array A must contain the matrix A.
           Unchanged on exit.
[in]LDA
          LDA is INTEGER
           On entry, LDA specifies the first dimension of A as declared
           in  the  calling  (sub)  program.   When  TRANS = 'N' or 'n'
           then  LDA must be at least  max( 1, n ), otherwise  LDA must
           be at least  max( 1, k ).
           Unchanged on exit.
[in]BETA
          BETA is DOUBLE PRECISION
           On entry, BETA specifies the scalar beta.
           Unchanged on exit.
[in,out]C
          C is DOUBLE PRECISION array, dimension (NT)
           NT = N*(N+1)/2. On entry, the symmetric matrix C in RFP
           Format. RFP Format is described by TRANSR, UPLO and N.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 164 of file dsfrk.f.

166*
167* -- LAPACK computational routine --
168* -- LAPACK is a software package provided by Univ. of Tennessee, --
169* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
170*
171* .. Scalar Arguments ..
172 DOUBLE PRECISION ALPHA, BETA
173 INTEGER K, LDA, N
174 CHARACTER TRANS, TRANSR, UPLO
175* ..
176* .. Array Arguments ..
177 DOUBLE PRECISION A( LDA, * ), C( * )
178* ..
179*
180* =====================================================================
181*
182* ..
183* .. Parameters ..
184 DOUBLE PRECISION ONE, ZERO
185 parameter( one = 1.0d+0, zero = 0.0d+0 )
186* ..
187* .. Local Scalars ..
188 LOGICAL LOWER, NORMALTRANSR, NISODD, NOTRANS
189 INTEGER INFO, NROWA, J, NK, N1, N2
190* ..
191* .. External Functions ..
192 LOGICAL LSAME
193 EXTERNAL lsame
194* ..
195* .. External Subroutines ..
196 EXTERNAL xerbla, dgemm, dsyrk
197* ..
198* .. Intrinsic Functions ..
199 INTRINSIC max
200* ..
201* .. Executable Statements ..
202*
203* Test the input parameters.
204*
205 info = 0
206 normaltransr = lsame( transr, 'N' )
207 lower = lsame( uplo, 'L' )
208 notrans = lsame( trans, 'N' )
209*
210 IF( notrans ) THEN
211 nrowa = n
212 ELSE
213 nrowa = k
214 END IF
215*
216 IF( .NOT.normaltransr .AND. .NOT.lsame( transr, 'T' ) ) THEN
217 info = -1
218 ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo, 'U' ) ) THEN
219 info = -2
220 ELSE IF( .NOT.notrans .AND. .NOT.lsame( trans, 'T' ) ) THEN
221 info = -3
222 ELSE IF( n.LT.0 ) THEN
223 info = -4
224 ELSE IF( k.LT.0 ) THEN
225 info = -5
226 ELSE IF( lda.LT.max( 1, nrowa ) ) THEN
227 info = -8
228 END IF
229 IF( info.NE.0 ) THEN
230 CALL xerbla( 'DSFRK ', -info )
231 RETURN
232 END IF
233*
234* Quick return if possible.
235*
236* The quick return case: ((ALPHA.EQ.0).AND.(BETA.NE.ZERO)) is not
237* done (it is in DSYRK for example) and left in the general case.
238*
239 IF( ( n.EQ.0 ) .OR. ( ( ( alpha.EQ.zero ) .OR. ( k.EQ.0 ) ) .AND.
240 $ ( beta.EQ.one ) ) )RETURN
241*
242 IF( ( alpha.EQ.zero ) .AND. ( beta.EQ.zero ) ) THEN
243 DO j = 1, ( ( n*( n+1 ) ) / 2 )
244 c( j ) = zero
245 END DO
246 RETURN
247 END IF
248*
249* C is N-by-N.
250* If N is odd, set NISODD = .TRUE., and N1 and N2.
251* If N is even, NISODD = .FALSE., and NK.
252*
253 IF( mod( n, 2 ).EQ.0 ) THEN
254 nisodd = .false.
255 nk = n / 2
256 ELSE
257 nisodd = .true.
258 IF( lower ) THEN
259 n2 = n / 2
260 n1 = n - n2
261 ELSE
262 n1 = n / 2
263 n2 = n - n1
264 END IF
265 END IF
266*
267 IF( nisodd ) THEN
268*
269* N is odd
270*
271 IF( normaltransr ) THEN
272*
273* N is odd and TRANSR = 'N'
274*
275 IF( lower ) THEN
276*
277* N is odd, TRANSR = 'N', and UPLO = 'L'
278*
279 IF( notrans ) THEN
280*
281* N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
282*
283 CALL dsyrk( 'L', 'N', n1, k, alpha, a( 1, 1 ), lda,
284 $ beta, c( 1 ), n )
285 CALL dsyrk( 'U', 'N', n2, k, alpha, a( n1+1, 1 ), lda,
286 $ beta, c( n+1 ), n )
287 CALL dgemm( 'N', 'T', n2, n1, k, alpha, a( n1+1, 1 ),
288 $ lda, a( 1, 1 ), lda, beta, c( n1+1 ), n )
289*
290 ELSE
291*
292* N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'T'
293*
294 CALL dsyrk( 'L', 'T', n1, k, alpha, a( 1, 1 ), lda,
295 $ beta, c( 1 ), n )
296 CALL dsyrk( 'U', 'T', n2, k, alpha, a( 1, n1+1 ), lda,
297 $ beta, c( n+1 ), n )
298 CALL dgemm( 'T', 'N', n2, n1, k, alpha, a( 1, n1+1 ),
299 $ lda, a( 1, 1 ), lda, beta, c( n1+1 ), n )
300*
301 END IF
302*
303 ELSE
304*
305* N is odd, TRANSR = 'N', and UPLO = 'U'
306*
307 IF( notrans ) THEN
308*
309* N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
310*
311 CALL dsyrk( 'L', 'N', n1, k, alpha, a( 1, 1 ), lda,
312 $ beta, c( n2+1 ), n )
313 CALL dsyrk( 'U', 'N', n2, k, alpha, a( n2, 1 ), lda,
314 $ beta, c( n1+1 ), n )
315 CALL dgemm( 'N', 'T', n1, n2, k, alpha, a( 1, 1 ),
316 $ lda, a( n2, 1 ), lda, beta, c( 1 ), n )
317*
318 ELSE
319*
320* N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'T'
321*
322 CALL dsyrk( 'L', 'T', n1, k, alpha, a( 1, 1 ), lda,
323 $ beta, c( n2+1 ), n )
324 CALL dsyrk( 'U', 'T', n2, k, alpha, a( 1, n2 ), lda,
325 $ beta, c( n1+1 ), n )
326 CALL dgemm( 'T', 'N', n1, n2, k, alpha, a( 1, 1 ),
327 $ lda, a( 1, n2 ), lda, beta, c( 1 ), n )
328*
329 END IF
330*
331 END IF
332*
333 ELSE
334*
335* N is odd, and TRANSR = 'T'
336*
337 IF( lower ) THEN
338*
339* N is odd, TRANSR = 'T', and UPLO = 'L'
340*
341 IF( notrans ) THEN
342*
343* N is odd, TRANSR = 'T', UPLO = 'L', and TRANS = 'N'
344*
345 CALL dsyrk( 'U', 'N', n1, k, alpha, a( 1, 1 ), lda,
346 $ beta, c( 1 ), n1 )
347 CALL dsyrk( 'L', 'N', n2, k, alpha, a( n1+1, 1 ), lda,
348 $ beta, c( 2 ), n1 )
349 CALL dgemm( 'N', 'T', n1, n2, k, alpha, a( 1, 1 ),
350 $ lda, a( n1+1, 1 ), lda, beta,
351 $ c( n1*n1+1 ), n1 )
352*
353 ELSE
354*
355* N is odd, TRANSR = 'T', UPLO = 'L', and TRANS = 'T'
356*
357 CALL dsyrk( 'U', 'T', n1, k, alpha, a( 1, 1 ), lda,
358 $ beta, c( 1 ), n1 )
359 CALL dsyrk( 'L', 'T', n2, k, alpha, a( 1, n1+1 ), lda,
360 $ beta, c( 2 ), n1 )
361 CALL dgemm( 'T', 'N', n1, n2, k, alpha, a( 1, 1 ),
362 $ lda, a( 1, n1+1 ), lda, beta,
363 $ c( n1*n1+1 ), n1 )
364*
365 END IF
366*
367 ELSE
368*
369* N is odd, TRANSR = 'T', and UPLO = 'U'
370*
371 IF( notrans ) THEN
372*
373* N is odd, TRANSR = 'T', UPLO = 'U', and TRANS = 'N'
374*
375 CALL dsyrk( 'U', 'N', n1, k, alpha, a( 1, 1 ), lda,
376 $ beta, c( n2*n2+1 ), n2 )
377 CALL dsyrk( 'L', 'N', n2, k, alpha, a( n1+1, 1 ), lda,
378 $ beta, c( n1*n2+1 ), n2 )
379 CALL dgemm( 'N', 'T', n2, n1, k, alpha, a( n1+1, 1 ),
380 $ lda, a( 1, 1 ), lda, beta, c( 1 ), n2 )
381*
382 ELSE
383*
384* N is odd, TRANSR = 'T', UPLO = 'U', and TRANS = 'T'
385*
386 CALL dsyrk( 'U', 'T', n1, k, alpha, a( 1, 1 ), lda,
387 $ beta, c( n2*n2+1 ), n2 )
388 CALL dsyrk( 'L', 'T', n2, k, alpha, a( 1, n1+1 ), lda,
389 $ beta, c( n1*n2+1 ), n2 )
390 CALL dgemm( 'T', 'N', n2, n1, k, alpha, a( 1, n1+1 ),
391 $ lda, a( 1, 1 ), lda, beta, c( 1 ), n2 )
392*
393 END IF
394*
395 END IF
396*
397 END IF
398*
399 ELSE
400*
401* N is even
402*
403 IF( normaltransr ) THEN
404*
405* N is even and TRANSR = 'N'
406*
407 IF( lower ) THEN
408*
409* N is even, TRANSR = 'N', and UPLO = 'L'
410*
411 IF( notrans ) THEN
412*
413* N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
414*
415 CALL dsyrk( 'L', 'N', nk, k, alpha, a( 1, 1 ), lda,
416 $ beta, c( 2 ), n+1 )
417 CALL dsyrk( 'U', 'N', nk, k, alpha, a( nk+1, 1 ), lda,
418 $ beta, c( 1 ), n+1 )
419 CALL dgemm( 'N', 'T', nk, nk, k, alpha, a( nk+1, 1 ),
420 $ lda, a( 1, 1 ), lda, beta, c( nk+2 ),
421 $ n+1 )
422*
423 ELSE
424*
425* N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'T'
426*
427 CALL dsyrk( 'L', 'T', nk, k, alpha, a( 1, 1 ), lda,
428 $ beta, c( 2 ), n+1 )
429 CALL dsyrk( 'U', 'T', nk, k, alpha, a( 1, nk+1 ), lda,
430 $ beta, c( 1 ), n+1 )
431 CALL dgemm( 'T', 'N', nk, nk, k, alpha, a( 1, nk+1 ),
432 $ lda, a( 1, 1 ), lda, beta, c( nk+2 ),
433 $ n+1 )
434*
435 END IF
436*
437 ELSE
438*
439* N is even, TRANSR = 'N', and UPLO = 'U'
440*
441 IF( notrans ) THEN
442*
443* N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
444*
445 CALL dsyrk( 'L', 'N', nk, k, alpha, a( 1, 1 ), lda,
446 $ beta, c( nk+2 ), n+1 )
447 CALL dsyrk( 'U', 'N', nk, k, alpha, a( nk+1, 1 ), lda,
448 $ beta, c( nk+1 ), n+1 )
449 CALL dgemm( 'N', 'T', nk, nk, k, alpha, a( 1, 1 ),
450 $ lda, a( nk+1, 1 ), lda, beta, c( 1 ),
451 $ n+1 )
452*
453 ELSE
454*
455* N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'T'
456*
457 CALL dsyrk( 'L', 'T', nk, k, alpha, a( 1, 1 ), lda,
458 $ beta, c( nk+2 ), n+1 )
459 CALL dsyrk( 'U', 'T', nk, k, alpha, a( 1, nk+1 ), lda,
460 $ beta, c( nk+1 ), n+1 )
461 CALL dgemm( 'T', 'N', nk, nk, k, alpha, a( 1, 1 ),
462 $ lda, a( 1, nk+1 ), lda, beta, c( 1 ),
463 $ n+1 )
464*
465 END IF
466*
467 END IF
468*
469 ELSE
470*
471* N is even, and TRANSR = 'T'
472*
473 IF( lower ) THEN
474*
475* N is even, TRANSR = 'T', and UPLO = 'L'
476*
477 IF( notrans ) THEN
478*
479* N is even, TRANSR = 'T', UPLO = 'L', and TRANS = 'N'
480*
481 CALL dsyrk( 'U', 'N', nk, k, alpha, a( 1, 1 ), lda,
482 $ beta, c( nk+1 ), nk )
483 CALL dsyrk( 'L', 'N', nk, k, alpha, a( nk+1, 1 ), lda,
484 $ beta, c( 1 ), nk )
485 CALL dgemm( 'N', 'T', nk, nk, k, alpha, a( 1, 1 ),
486 $ lda, a( nk+1, 1 ), lda, beta,
487 $ c( ( ( nk+1 )*nk )+1 ), nk )
488*
489 ELSE
490*
491* N is even, TRANSR = 'T', UPLO = 'L', and TRANS = 'T'
492*
493 CALL dsyrk( 'U', 'T', nk, k, alpha, a( 1, 1 ), lda,
494 $ beta, c( nk+1 ), nk )
495 CALL dsyrk( 'L', 'T', nk, k, alpha, a( 1, nk+1 ), lda,
496 $ beta, c( 1 ), nk )
497 CALL dgemm( 'T', 'N', nk, nk, k, alpha, a( 1, 1 ),
498 $ lda, a( 1, nk+1 ), lda, beta,
499 $ c( ( ( nk+1 )*nk )+1 ), nk )
500*
501 END IF
502*
503 ELSE
504*
505* N is even, TRANSR = 'T', and UPLO = 'U'
506*
507 IF( notrans ) THEN
508*
509* N is even, TRANSR = 'T', UPLO = 'U', and TRANS = 'N'
510*
511 CALL dsyrk( 'U', 'N', nk, k, alpha, a( 1, 1 ), lda,
512 $ beta, c( nk*( nk+1 )+1 ), nk )
513 CALL dsyrk( 'L', 'N', nk, k, alpha, a( nk+1, 1 ), lda,
514 $ beta, c( nk*nk+1 ), nk )
515 CALL dgemm( 'N', 'T', nk, nk, k, alpha, a( nk+1, 1 ),
516 $ lda, a( 1, 1 ), lda, beta, c( 1 ), nk )
517*
518 ELSE
519*
520* N is even, TRANSR = 'T', UPLO = 'U', and TRANS = 'T'
521*
522 CALL dsyrk( 'U', 'T', nk, k, alpha, a( 1, 1 ), lda,
523 $ beta, c( nk*( nk+1 )+1 ), nk )
524 CALL dsyrk( 'L', 'T', nk, k, alpha, a( 1, nk+1 ), lda,
525 $ beta, c( nk*nk+1 ), nk )
526 CALL dgemm( 'T', 'N', nk, nk, k, alpha, a( 1, nk+1 ),
527 $ lda, a( 1, 1 ), lda, beta, c( 1 ), nk )
528*
529 END IF
530*
531 END IF
532*
533 END IF
534*
535 END IF
536*
537 RETURN
538*
539* End of DSFRK
540*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
dgemm
subroutine dgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
DGEMM
Definition dgemm.f:188
dsyrk
subroutine dsyrk(uplo, trans, n, k, alpha, a, lda, beta, c, ldc)
DSYRK
Definition dsyrk.f:169
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the call graph for this function:
Here is the caller graph for this function: