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