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