LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
ztrttf.f
Go to the documentation of this file.
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 *> \date September 2012
105 *
106 *> \ingroup complex16OTHERcomputational
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 ztrttf( TRANSR, UPLO, N, A, LDA, ARF, 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*16 a( 0: lda-1, 0: * ), arf( 0: * )
230 * ..
231 *
232 * =====================================================================
233 *
234 * .. Parameters ..
235 * ..
236 * .. Local Scalars ..
237  LOGICAL lower, nisodd, normaltransr
238  INTEGER i, ij, j, k, l, n1, n2, nt, nx2, np1x2
239 * ..
240 * .. External Functions ..
241  LOGICAL lsame
242  EXTERNAL lsame
243 * ..
244 * .. External Subroutines ..
245  EXTERNAL xerbla
246 * ..
247 * .. Intrinsic Functions ..
248  INTRINSIC dconjg, max, mod
249 * ..
250 * .. Executable Statements ..
251 *
252 * Test the input parameters.
253 *
254  info = 0
255  normaltransr = lsame( transr, 'N' )
256  lower = lsame( uplo, 'L' )
257  IF( .NOT.normaltransr .AND. .NOT.lsame( transr, 'C' ) ) THEN
258  info = -1
259  ELSE IF( .NOT.lower .AND. .NOT.lsame( uplo, 'U' ) ) THEN
260  info = -2
261  ELSE IF( n.LT.0 ) THEN
262  info = -3
263  ELSE IF( lda.LT.max( 1, n ) ) THEN
264  info = -5
265  END IF
266  IF( info.NE.0 ) THEN
267  CALL xerbla( 'ZTRTTF', -info )
268  return
269  END IF
270 *
271 * Quick return if possible
272 *
273  IF( n.LE.1 ) THEN
274  IF( n.EQ.1 ) THEN
275  IF( normaltransr ) THEN
276  arf( 0 ) = a( 0, 0 )
277  ELSE
278  arf( 0 ) = dconjg( a( 0, 0 ) )
279  END IF
280  END IF
281  return
282  END IF
283 *
284 * Size of array ARF(1:2,0:nt-1)
285 *
286  nt = n*( n+1 ) / 2
287 *
288 * set N1 and N2 depending on LOWER: for N even N1=N2=K
289 *
290  IF( lower ) THEN
291  n2 = n / 2
292  n1 = n - n2
293  ELSE
294  n1 = n / 2
295  n2 = n - n1
296  END IF
297 *
298 * If N is odd, set NISODD = .TRUE., LDA=N+1 and A is (N+1)--by--K2.
299 * If N is even, set K = N/2 and NISODD = .FALSE., LDA=N and A is
300 * N--by--(N+1)/2.
301 *
302  IF( mod( n, 2 ).EQ.0 ) THEN
303  k = n / 2
304  nisodd = .false.
305  IF( .NOT.lower )
306  $ np1x2 = n + n + 2
307  ELSE
308  nisodd = .true.
309  IF( .NOT.lower )
310  $ nx2 = n + n
311  END IF
312 *
313  IF( nisodd ) THEN
314 *
315 * N is odd
316 *
317  IF( normaltransr ) THEN
318 *
319 * N is odd and TRANSR = 'N'
320 *
321  IF( lower ) THEN
322 *
323 * SRPA for LOWER, NORMAL and N is odd ( a(0:n-1,0:n1-1) )
324 * T1 -> a(0,0), T2 -> a(0,1), S -> a(n1,0)
325 * T1 -> a(0), T2 -> a(n), S -> a(n1); lda=n
326 *
327  ij = 0
328  DO j = 0, n2
329  DO i = n1, n2 + j
330  arf( ij ) = dconjg( a( n2+j, i ) )
331  ij = ij + 1
332  END DO
333  DO i = j, n - 1
334  arf( ij ) = a( i, j )
335  ij = ij + 1
336  END DO
337  END DO
338 *
339  ELSE
340 *
341 * SRPA for UPPER, NORMAL and N is odd ( a(0:n-1,0:n2-1)
342 * T1 -> a(n1+1,0), T2 -> a(n1,0), S -> a(0,0)
343 * T1 -> a(n2), T2 -> a(n1), S -> a(0); lda=n
344 *
345  ij = nt - n
346  DO j = n - 1, n1, -1
347  DO i = 0, j
348  arf( ij ) = a( i, j )
349  ij = ij + 1
350  END DO
351  DO l = j - n1, n1 - 1
352  arf( ij ) = dconjg( a( j-n1, l ) )
353  ij = ij + 1
354  END DO
355  ij = ij - nx2
356  END DO
357 *
358  END IF
359 *
360  ELSE
361 *
362 * N is odd and TRANSR = 'C'
363 *
364  IF( lower ) THEN
365 *
366 * SRPA for LOWER, TRANSPOSE and N is odd
367 * T1 -> A(0,0) , T2 -> A(1,0) , S -> A(0,n1)
368 * T1 -> A(0+0) , T2 -> A(1+0) , S -> A(0+n1*n1); lda=n1
369 *
370  ij = 0
371  DO j = 0, n2 - 1
372  DO i = 0, j
373  arf( ij ) = dconjg( a( j, i ) )
374  ij = ij + 1
375  END DO
376  DO i = n1 + j, n - 1
377  arf( ij ) = a( i, n1+j )
378  ij = ij + 1
379  END DO
380  END DO
381  DO j = n2, n - 1
382  DO i = 0, n1 - 1
383  arf( ij ) = dconjg( a( j, i ) )
384  ij = ij + 1
385  END DO
386  END DO
387 *
388  ELSE
389 *
390 * SRPA for UPPER, TRANSPOSE and N is odd
391 * T1 -> A(0,n1+1), T2 -> A(0,n1), S -> A(0,0)
392 * T1 -> A(n2*n2), T2 -> A(n1*n2), S -> A(0); lda=n2
393 *
394  ij = 0
395  DO j = 0, n1
396  DO i = n1, n - 1
397  arf( ij ) = dconjg( a( j, i ) )
398  ij = ij + 1
399  END DO
400  END DO
401  DO j = 0, n1 - 1
402  DO i = 0, j
403  arf( ij ) = a( i, j )
404  ij = ij + 1
405  END DO
406  DO l = n2 + j, n - 1
407  arf( ij ) = dconjg( a( n2+j, l ) )
408  ij = ij + 1
409  END DO
410  END DO
411 *
412  END IF
413 *
414  END IF
415 *
416  ELSE
417 *
418 * N is even
419 *
420  IF( normaltransr ) THEN
421 *
422 * N is even and TRANSR = 'N'
423 *
424  IF( lower ) THEN
425 *
426 * SRPA for LOWER, NORMAL, and N is even ( a(0:n,0:k-1) )
427 * T1 -> a(1,0), T2 -> a(0,0), S -> a(k+1,0)
428 * T1 -> a(1), T2 -> a(0), S -> a(k+1); lda=n+1
429 *
430  ij = 0
431  DO j = 0, k - 1
432  DO i = k, k + j
433  arf( ij ) = dconjg( a( k+j, i ) )
434  ij = ij + 1
435  END DO
436  DO i = j, n - 1
437  arf( ij ) = a( i, j )
438  ij = ij + 1
439  END DO
440  END DO
441 *
442  ELSE
443 *
444 * SRPA for UPPER, NORMAL, and N is even ( a(0:n,0:k-1) )
445 * T1 -> a(k+1,0) , T2 -> a(k,0), S -> a(0,0)
446 * T1 -> a(k+1), T2 -> a(k), S -> a(0); lda=n+1
447 *
448  ij = nt - n - 1
449  DO j = n - 1, k, -1
450  DO i = 0, j
451  arf( ij ) = a( i, j )
452  ij = ij + 1
453  END DO
454  DO l = j - k, k - 1
455  arf( ij ) = dconjg( a( j-k, l ) )
456  ij = ij + 1
457  END DO
458  ij = ij - np1x2
459  END DO
460 *
461  END IF
462 *
463  ELSE
464 *
465 * N is even and TRANSR = 'C'
466 *
467  IF( lower ) THEN
468 *
469 * SRPA for LOWER, TRANSPOSE and N is even (see paper, A=B)
470 * T1 -> A(0,1) , T2 -> A(0,0) , S -> A(0,k+1) :
471 * T1 -> A(0+k) , T2 -> A(0+0) , S -> A(0+k*(k+1)); lda=k
472 *
473  ij = 0
474  j = k
475  DO i = k, n - 1
476  arf( ij ) = a( i, j )
477  ij = ij + 1
478  END DO
479  DO j = 0, k - 2
480  DO i = 0, j
481  arf( ij ) = dconjg( a( j, i ) )
482  ij = ij + 1
483  END DO
484  DO i = k + 1 + j, n - 1
485  arf( ij ) = a( i, k+1+j )
486  ij = ij + 1
487  END DO
488  END DO
489  DO j = k - 1, n - 1
490  DO i = 0, k - 1
491  arf( ij ) = dconjg( a( j, i ) )
492  ij = ij + 1
493  END DO
494  END DO
495 *
496  ELSE
497 *
498 * SRPA for UPPER, TRANSPOSE and N is even (see paper, A=B)
499 * T1 -> A(0,k+1) , T2 -> A(0,k) , S -> A(0,0)
500 * T1 -> A(0+k*(k+1)) , T2 -> A(0+k*k) , S -> A(0+0)); lda=k
501 *
502  ij = 0
503  DO j = 0, k
504  DO i = k, n - 1
505  arf( ij ) = dconjg( a( j, i ) )
506  ij = ij + 1
507  END DO
508  END DO
509  DO j = 0, k - 2
510  DO i = 0, j
511  arf( ij ) = a( i, j )
512  ij = ij + 1
513  END DO
514  DO l = k + 1 + j, n - 1
515  arf( ij ) = dconjg( a( k+1+j, l ) )
516  ij = ij + 1
517  END DO
518  END DO
519 *
520 * Note that here J = K-1
521 *
522  DO i = 0, j
523  arf( ij ) = a( i, j )
524  ij = ij + 1
525  END DO
526 *
527  END IF
528 *
529  END IF
530 *
531  END IF
532 *
533  return
534 *
535 * End of ZTRTTF
536 *
537  END