1       SUBROUTINE DTFTTR( TRANSR, UPLO, N, ARF, A, LDA, INFO )
  2 *
  3 *  -- LAPACK routine (version 3.3.1)                                    --
  4 *
  5 *  -- Contributed by Fred Gustavson of the IBM Watson Research Center --
  6 *  -- April 2011                                                      --
  7 *
  8 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  9 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 10 *
 11 *     .. Scalar Arguments ..
 12       CHARACTER          TRANSR, UPLO
 13       INTEGER            INFO, N, LDA
 14 *     ..
 15 *     .. Array Arguments ..
 16       DOUBLE PRECISION   A( 0: LDA-10* ), ARF( 0* )
 17 *     ..
 18 *
 19 *  Purpose
 20 *  =======
 21 *
 22 *  DTFTTR copies a triangular matrix A from rectangular full packed
 23 *  format (TF) to standard full format (TR).
 24 *
 25 *  Arguments
 26 *  =========
 27 *
 28 *  TRANSR  (input) CHARACTER*1
 29 *          = 'N':  ARF is in Normal format;
 30 *          = 'T':  ARF is in Transpose format.
 31 *
 32 *  UPLO    (input) CHARACTER*1
 33 *          = 'U':  A is upper triangular;
 34 *          = 'L':  A is lower triangular.
 35 *
 36 *  N       (input) INTEGER
 37 *          The order of the matrices ARF and A. N >= 0.
 38 *
 39 *  ARF     (input) DOUBLE PRECISION array, dimension (N*(N+1)/2).
 40 *          On entry, the upper (if UPLO = 'U') or lower (if UPLO = 'L')
 41 *          matrix A in RFP format. See the "Notes" below for more
 42 *          details.
 43 *
 44 *  A       (output) DOUBLE PRECISION array, dimension (LDA,N)
 45 *          On exit, the triangular matrix A.  If UPLO = 'U', the
 46 *          leading N-by-N upper triangular part of the array A contains
 47 *          the upper triangular matrix, and the strictly lower
 48 *          triangular part of A is not referenced.  If UPLO = 'L', the
 49 *          leading N-by-N lower triangular part of the array A contains
 50 *          the lower triangular matrix, and the strictly upper
 51 *          triangular part of A is not referenced.
 52 *
 53 *  LDA     (input) INTEGER
 54 *          The leading dimension of the array A.  LDA >= max(1,N).
 55 *
 56 *  INFO    (output) INTEGER
 57 *          = 0:  successful exit
 58 *          < 0:  if INFO = -i, the i-th argument had an illegal value
 59 *
 60 *  Further Details
 61 *  ===============
 62 *
 63 *  We first consider Rectangular Full Packed (RFP) Format when N is
 64 *  even. We give an example where N = 6.
 65 *
 66 *      AP is Upper             AP is Lower
 67 *
 68 *   00 01 02 03 04 05       00
 69 *      11 12 13 14 15       10 11
 70 *         22 23 24 25       20 21 22
 71 *            33 34 35       30 31 32 33
 72 *               44 45       40 41 42 43 44
 73 *                  55       50 51 52 53 54 55
 74 *
 75 *
 76 *  Let TRANSR = 'N'. RFP holds AP as follows:
 77 *  For UPLO = 'U' the upper trapezoid A(0:5,0:2) consists of the last
 78 *  three columns of AP upper. The lower triangle A(4:6,0:2) consists of
 79 *  the transpose of the first three columns of AP upper.
 80 *  For UPLO = 'L' the lower trapezoid A(1:6,0:2) consists of the first
 81 *  three columns of AP lower. The upper triangle A(0:2,0:2) consists of
 82 *  the transpose of the last three columns of AP lower.
 83 *  This covers the case N even and TRANSR = 'N'.
 84 *
 85 *         RFP A                   RFP A
 86 *
 87 *        03 04 05                33 43 53
 88 *        13 14 15                00 44 54
 89 *        23 24 25                10 11 55
 90 *        33 34 35                20 21 22
 91 *        00 44 45                30 31 32
 92 *        01 11 55                40 41 42
 93 *        02 12 22                50 51 52
 94 *
 95 *  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
 96 *  transpose of RFP A above. One therefore gets:
 97 *
 98 *
 99 *           RFP A                   RFP A
