SRC/dtrttf.f

   1 *> \brief \b DTRTTF copies a triangular matrix from the standard full format (TR) to the rectangular full packed format (TF).
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download DTRTTF + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtrttf.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtrttf.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtrttf.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE DTRTTF( TRANSR, UPLO, N, A, LDA, ARF, INFO )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       CHARACTER          TRANSR, UPLO
  25 *       INTEGER            INFO, N, LDA
  26 *       ..
  27 *       .. Array Arguments ..
  28 *       DOUBLE PRECISION   A( 0: LDA-1, 0: * ), ARF( 0: * )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *> DTRTTF copies a triangular matrix A from standard full format (TR)
  38 *> to rectangular full packed format (TF) .
  39 *> \endverbatim
  40 *
  41 *  Arguments:
  42 *  ==========
  43 *
  44 *> \param[in] TRANSR
  45 *> \verbatim
  46 *>          TRANSR is CHARACTER*1
  47 *>          = 'N':  ARF in Normal form is wanted;
  48 *>          = 'T':  ARF in Transpose form is wanted.
  49 *> \endverbatim
  50 *>
  51 *> \param[in] UPLO
  52 *> \verbatim
  53 *>          UPLO is CHARACTER*1
  54 *>          = 'U':  Upper triangle of A is stored;
  55 *>          = 'L':  Lower triangle of A is stored.
  56 *> \endverbatim
  57 *>
  58 *> \param[in] N
  59 *> \verbatim
  60 *>          N is INTEGER
  61 *>          The order of the matrix A. N >= 0.
  62 *> \endverbatim
  63 *>
  64 *> \param[in] A
  65 *> \verbatim
  66 *>          A is DOUBLE PRECISION array, dimension (LDA,N).
  67 *>          On entry, the triangular matrix A.  If UPLO = 'U', the
  68 *>          leading N-by-N upper triangular part of the array A contains
  69 *>          the upper triangular matrix, and the strictly lower
  70 *>          triangular part of A is not referenced.  If UPLO = 'L', the
  71 *>          leading N-by-N lower triangular part of the array A contains
  72 *>          the lower triangular matrix, and the strictly upper
  73 *>          triangular part of A is not referenced.
  74 *> \endverbatim
  75 *>
  76 *> \param[in] LDA
  77 *> \verbatim
  78 *>          LDA is INTEGER
  79 *>          The leading dimension of the matrix A. LDA >= max(1,N).
  80 *> \endverbatim
  81 *>
  82 *> \param[out] ARF
  83 *> \verbatim
  84 *>          ARF is DOUBLE PRECISION array, dimension (NT).
  85 *>          NT=N*(N+1)/2. On exit, the triangular matrix A in RFP format.
  86 *> \endverbatim
  87 *>
  88 *> \param[out] INFO
  89 *> \verbatim
  90 *>          INFO is INTEGER
  91 *>          = 0:  successful exit
  92 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
  93 *> \endverbatim
  94 *
  95 *  Authors:
  96 *  ========
  97 *
  98 *> \author Univ. of Tennessee
  99 *> \author Univ. of California Berkeley
 100 *> \author Univ. of Colorado Denver
 101 *> \author NAG Ltd.
 102 *
 103 *> \date September 2012
 104 *
 105 *> \ingroup doubleOTHERcomputational
 106 *
 107 *> \par Further Details:
 108 *  =====================
 109 *>
 110 *> \verbatim
 111 *>
 112 *>  We first consider Rectangular Full Packed (RFP) Format when N is
 113 *>  even. We give an example where N = 6.
 114 *>
 115 *>      AP is Upper             AP is Lower
 116 *>
 117 *>   00 01 02 03 04 05       00
 118 *>      11 12 13 14 15       10 11
 119 *>         22 23 24 25       20 21 22
 120 *>            33 34 35       30 31 32 33
 121 *>               44 45       40 41 42 43 44
 122 *>                  55       50 51 52 53 54 55
 123 *>
 124 *>
 125 *>  Let TRANSR = 'N'. RFP holds AP as follows:
 126 *>  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
 127 *>  three columns of AP upper. The lower triangle A(4:6,0:2) consists of
 128 *>  the transpose of the first three columns of AP upper.
 129 *>  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
 130 *>  three columns of AP lower. The upper triangle A(0:2,0:2) consists of
 131 *>  the transpose of the last three columns of AP lower.
 132 *>  This covers the case N even and TRANSR = 'N'.
 133 *>
 134 *>         RFP A                   RFP A
 135 *>
 136 *>        03 04 05                33 43 53
 137 *>        13 14 15                00 44 54
 138 *>        23 24 25                10 11 55
 139 *>        33 34 35                20 21 22
 140 *>        00 44 45                30 31 32
 141 *>        01 11 55                40 41 42
 142 *>        02 12 22                50 51 52
 143 *>
 144 *>  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
 145 *>  transpose of RFP A above. One therefore gets:
 146 *>
 147 *>
 148 *>           RFP A                   RFP A
 149 *>
 150 *>     03 13 23 33 00 01 02    33 00 10 20 30 40 50
 151 *>     04 14 24 34 44 11 12    43 44 11 21 31 41 51
 152 *>     05 15 25 35 45 55 22    53 54 55 22 32 42 52
 153 *>
 154 *>
 155 *>  We then consider Rectangular Full Packed (RFP) Format when N is
 156 *>  odd. We give an example where N = 5.
 157 *>
 158 *>     AP is Upper                 AP is Lower
 159 *>
 160 *>   00 01 02 03 04              00
 161 *>      11 12 13 14              10 11
 162 *>         22 23 24              20 21 22
 163 *>            33 34              30 31 32 33
 164 *>               44              40 41 42 43 44
 165 *>
 166 *>
 167 *>  Let TRANSR = 'N'. RFP holds AP as follows:
 168 *>  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
 169 *>  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
 170 *>  the transpose of the first two columns of AP upper.
 171 *>  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
 172 *>  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
 173 *>  the transpose of the last two columns of AP lower.
 174 *>  This covers the case N odd and TRANSR = 'N'.
 175 *>
 176 *>         RFP A                   RFP A
 177 *>
 178 *>        02 03 04                00 33 43
 179 *>        12 13 14                10 11 44
 180 *>        22 23 24                20 21 22
 181 *>        00 33 34                30 31 32
 182 *>        01 11 44                40 41 42
 183 *>
 184 *>  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
 185 *>  transpose of RFP A above. One therefore gets:
 186 *>
 187 *>           RFP A                   RFP A
 188 *>
 189 *>     02 12 22 00 01             00 10 20 30 40 50
 190 *>     03 13 23 33 11             33 11 21 31 41 51
 191 *>     04 14 24 34 44             43 44 22 32 42 52
 192 *> \endverbatim
 193 *
 194 *  =====================================================================
 195       SUBROUTINE DTRTTF( TRANSR, UPLO, N, A, LDA, ARF, INFO )
 196 *
 197 *  -- LAPACK computational routine (version 3.4.2) --
 198 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 199 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 200 *     September 2012
 201 *
 202 *     .. Scalar Arguments ..
 203       CHARACTER          TRANSR, UPLO
 204       INTEGER            INFO, N, LDA
 205 *     ..
 206 *     .. Array Arguments ..
 207       DOUBLE PRECISION   A( 0: LDA-1, 0: * ), ARF( 0: * )
 208 *     ..
 209 *
 210 *  =====================================================================
 211 *
 212 *     ..
 213 *     .. Local Scalars ..
 214       LOGICAL            LOWER, NISODD, NORMALTRANSR
 215       INTEGER            I, IJ, J, K, L, N1, N2, NT, NX2, NP1X2
 216 *     ..
 217 *     .. External Functions ..
 218       LOGICAL            LSAME
 219       EXTERNAL           LSAME
 220 *     ..
 221 *     .. External Subroutines ..
 222       EXTERNAL           XERBLA
 223 *     ..
 224 *     .. Intrinsic Functions ..
 225       INTRINSIC          MAX, MOD
 226 *     ..
 227 *     .. Executable Statements ..
 228 *
 229 *     Test the input parameters.
 230 *
 231       INFO = 0
 232       NORMALTRANSR = LSAME( TRANSR, 'N' )
 233       LOWER = LSAME( UPLO, 'L' )
 234       IF( .NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN
 235          INFO = -1
 236       ELSE IF( .NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
 237          INFO = -2
 238       ELSE IF( N.LT.0 ) THEN
 239          INFO = -3
 240       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
 241          INFO = -5
 242       END IF
 243       IF( INFO.NE.0 ) THEN
 244          CALL XERBLA( 'DTRTTF', -INFO )
 245          RETURN
 246       END IF
 247 *
 248 *     Quick return if possible
 249 *
 250       IF( N.LE.1 ) THEN
 251          IF( N.EQ.1 ) THEN
 252             ARF( 0 ) = A( 0, 0 )
 253          END IF
 254          RETURN
 255       END IF
 256 *
 257 *     Size of array ARF(0:nt-1)
 258 *
 259       NT = N*( N+1 ) / 2
 260 *
 261 *     Set N1 and N2 depending on LOWER: for N even N1=N2=K
 262 *
 263       IF( LOWER ) THEN
 264          N2 = N / 2
 265          N1 = N - N2
 266       ELSE
 267          N1 = N / 2
 268          N2 = N - N1
 269       END IF
 270 *
 271 *     If N is odd, set NISODD = .TRUE., LDA=N+1 and A is (N+1)--by--K2.
 272 *     If N is even, set K = N/2 and NISODD = .FALSE., LDA=N and A is
 273 *     N--by--(N+1)/2.
 274 *
 275       IF( MOD( N, 2 ).EQ.0 ) THEN
 276          K = N / 2
 277          NISODD = .FALSE.
 278          IF( .NOT.LOWER )
 279      $      NP1X2 = N + N + 2
 280       ELSE
 281          NISODD = .TRUE.
 282          IF( .NOT.LOWER )
 283      $      NX2 = N + N
 284       END IF
 285 *
 286       IF( NISODD ) THEN
 287 *
 288 *        N is odd
 289 *
 290          IF( NORMALTRANSR ) THEN
 291 *
 292 *           N is odd and TRANSR = 'N'
 293 *
 294             IF( LOWER ) THEN
 295 *
 296 *              N is odd, TRANSR = 'N', and UPLO = 'L'
 297 *
 298                IJ = 0
 299                DO J = 0, N2
 300                   DO I = N1, N2 + J
 301                      ARF( IJ ) = A( N2+J, I )
 302                      IJ = IJ + 1
 303                   END DO
 304                   DO I = J, N - 1
 305                      ARF( IJ ) = A( I, J )
 306                      IJ = IJ + 1
 307                   END DO
 308                END DO
 309 *
 310             ELSE
 311 *
 312 *              N is odd, TRANSR = 'N', and UPLO = 'U'
 313 *
 314                IJ = NT - N
 315                DO J = N - 1, N1, -1
 316                   DO I = 0, J
 317                      ARF( IJ ) = A( I, J )
 318                      IJ = IJ + 1
 319                   END DO
 320                   DO L = J - N1, N1 - 1
 321                      ARF( IJ ) = A( J-N1, L )
 322                      IJ = IJ + 1
 323                   END DO
 324                   IJ = IJ - NX2
 325                END DO
 326 *
 327             END IF
 328 *
 329          ELSE
 330 *
 331 *           N is odd and TRANSR = 'T'
 332 *
 333             IF( LOWER ) THEN
 334 *
 335 *              N is odd, TRANSR = 'T', and UPLO = 'L'
 336 *
 337                IJ = 0
 338                DO J = 0, N2 - 1
 339                   DO I = 0, J
 340                      ARF( IJ ) = A( J, I )
 341                      IJ = IJ + 1
 342                   END DO
 343                   DO I = N1 + J, N - 1
 344                      ARF( IJ ) = A( I, N1+J )
 345                      IJ = IJ + 1
 346                   END DO
 347                END DO
 348                DO J = N2, N - 1
 349                   DO I = 0, N1 - 1
 350                      ARF( IJ ) = A( J, I )
 351                      IJ = IJ + 1
 352                   END DO
 353                END DO
 354 *
 355             ELSE
 356 *
 357 *              N is odd, TRANSR = 'T', and UPLO = 'U'
 358 *
 359                IJ = 0
 360                DO J = 0, N1
 361                   DO I = N1, N - 1
 362                      ARF( IJ ) = A( J, I )
 363                      IJ = IJ + 1
 364                   END DO
 365                END DO
 366                DO J = 0, N1 - 1
 367                   DO I = 0, J
 368                      ARF( IJ ) = A( I, J )
 369                      IJ = IJ + 1
 370                   END DO
 371                   DO L = N2 + J, N - 1
 372                      ARF( IJ ) = A( N2+J, L )
 373                      IJ = IJ + 1
 374                   END DO
 375                END DO
 376 *
 377             END IF
 378 *
 379          END IF
 380 *
 381       ELSE
 382 *
 383 *        N is even
 384 *
 385          IF( NORMALTRANSR ) THEN
 386 *
 387 *           N is even and TRANSR = 'N'
 388 *
 389             IF( LOWER ) THEN
 390 *
 391 *              N is even, TRANSR = 'N', and UPLO = 'L'
 392 *
 393                IJ = 0
 394                DO J = 0, K - 1
 395                   DO I = K, K + J
 396                      ARF( IJ ) = A( K+J, I )
 397                      IJ = IJ + 1
 398                   END DO
 399                   DO I = J, N - 1
 400                      ARF( IJ ) = A( I, J )
 401                      IJ = IJ + 1
 402                   END DO
 403                END DO
 404 *
 405             ELSE
 406 *
 407 *              N is even, TRANSR = 'N', and UPLO = 'U'
 408 *
 409                IJ = NT - N - 1
 410                DO J = N - 1, K, -1
 411                   DO I = 0, J
 412                      ARF( IJ ) = A( I, J )
 413                      IJ = IJ + 1
 414                   END DO
 415                   DO L = J - K, K - 1
 416                      ARF( IJ ) = A( J-K, L )
 417                      IJ = IJ + 1
 418                   END DO
 419                   IJ = IJ - NP1X2
 420                END DO
 421 *
 422             END IF
 423 *
 424          ELSE
 425 *
 426 *           N is even and TRANSR = 'T'
 427 *
 428             IF( LOWER ) THEN
 429 *
 430 *              N is even, TRANSR = 'T', and UPLO = 'L'
 431 *
 432                IJ = 0
 433                J = K
 434                DO I = K, N - 1
 435                   ARF( IJ ) = A( I, J )
 436                   IJ = IJ + 1
 437                END DO
 438                DO J = 0, K - 2
 439                   DO I = 0, J
 440                      ARF( IJ ) = A( J, I )
 441                      IJ = IJ + 1
 442                   END DO
 443                   DO I = K + 1 + J, N - 1
 444                      ARF( IJ ) = A( I, K+1+J )
 445                      IJ = IJ + 1
 446                   END DO
 447                END DO
 448                DO J = K - 1, N - 1
 449                   DO I = 0, K - 1
 450                      ARF( IJ ) = A( J, I )
 451                      IJ = IJ + 1
 452                   END DO
 453                END DO
 454 *
 455             ELSE
 456 *
 457 *              N is even, TRANSR = 'T', and UPLO = 'U'
 458 *
 459                IJ = 0
 460                DO J = 0, K
 461                   DO I = K, N - 1
 462                      ARF( IJ ) = A( J, I )
 463                      IJ = IJ + 1
 464                   END DO
 465                END DO
 466                DO J = 0, K - 2
 467                   DO I = 0, J
 468                      ARF( IJ ) = A( I, J )
 469                      IJ = IJ + 1
 470                   END DO
 471                   DO L = K + 1 + J, N - 1
 472                      ARF( IJ ) = A( K+1+J, L )
 473                      IJ = IJ + 1
 474                   END DO
 475                END DO
 476 *              Note that here, on exit of the loop, J = K-1
 477                DO I = 0, J
 478                   ARF( IJ ) = A( I, J )
 479                   IJ = IJ + 1
 480                END DO
 481 *
 482             END IF
 483 *
 484          END IF
 485 *
 486       END IF
 487 *
 488       RETURN
 489 *
 490 *     End of DTRTTF
 491 *
 492       END