SRC/dsfrk.f

   1 *> \brief \b DSFRK performs a symmetric rank-k operation for matrix in RFP format.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download DSFRK + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsfrk.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsfrk.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsfrk.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE DSFRK( TRANSR, UPLO, TRANS, N, K, ALPHA, A, LDA, BETA,
  22 *                         C )
  23 *
  24 *       .. Scalar Arguments ..
  25 *       DOUBLE PRECISION   ALPHA, BETA
  26 *       INTEGER            K, LDA, N
  27 *       CHARACTER          TRANS, TRANSR, UPLO
  28 *       ..
  29 *       .. Array Arguments ..
  30 *       DOUBLE PRECISION   A( LDA, * ), C( * )
  31 *       ..
  32 *
  33 *
  34 *> \par Purpose:
  35 *  =============
  36 *>
  37 *> \verbatim
  38 *>
  39 *> Level 3 BLAS like routine for C in RFP Format.
  40 *>
  41 *> DSFRK performs one of the symmetric rank--k operations
  42 *>
  43 *>    C := alpha*A*A**T + beta*C,
  44 *>
  45 *> or
  46 *>
  47 *>    C := alpha*A**T*A + beta*C,
  48 *>
  49 *> where alpha and beta are real scalars, C is an n--by--n symmetric
  50 *> matrix and A is an n--by--k matrix in the first case and a k--by--n
  51 *> matrix in the second case.
  52 *> \endverbatim
  53 *
  54 *  Arguments:
  55 *  ==========
  56 *
  57 *> \param[in] TRANSR
  58 *> \verbatim
  59 *>          TRANSR is CHARACTER*1
  60 *>          = 'N':  The Normal Form of RFP A is stored;
  61 *>          = 'T':  The Transpose Form of RFP A is stored.
  62 *> \endverbatim
  63 *>
  64 *> \param[in] UPLO
  65 *> \verbatim
  66 *>          UPLO is CHARACTER*1
  67 *>           On  entry, UPLO specifies whether the upper or lower
  68 *>           triangular part of the array C is to be referenced as
  69 *>           follows:
  70 *>
  71 *>              UPLO = 'U' or 'u'   Only the upper triangular part of C
  72 *>                                  is to be referenced.
  73 *>
  74 *>              UPLO = 'L' or 'l'   Only the lower triangular part of C
  75 *>                                  is to be referenced.
  76 *>
  77 *>           Unchanged on exit.
  78 *> \endverbatim
  79 *>
  80 *> \param[in] TRANS
  81 *> \verbatim
  82 *>          TRANS is CHARACTER*1
  83 *>           On entry, TRANS specifies the operation to be performed as
  84 *>           follows:
  85 *>
  86 *>              TRANS = 'N' or 'n'   C := alpha*A*A**T + beta*C.
  87 *>
  88 *>              TRANS = 'T' or 't'   C := alpha*A**T*A + beta*C.
  89 *>
  90 *>           Unchanged on exit.
  91 *> \endverbatim
  92 *>
  93 *> \param[in] N
  94 *> \verbatim
  95 *>          N is INTEGER
  96 *>           On entry, N specifies the order of the matrix C. N must be
  97 *>           at least zero.
  98 *>           Unchanged on exit.
  99 *> \endverbatim
 100 *>
 101 *> \param[in] K
 102 *> \verbatim
 103 *>          K is INTEGER
 104 *>           On entry with TRANS = 'N' or 'n', K specifies the number
 105 *>           of  columns of the matrix A, and on entry with TRANS = 'T'
 106 *>           or 't', K specifies the number of rows of the matrix A. K
 107 *>           must be at least zero.
 108 *>           Unchanged on exit.
 109 *> \endverbatim
 110 *>
 111 *> \param[in] ALPHA
 112 *> \verbatim
 113 *>          ALPHA is DOUBLE PRECISION
 114 *>           On entry, ALPHA specifies the scalar alpha.
 115 *>           Unchanged on exit.
 116 *> \endverbatim
 117 *>
 118 *> \param[in] A
 119 *> \verbatim
 120 *>          A is DOUBLE PRECISION array, dimension (LDA,ka)
 121 *>           where KA
 122 *>           is K  when TRANS = 'N' or 'n', and is N otherwise. Before
 123 *>           entry with TRANS = 'N' or 'n', the leading N--by--K part of
 124 *>           the array A must contain the matrix A, otherwise the leading
 125 *>           K--by--N part of the array A must contain the matrix A.
 126 *>           Unchanged on exit.
 127 *> \endverbatim
 128 *>
 129 *> \param[in] LDA
 130 *> \verbatim
 131 *>          LDA is INTEGER
 132 *>           On entry, LDA specifies the first dimension of A as declared
 133 *>           in  the  calling  (sub)  program.   When  TRANS = 'N' or 'n'
 134 *>           then  LDA must be at least  max( 1, n ), otherwise  LDA must
 135 *>           be at least  max( 1, k ).
 136 *>           Unchanged on exit.
 137 *> \endverbatim
 138 *>
 139 *> \param[in] BETA
 140 *> \verbatim
 141 *>          BETA is DOUBLE PRECISION
 142 *>           On entry, BETA specifies the scalar beta.
 143 *>           Unchanged on exit.
 144 *> \endverbatim
 145 *>
 146 *> \param[in,out] C
 147 *> \verbatim
 148 *>          C is DOUBLE PRECISION array, dimension (NT)
 149 *>           NT = N*(N+1)/2. On entry, the symmetric matrix C in RFP
 150 *>           Format. RFP Format is described by TRANSR, UPLO and N.
 151 *> \endverbatim
 152 *
 153 *  Authors:
 154 *  ========
 155 *
 156 *> \author Univ. of Tennessee
 157 *> \author Univ. of California Berkeley
 158 *> \author Univ. of Colorado Denver
 159 *> \author NAG Ltd.
 160 *
 161 *> \date September 2012
 162 *
 163 *> \ingroup doubleOTHERcomputational
 164 *
 165 *  =====================================================================
 166       SUBROUTINE DSFRK( TRANSR, UPLO, TRANS, N, K, ALPHA, A, LDA, BETA,
 167      $                  C )
 168 *
 169 *  -- LAPACK computational routine (version 3.4.2) --
 170 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 171 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 172 *     September 2012
 173 *
 174 *     .. Scalar Arguments ..
 175       DOUBLE PRECISION   ALPHA, BETA
 176       INTEGER            K, LDA, N
 177       CHARACTER          TRANS, TRANSR, UPLO
 178 *     ..
 179 *     .. Array Arguments ..
 180       DOUBLE PRECISION   A( LDA, * ), C( * )
 181 *     ..
 182 *
 183 *  =====================================================================
 184 *
 185 *     ..
 186 *     .. Parameters ..
 187       DOUBLE PRECISION   ONE, ZERO
 188       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
 189 *     ..
 190 *     .. Local Scalars ..
 191       LOGICAL            LOWER, NORMALTRANSR, NISODD, NOTRANS
 192       INTEGER            INFO, NROWA, J, NK, N1, N2
 193 *     ..
 194 *     .. External Functions ..
 195       LOGICAL            LSAME
 196       EXTERNAL           LSAME
 197 *     ..
 198 *     .. External Subroutines ..
 199       EXTERNAL           XERBLA, DGEMM, DSYRK
 200 *     ..
 201 *     .. Intrinsic Functions ..
 202       INTRINSIC          MAX
 203 *     ..
 204 *     .. Executable Statements ..
 205 *
 206 *     Test the input parameters.
 207 *
 208       INFO = 0
 209       NORMALTRANSR = LSAME( TRANSR, 'N' )
 210       LOWER = LSAME( UPLO, 'L' )
 211       NOTRANS = LSAME( TRANS, 'N' )
 212 *
 213       IF( NOTRANS ) THEN
 214          NROWA = N
 215       ELSE
 216          NROWA = K
 217       END IF
 218 *
 219       IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN
 220          INFO = -1
 221       ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
 222          INFO = -2
 223       ELSE IF( .NOT.NOTRANS .AND. .NOT.LSAME( TRANS, 'T' ) ) THEN
 224          INFO = -3
 225       ELSE IF( N.LT.0 ) THEN
 226          INFO = -4
 227       ELSE IF( K.LT.0 ) THEN
 228          INFO = -5
 229       ELSE IF( LDA.LT.MAX( 1, NROWA ) ) THEN
 230          INFO = -8
 231       END IF
 232       IF( INFO.NE.0 ) THEN
 233          CALL XERBLA( 'DSFRK ', -INFO )
 234          RETURN
 235       END IF
 236 *
 237 *     Quick return if possible.
 238 *
 239 *     The quick return case: ((ALPHA.EQ.0).AND.(BETA.NE.ZERO)) is not
 240 *     done (it is in DSYRK for example) and left in the general case.
 241 *
 242       IF( ( N.EQ.0 ) .OR. ( ( ( ALPHA.EQ.ZERO ) .OR. ( K.EQ.0 ) ) .AND.
 243      $    ( BETA.EQ.ONE ) ) )RETURN
 244 *
 245       IF( ( ALPHA.EQ.ZERO ) .AND. ( BETA.EQ.ZERO ) ) THEN
 246          DO J = 1, ( ( N*( N+1 ) ) / 2 )
 247             C( J ) = ZERO
 248          END DO
 249          RETURN
 250       END IF
 251 *
 252 *     C is N-by-N.
 253 *     If N is odd, set NISODD = .TRUE., and N1 and N2.
 254 *     If N is even, NISODD = .FALSE., and NK.
 255 *
 256       IF( MOD( N, 2 ).EQ.0 ) THEN
 257          NISODD = .FALSE.
 258          NK = N / 2
 259       ELSE
 260          NISODD = .TRUE.
 261          IF( LOWER ) THEN
 262             N2 = N / 2
 263             N1 = N - N2
 264          ELSE
 265             N1 = N / 2
 266             N2 = N - N1
 267          END IF
 268       END IF
 269 *
 270       IF( NISODD ) THEN
 271 *
 272 *        N is odd
 273 *
 274          IF( NORMALTRANSR ) THEN
 275 *
 276 *           N is odd and TRANSR = 'N'
 277 *
 278             IF( LOWER ) THEN
 279 *
 280 *              N is odd, TRANSR = 'N', and UPLO = 'L'
 281 *
 282                IF( NOTRANS ) THEN
 283 *
 284 *                 N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
 285 *
 286                   CALL DSYRK( 'L', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
 287      $                        BETA, C( 1 ), N )
 288                   CALL DSYRK( 'U', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
 289      $                        BETA, C( N+1 ), N )
 290                   CALL DGEMM( 'N', 'T', N2, N1, K, ALPHA, A( N1+1, 1 ),
 291      $                        LDA, A( 1, 1 ), LDA, BETA, C( N1+1 ), N )
 292 *
 293                ELSE
 294 *
 295 *                 N is odd, TRANSR = 'N', UPLO = 'L', and TRANS = 'T'
 296 *
 297                   CALL DSYRK( 'L', 'T', N1, K, ALPHA, A( 1, 1 ), LDA,
 298      $                        BETA, C( 1 ), N )
 299                   CALL DSYRK( 'U', 'T', N2, K, ALPHA, A( 1, N1+1 ), LDA,
 300      $                        BETA, C( N+1 ), N )
 301                   CALL DGEMM( 'T', 'N', N2, N1, K, ALPHA, A( 1, N1+1 ),
 302      $                        LDA, A( 1, 1 ), LDA, BETA, C( N1+1 ), N )
 303 *
 304                END IF
 305 *
 306             ELSE
 307 *
 308 *              N is odd, TRANSR = 'N', and UPLO = 'U'
 309 *
 310                IF( NOTRANS ) THEN
 311 *
 312 *                 N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
 313 *
 314                   CALL DSYRK( 'L', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
 315      $                        BETA, C( N2+1 ), N )
 316                   CALL DSYRK( 'U', 'N', N2, K, ALPHA, A( N2, 1 ), LDA,
 317      $                        BETA, C( N1+1 ), N )
 318                   CALL DGEMM( 'N', 'T', N1, N2, K, ALPHA, A( 1, 1 ),
 319      $                        LDA, A( N2, 1 ), LDA, BETA, C( 1 ), N )
 320 *
 321                ELSE
 322 *
 323 *                 N is odd, TRANSR = 'N', UPLO = 'U', and TRANS = 'T'
 324 *
 325                   CALL DSYRK( 'L', 'T', N1, K, ALPHA, A( 1, 1 ), LDA,
 326      $                        BETA, C( N2+1 ), N )
 327                   CALL DSYRK( 'U', 'T', N2, K, ALPHA, A( 1, N2 ), LDA,
 328      $                        BETA, C( N1+1 ), N )
 329                   CALL DGEMM( 'T', 'N', N1, N2, K, ALPHA, A( 1, 1 ),
 330      $                        LDA, A( 1, N2 ), LDA, BETA, C( 1 ), N )
 331 *
 332                END IF
 333 *
 334             END IF
 335 *
 336          ELSE
 337 *
 338 *           N is odd, and TRANSR = 'T'
 339 *
 340             IF( LOWER ) THEN
 341 *
 342 *              N is odd, TRANSR = 'T', and UPLO = 'L'
 343 *
 344                IF( NOTRANS ) THEN
 345 *
 346 *                 N is odd, TRANSR = 'T', UPLO = 'L', and TRANS = 'N'
 347 *
 348                   CALL DSYRK( 'U', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
 349      $                        BETA, C( 1 ), N1 )
 350                   CALL DSYRK( 'L', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
 351      $                        BETA, C( 2 ), N1 )
 352                   CALL DGEMM( 'N', 'T', N1, N2, K, ALPHA, A( 1, 1 ),
 353      $                        LDA, A( N1+1, 1 ), LDA, BETA,
 354      $                        C( N1*N1+1 ), N1 )
 355 *
 356                ELSE
 357 *
 358 *                 N is odd, TRANSR = 'T', UPLO = 'L', and TRANS = 'T'
 359 *
 360                   CALL DSYRK( 'U', 'T', N1, K, ALPHA, A( 1, 1 ), LDA,
 361      $                        BETA, C( 1 ), N1 )
 362                   CALL DSYRK( 'L', 'T', N2, K, ALPHA, A( 1, N1+1 ), LDA,
 363      $                        BETA, C( 2 ), N1 )
 364                   CALL DGEMM( 'T', 'N', N1, N2, K, ALPHA, A( 1, 1 ),
 365      $                        LDA, A( 1, N1+1 ), LDA, BETA,
 366      $                        C( N1*N1+1 ), N1 )
 367 *
 368                END IF
 369 *
 370             ELSE
 371 *
 372 *              N is odd, TRANSR = 'T', and UPLO = 'U'
 373 *
 374                IF( NOTRANS ) THEN
 375 *
 376 *                 N is odd, TRANSR = 'T', UPLO = 'U', and TRANS = 'N'
 377 *
 378                   CALL DSYRK( 'U', 'N', N1, K, ALPHA, A( 1, 1 ), LDA,
 379      $                        BETA, C( N2*N2+1 ), N2 )
 380                   CALL DSYRK( 'L', 'N', N2, K, ALPHA, A( N1+1, 1 ), LDA,
 381      $                        BETA, C( N1*N2+1 ), N2 )
 382                   CALL DGEMM( 'N', 'T', N2, N1, K, ALPHA, A( N1+1, 1 ),
 383      $                        LDA, A( 1, 1 ), LDA, BETA, C( 1 ), N2 )
 384 *
 385                ELSE
 386 *
 387 *                 N is odd, TRANSR = 'T', UPLO = 'U', and TRANS = 'T'
 388 *
 389                   CALL DSYRK( 'U', 'T', N1, K, ALPHA, A( 1, 1 ), LDA,
 390      $                        BETA, C( N2*N2+1 ), N2 )
 391                   CALL DSYRK( 'L', 'T', N2, K, ALPHA, A( 1, N1+1 ), LDA,
 392      $                        BETA, C( N1*N2+1 ), N2 )
 393                   CALL DGEMM( 'T', 'N', N2, N1, K, ALPHA, A( 1, N1+1 ),
 394      $                        LDA, A( 1, 1 ), LDA, BETA, C( 1 ), N2 )
 395 *
 396                END IF
 397 *
 398             END IF
 399 *
 400          END IF
 401 *
 402       ELSE
 403 *
 404 *        N is even
 405 *
 406          IF( NORMALTRANSR ) THEN
 407 *
 408 *           N is even and TRANSR = 'N'
 409 *
 410             IF( LOWER ) THEN
 411 *
 412 *              N is even, TRANSR = 'N', and UPLO = 'L'
 413 *
 414                IF( NOTRANS ) THEN
 415 *
 416 *                 N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'N'
 417 *
 418                   CALL DSYRK( 'L', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
 419      $                        BETA, C( 2 ), N+1 )
 420                   CALL DSYRK( 'U', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
 421      $                        BETA, C( 1 ), N+1 )
 422                   CALL DGEMM( 'N', 'T', NK, NK, K, ALPHA, A( NK+1, 1 ),
 423      $                        LDA, A( 1, 1 ), LDA, BETA, C( NK+2 ),
 424      $                        N+1 )
 425 *
 426                ELSE
 427 *
 428 *                 N is even, TRANSR = 'N', UPLO = 'L', and TRANS = 'T'
 429 *
 430                   CALL DSYRK( 'L', 'T', NK, K, ALPHA, A( 1, 1 ), LDA,
 431      $                        BETA, C( 2 ), N+1 )
 432                   CALL DSYRK( 'U', 'T', NK, K, ALPHA, A( 1, NK+1 ), LDA,
 433      $                        BETA, C( 1 ), N+1 )
 434                   CALL DGEMM( 'T', 'N', NK, NK, K, ALPHA, A( 1, NK+1 ),
 435      $                        LDA, A( 1, 1 ), LDA, BETA, C( NK+2 ),
 436      $                        N+1 )
 437 *
 438                END IF
 439 *
 440             ELSE
 441 *
 442 *              N is even, TRANSR = 'N', and UPLO = 'U'
 443 *
 444                IF( NOTRANS ) THEN
 445 *
 446 *                 N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'N'
 447 *
 448                   CALL DSYRK( 'L', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
 449      $                        BETA, C( NK+2 ), N+1 )
 450                   CALL DSYRK( 'U', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
 451      $                        BETA, C( NK+1 ), N+1 )
 452                   CALL DGEMM( 'N', 'T', NK, NK, K, ALPHA, A( 1, 1 ),
 453      $                        LDA, A( NK+1, 1 ), LDA, BETA, C( 1 ),
 454      $                        N+1 )
 455 *
 456                ELSE
 457 *
 458 *                 N is even, TRANSR = 'N', UPLO = 'U', and TRANS = 'T'
 459 *
 460                   CALL DSYRK( 'L', 'T', NK, K, ALPHA, A( 1, 1 ), LDA,
 461      $                        BETA, C( NK+2 ), N+1 )
 462                   CALL DSYRK( 'U', 'T', NK, K, ALPHA, A( 1, NK+1 ), LDA,
 463      $                        BETA, C( NK+1 ), N+1 )
 464                   CALL DGEMM( 'T', 'N', NK, NK, K, ALPHA, A( 1, 1 ),
 465      $                        LDA, A( 1, NK+1 ), LDA, BETA, C( 1 ),
 466      $                        N+1 )
 467 *
 468                END IF
 469 *
 470             END IF
 471 *
 472          ELSE
 473 *
 474 *           N is even, and TRANSR = 'T'
 475 *
 476             IF( LOWER ) THEN
 477 *
 478 *              N is even, TRANSR = 'T', and UPLO = 'L'
 479 *
 480                IF( NOTRANS ) THEN
 481 *
 482 *                 N is even, TRANSR = 'T', UPLO = 'L', and TRANS = 'N'
 483 *
 484                   CALL DSYRK( 'U', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
 485      $                        BETA, C( NK+1 ), NK )
 486                   CALL DSYRK( 'L', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
 487      $                        BETA, C( 1 ), NK )
 488                   CALL DGEMM( 'N', 'T', NK, NK, K, ALPHA, A( 1, 1 ),
 489      $                        LDA, A( NK+1, 1 ), LDA, BETA,
 490      $                        C( ( ( NK+1 )*NK )+1 ), NK )
 491 *
 492                ELSE
 493 *
 494 *                 N is even, TRANSR = 'T', UPLO = 'L', and TRANS = 'T'
 495 *
 496                   CALL DSYRK( 'U', 'T', NK, K, ALPHA, A( 1, 1 ), LDA,
 497      $                        BETA, C( NK+1 ), NK )
 498                   CALL DSYRK( 'L', 'T', NK, K, ALPHA, A( 1, NK+1 ), LDA,
 499      $                        BETA, C( 1 ), NK )
 500                   CALL DGEMM( 'T', 'N', NK, NK, K, ALPHA, A( 1, 1 ),
 501      $                        LDA, A( 1, NK+1 ), LDA, BETA,
 502      $                        C( ( ( NK+1 )*NK )+1 ), NK )
 503 *
 504                END IF
 505 *
 506             ELSE
 507 *
 508 *              N is even, TRANSR = 'T', and UPLO = 'U'
 509 *
 510                IF( NOTRANS ) THEN
 511 *
 512 *                 N is even, TRANSR = 'T', UPLO = 'U', and TRANS = 'N'
 513 *
 514                   CALL DSYRK( 'U', 'N', NK, K, ALPHA, A( 1, 1 ), LDA,
 515      $                        BETA, C( NK*( NK+1 )+1 ), NK )
 516                   CALL DSYRK( 'L', 'N', NK, K, ALPHA, A( NK+1, 1 ), LDA,
 517      $                        BETA, C( NK*NK+1 ), NK )
 518                   CALL DGEMM( 'N', 'T', NK, NK, K, ALPHA, A( NK+1, 1 ),
 519      $                        LDA, A( 1, 1 ), LDA, BETA, C( 1 ), NK )
 520 *
 521                ELSE
 522 *
 523 *                 N is even, TRANSR = 'T', UPLO = 'U', and TRANS = 'T'
 524 *
 525                   CALL DSYRK( 'U', 'T', NK, K, ALPHA, A( 1, 1 ), LDA,
 526      $                        BETA, C( NK*( NK+1 )+1 ), NK )
 527                   CALL DSYRK( 'L', 'T', NK, K, ALPHA, A( 1, NK+1 ), LDA,
 528      $                        BETA, C( NK*NK+1 ), NK )
 529                   CALL DGEMM( 'T', 'N', NK, NK, K, ALPHA, A( 1, NK+1 ),
 530      $                        LDA, A( 1, 1 ), LDA, BETA, C( 1 ), NK )
 531 *
 532                END IF
 533 *
 534             END IF
 535 *
 536          END IF
 537 *
 538       END IF
 539 *
 540       RETURN
 541 *
 542 *     End of DSFRK
 543 *
 544       END