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