File:  [local] / rpl / lapack / lapack / dtfttr.f
Revision 1.17: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:39:12 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 DTFTTR copies a triangular matrix from the rectangular full packed format (TF) to the standard full format (TR).
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *> \htmlonly
    9: *> Download DTFTTR + dependencies
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtfttr.f">
   11: *> [TGZ]</a>
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtfttr.f">
   13: *> [ZIP]</a>
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtfttr.f">
   15: *> [TXT]</a>
   16: *> \endhtmlonly
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DTFTTR( TRANSR, UPLO, N, ARF, A, LDA, 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: *> DTFTTR copies a triangular matrix A from rectangular full packed
   38: *> format (TF) to standard full format (TR).
   39: *> \endverbatim
   40: *
   41: *  Arguments:
   42: *  ==========
   43: *
   44: *> \param[in] TRANSR
   45: *> \verbatim
   46: *>          TRANSR is CHARACTER*1
   47: *>          = 'N':  ARF is in Normal format;
   48: *>          = 'T':  ARF is in Transpose format.
   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 matrices ARF and A. N >= 0.
   62: *> \endverbatim
   63: *>
   64: *> \param[in] ARF
   65: *> \verbatim
   66: *>          ARF is DOUBLE PRECISION array, dimension (N*(N+1)/2).
   67: *>          On entry, the upper (if UPLO = 'U') or lower (if UPLO = 'L')
   68: *>          matrix A in RFP format. See the "Notes" below for more
   69: *>          details.
   70: *> \endverbatim
   71: *>
   72: *> \param[out] A
   73: *> \verbatim
   74: *>          A is DOUBLE PRECISION array, dimension (LDA,N)
   75: *>          On exit, the triangular matrix A.  If UPLO = 'U', the
   76: *>          leading N-by-N upper triangular part of the array A contains
   77: *>          the upper triangular matrix, and the strictly lower
   78: *>          triangular part of A is not referenced.  If UPLO = 'L', the
   79: *>          leading N-by-N lower triangular part of the array A contains
   80: *>          the lower triangular matrix, and the strictly upper
   81: *>          triangular part of A is not referenced.
   82: *> \endverbatim
   83: *>
   84: *> \param[in] LDA
   85: *> \verbatim
   86: *>          LDA is INTEGER
   87: *>          The leading dimension of the array A.  LDA >= max(1,N).
   88: *> \endverbatim
   89: *>
   90: *> \param[out] INFO
   91: *> \verbatim
   92: *>          INFO is INTEGER
   93: *>          = 0:  successful exit
   94: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
   95: *> \endverbatim
   96: *
   97: *  Authors:
   98: *  ========
   99: *
  100: *> \author Univ. of Tennessee
  101: *> \author Univ. of California Berkeley
  102: *> \author Univ. of Colorado Denver
  103: *> \author NAG Ltd.
  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 DTFTTR( TRANSR, UPLO, N, ARF, A, LDA, INFO )
  196: *
  197: *  -- LAPACK computational routine --
  198: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  199: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  200: *
  201: *     .. Scalar Arguments ..
  202:       CHARACTER          TRANSR, UPLO
  203:       INTEGER            INFO, N, LDA
  204: *     ..
  205: *     .. Array Arguments ..
  206:       DOUBLE PRECISION   A( 0: LDA-1, 0: * ), ARF( 0: * )
  207: *     ..
  208: *
  209: *  =====================================================================
  210: *
  211: *     ..
  212: *     .. Local Scalars ..
  213:       LOGICAL            LOWER, NISODD, NORMALTRANSR
  214:       INTEGER            N1, N2, K, NT, NX2, NP1X2
  215:       INTEGER            I, J, L, IJ
  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 = -6
  242:       END IF
  243:       IF( INFO.NE.0 ) THEN
  244:          CALL XERBLA( 'DTFTTR', -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:             A( 0, 0 ) = ARF( 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:                      A( N2+J, I ) = ARF( IJ )
  302:                      IJ = IJ + 1
  303:                   END DO
  304:                   DO I = J, N - 1
  305:                      A( I, J ) = ARF( IJ )
  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:                      A( I, J ) = ARF( IJ )
  318:                      IJ = IJ + 1
  319:                   END DO
  320:                   DO L = J - N1, N1 - 1
  321:                      A( J-N1, L ) = ARF( IJ )
  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:                      A( J, I ) = ARF( IJ )
  341:                      IJ = IJ + 1
  342:                   END DO
  343:                   DO I = N1 + J, N - 1
  344:                      A( I, N1+J ) = ARF( IJ )
  345:                      IJ = IJ + 1
  346:                   END DO
  347:                END DO
  348:                DO J = N2, N - 1
  349:                   DO I = 0, N1 - 1
  350:                      A( J, I ) = ARF( IJ )
  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:                      A( J, I ) = ARF( IJ )
  363:                      IJ = IJ + 1
  364:                   END DO
  365:                END DO
  366:                DO J = 0, N1 - 1
  367:                   DO I = 0, J
  368:                      A( I, J ) = ARF( IJ )
  369:                      IJ = IJ + 1
  370:                   END DO
  371:                   DO L = N2 + J, N - 1
  372:                      A( N2+J, L ) = ARF( IJ )
  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:                      A( K+J, I ) = ARF( IJ )
  397:                      IJ = IJ + 1
  398:                   END DO
  399:                   DO I = J, N - 1
  400:                      A( I, J ) = ARF( IJ )
  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:                      A( I, J ) = ARF( IJ )
  413:                      IJ = IJ + 1
  414:                   END DO
  415:                   DO L = J - K, K - 1
  416:                      A( J-K, L ) = ARF( IJ )
  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:                   A( I, J ) = ARF( IJ )
  436:                   IJ = IJ + 1
  437:                END DO
  438:                DO J = 0, K - 2
  439:                   DO I = 0, J
  440:                      A( J, I ) = ARF( IJ )
  441:                      IJ = IJ + 1
  442:                   END DO
  443:                   DO I = K + 1 + J, N - 1
  444:                      A( I, K+1+J ) = ARF( IJ )
  445:                      IJ = IJ + 1
  446:                   END DO
  447:                END DO
  448:                DO J = K - 1, N - 1
  449:                   DO I = 0, K - 1
  450:                      A( J, I ) = ARF( IJ )
  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:                      A( J, I ) = ARF( IJ )
  463:                      IJ = IJ + 1
  464:                   END DO
  465:                END DO
  466:                DO J = 0, K - 2
  467:                   DO I = 0, J
  468:                      A( I, J ) = ARF( IJ )
  469:                      IJ = IJ + 1
  470:                   END DO
  471:                   DO L = K + 1 + J, N - 1
  472:                      A( K+1+J, L ) = ARF( IJ )
  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:                   A( I, J ) = ARF( IJ )
  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 DTFTTR
  491: *
  492:       END

CVSweb interface <joel.bertrand@systella.fr>