100 *
101 *     03 13 23 33 00 01 02    33 00 10 20 30 40 50
102 *     04 14 24 34 44 11 12    43 44 11 21 31 41 51
103 *     05 15 25 35 45 55 22    53 54 55 22 32 42 52
104 *
105 *
106 *  We then consider Rectangular Full Packed (RFP) Format when N is
107 *  odd. We give an example where N = 5.
108 *
109 *     AP is Upper                 AP is Lower
110 *
111 *   00 01 02 03 04              00
112 *      11 12 13 14              10 11
113 *         22 23 24              20 21 22
114 *            33 34              30 31 32 33
115 *               44              40 41 42 43 44
116 *
117 *
118 *  Let TRANSR = 'N'. RFP holds AP as follows:
119 *  For UPLO = 'U' the upper trapezoid A(0:4,0:2) consists of the last
120 *  three columns of AP upper. The lower triangle A(3:4,0:1) consists of
121 *  the transpose of the first two columns of AP upper.
122 *  For UPLO = 'L' the lower trapezoid A(0:4,0:2) consists of the first
123 *  three columns of AP lower. The upper triangle A(0:1,1:2) consists of
124 *  the transpose of the last two columns of AP lower.
125 *  This covers the case N odd and TRANSR = 'N'.
126 *
127 *         RFP A                   RFP A
128 *
129 *        02 03 04                00 33 43
130 *        12 13 14                10 11 44
131 *        22 23 24                20 21 22
132 *        00 33 34                30 31 32
133 *        01 11 44                40 41 42
134 *
135 *  Now let TRANSR = 'T'. RFP A in both UPLO cases is just the
136 *  transpose of RFP A above. One therefore gets:
137 *
138 *           RFP A                   RFP A
139 *
140 *     02 12 22 00 01             00 10 20 30 40 50
141 *     03 13 23 33 11             33 11 21 31 41 51
142 *     04 14 24 34 44             43 44 22 32 42 52
143 *
144 *  Reference
145 *  =========
146 *
147 *  =====================================================================
148 *
149 *     ..
150 *     .. Local Scalars ..
151       LOGICAL            LOWER, NISODD, NORMALTRANSR
152       INTEGER            N1, N2, K, NT, NX2, NP1X2
153       INTEGER            I, J, L, IJ
154 *     ..
155 *     .. External Functions ..
156       LOGICAL            LSAME
157       EXTERNAL           LSAME
158 *     ..
159 *     .. External Subroutines ..
160       EXTERNAL           XERBLA
161 *     ..
162 *     .. Intrinsic Functions ..
163       INTRINSIC          MAXMOD
164 *     ..
165 *     .. Executable Statements ..
166 *
167 *     Test the input parameters.
168 *
169       INFO = 0
170       NORMALTRANSR = LSAME( TRANSR, 'N' )
171       LOWER = LSAME( UPLO, 'L' )
172       IF.NOT.NORMALTRANSR .AND. .NOT.LSAME( TRANSR, 'T' ) ) THEN
173          INFO = -1
174       ELSE IF.NOT.LOWER .AND. .NOT.LSAME( UPLO, 'U' ) ) THEN
175          INFO = -2
176       ELSE IF( N.LT.0 ) THEN
177          INFO = -3
178       ELSE IF( LDA.LT.MAX1, N ) ) THEN
179          INFO = -6
180       END IF
181       IF( INFO.NE.0 ) THEN
182          CALL XERBLA( 'DTFTTR'-INFO )
183          RETURN
184       END IF
185 *
186 *     Quick return if possible
187 *
188       IF( N.LE.1 ) THEN
189          IF( N.EQ.1 ) THEN
190             A( 00 ) = ARF( 0 )
191          END IF
192          RETURN
193       END IF
194 *
195 *     Size of array ARF(0:nt-1)
196 *
197       NT = N*( N+1 ) / 2
198 *
199 *     set N1 and N2 depending on LOWER: for N even N1=N2=K
200 *
201       IF( LOWER ) THEN
202          N2 = N / 2
203          N1 = N - N2
204       ELSE
205          N1 = N / 2
206          N2 = N - N1
207       END IF
208 *
209 *     If N is odd, set NISODD = .TRUE., LDA=N+1 and A is (N+1)--by--K2.
210 *     If N is even, set K = N/2 and NISODD = .FALSE., LDA=N and A is
211 *     N--by--(N+1)/2.
212 *
213       IFMOD( N, 2 ).EQ.0 ) THEN
214          K = N / 2
215          NISODD = .FALSE.
216          IF.NOT.LOWER )
217      $      NP1X2 = N + N + 2
218       ELSE
219          NISODD = .TRUE.
220          IF.NOT.LOWER )
221      $      NX2 = N + N
222       END IF
223 *
224       IF( NISODD ) THEN
225 *
226 *        N is odd
227 *
228          IF( NORMALTRANSR ) THEN
229 *
230 *           N is odd and TRANSR = 'N'
231 *
232             IF( LOWER ) THEN
233 *
234 *              N is odd, TRANSR = 'N', and UPLO = 'L'
235 *
236                IJ = 0
237                DO J = 0, N2
238                   DO I = N1, N2 + J
239                      A( N2+J, I ) = ARF( IJ )
240                      IJ = IJ + 1
241                   END DO
242                   DO I = J, N - 1
243                      A( I, J ) = ARF( IJ )
244                      IJ = IJ + 1
245                   END DO
246                END DO
247 *
248             ELSE
249 *
250 *              N is odd, TRANSR = 'N', and UPLO = 'U'
251 *
252                IJ = NT - N
253                DO J = N - 1, N1, -1
254                   DO I = 0, J
255                      A( I, J ) = ARF( IJ )
256                      IJ = IJ + 1
257                   END DO
258                   DO L = J - N1, N1 - 1
259                      A( J-N1, L ) = ARF( IJ )
260                      IJ = IJ + 1
261                   END DO
262                   IJ = IJ - NX2
263                END DO
264 *
265             END IF
266 *
267          ELSE
268 *
269 *           N is odd and TRANSR = 'T'
270 *
271             IF( LOWER ) THEN
272 *
273 *              N is odd, TRANSR = 'T', and UPLO = 'L'
274 *
275                IJ = 0
276                DO J = 0, N2 - 1
277                   DO I = 0, J
278                      A( J, I ) = ARF( IJ )
279                      IJ = IJ + 1
280                   END DO
281                   DO I = N1 + J, N - 1
282                      A( I, N1+J ) = ARF( IJ )
283                      IJ = IJ + 1
284                   END DO
285                END DO
286                DO J = N2, N - 1
287                   DO I = 0, N1 - 1
288                      A( J, I ) = ARF( IJ )
289                      IJ = IJ + 1
290                   END DO
291                END DO
292 *
293             ELSE
294 *
295 *              N is odd, TRANSR = 'T', and UPLO = 'U'
296 *
297                IJ = 0
298                DO J = 0, N1
299                   DO I = N1, N - 1
300                      A( J, I ) = ARF( IJ )
301                      IJ = IJ + 1
302                   END DO
303                END DO
304                DO J = 0, N1 - 1
305                   DO I = 0, J
306                      A( I, J ) = ARF( IJ )
307                      IJ = IJ + 1
308                   END DO
309                   DO L = N2 + J, N - 1
310                      A( N2+J, L ) = ARF( IJ )
311                      IJ = IJ + 1
312                   END DO
313                END DO
314 *
315             END IF
316 *
317          END IF
318 *
319       ELSE
320 *
321 *        N is even
322 *
323          IF( NORMALTRANSR ) THEN
324 *
325 *           N is even and TRANSR = 'N'
326 *
327             IF( LOWER ) THEN
328 *
329 *              N is even, TRANSR = 'N', and UPLO = 'L'
330 *
331                IJ = 0
332                DO J = 0, K - 1
333                   DO I = K, K + J
334                      A( K+J, I ) = ARF( IJ )
335                      IJ = IJ + 1
336                   END DO
337                   DO I = J, N - 1
338                      A( I, J ) = ARF( IJ )
339                      IJ = IJ + 1
340                   END DO
341                END DO
342 *
343             ELSE
344 *
345 *              N is even, TRANSR = 'N', and UPLO = 'U'
346 *
347                IJ = NT - N - 1
348                DO J = N - 1, K, -1
349                   DO I = 0, J
350                      A( I, J ) = ARF( IJ )
351                      IJ = IJ + 1
352                   END DO
353                   DO L = J - K, K - 1
354                      A( J-K, L ) = ARF( IJ )
355                      IJ = IJ + 1
356                   END DO
357                   IJ = IJ - NP1X2
358                END DO
359 *
360             END IF
361 *
362          ELSE
363 *
364 *           N is even and TRANSR = 'T'
365 *
366             IF( LOWER ) THEN
367 *
368 *              N is even, TRANSR = 'T', and UPLO = 'L'
369 *
370                IJ = 0
371                J = K
372                DO I = K, N - 1
373                   A( I, J ) = ARF( IJ )
374                   IJ = IJ + 1
375                END DO
376                DO J = 0, K - 2
377                   DO I = 0, J
378                      A( J, I ) = ARF( IJ )
379                      IJ = IJ + 1
380                   END DO
381                   DO I = K + 1 + J, N - 1
382                      A( I, K+1+J ) = ARF( IJ )
383                      IJ = IJ + 1
384                   END DO
385                END DO
386                DO J = K - 1, N - 1
387                   DO I = 0, K - 1
388                      A( J, I ) = ARF( IJ )
389                      IJ = IJ + 1
390                   END DO
391                END DO
392 *
393             ELSE
394 *
395 *              N is even, TRANSR = 'T', and UPLO = 'U'
396 *
397                IJ = 0
398                DO J = 0, K
399                   DO I = K, N - 1
400                      A( J, I ) = ARF( IJ )
401                      IJ = IJ + 1
402                   END DO
403                END DO
404                DO J = 0, K - 2
405                   DO I = 0, J
406                      A( I, J ) = ARF( IJ )
407                      IJ = IJ + 1
408                   END DO
409                   DO L = K + 1 + J, N - 1
410                      A( K+1+J, L ) = ARF( IJ )
411                      IJ = IJ + 1
412                   END DO
413                END DO
414 *              Note that here, on exit of the loop, J = K-1
415                DO I = 0, J
416                   A( I, J ) = ARF( IJ )
417                   IJ = IJ + 1
418                END DO
419 *
420             END IF
421 *
422          END IF
423 *
424       END IF
425 *
426       RETURN
427 *
428 *     End of DTFTTR
429 *
430       END