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