File:  [local] / rpl / lapack / lapack / ztrttf.f
Revision 1.17: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:39:42 2023 UTC (8 months, 3 weeks ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_35, rpl-4_1_34, HEAD
Première mise à jour de lapack et blas.

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

CVSweb interface <joel.bertrand@systella.fr>