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