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