1*> \brief \b DTRSV
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*  Definition:
9*  ===========
10*
11*       SUBROUTINE DTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
12*
13*       .. Scalar Arguments ..
14*       INTEGER INCX,LDA,N
15*       CHARACTER DIAG,TRANS,UPLO
16*       ..
17*       .. Array Arguments ..
18*       DOUBLE PRECISION A(LDA,*),X(*)
19*       ..
20*
21*
22*> \par Purpose:
23*  =============
24*>
25*> \verbatim
26*>
27*> DTRSV  solves one of the systems of equations
28*>
29*>    A*x = b,   or   A**T*x = b,
30*>
31*> where b and x are n element vectors and A is an n by n unit, or
32*> non-unit, upper or lower triangular matrix.
33*>
34*> No test for singularity or near-singularity is included in this
35*> routine. Such tests must be performed before calling this routine.
36*> \endverbatim
37*
38*  Arguments:
39*  ==========
40*
41*> \param[in] UPLO
42*> \verbatim
43*>          UPLO is CHARACTER*1
44*>           On entry, UPLO specifies whether the matrix is an upper or
45*>           lower triangular matrix as follows:
46*>
47*>              UPLO = 'U' or 'u'   A is an upper triangular matrix.
48*>
49*>              UPLO = 'L' or 'l'   A is a lower triangular matrix.
50*> \endverbatim
51*>
52*> \param[in] TRANS
53*> \verbatim
54*>          TRANS is CHARACTER*1
55*>           On entry, TRANS specifies the equations to be solved as
56*>           follows:
57*>
58*>              TRANS = 'N' or 'n'   A*x = b.
59*>
60*>              TRANS = 'T' or 't'   A**T*x = b.
61*>
62*>              TRANS = 'C' or 'c'   A**T*x = b.
63*> \endverbatim
64*>
65*> \param[in] DIAG
66*> \verbatim
67*>          DIAG is CHARACTER*1
68*>           On entry, DIAG specifies whether or not A is unit
69*>           triangular as follows:
70*>
71*>              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
72*>
73*>              DIAG = 'N' or 'n'   A is not assumed to be unit
74*>                                  triangular.
75*> \endverbatim
76*>
77*> \param[in] N
78*> \verbatim
79*>          N is INTEGER
80*>           On entry, N specifies the order of the matrix A.
81*>           N must be at least zero.
82*> \endverbatim
83*>
84*> \param[in] A
85*> \verbatim
86*>          A is DOUBLE PRECISION array, dimension ( LDA, N )
87*>           Before entry with  UPLO = 'U' or 'u', the leading n by n
88*>           upper triangular part of the array A must contain the upper
89*>           triangular matrix and the strictly lower triangular part of
90*>           A is not referenced.
91*>           Before entry with UPLO = 'L' or 'l', the leading n by n
92*>           lower triangular part of the array A must contain the lower
93*>           triangular matrix and the strictly upper triangular part of
94*>           A is not referenced.
95*>           Note that when  DIAG = 'U' or 'u', the diagonal elements of
96*>           A are not referenced either, but are assumed to be unity.
97*> \endverbatim
98*>
99*> \param[in] LDA
100*> \verbatim
101*>          LDA is INTEGER
102*>           On entry, LDA specifies the first dimension of A as declared
103*>           in the calling (sub) program. LDA must be at least
104*>           max( 1, n ).
105*> \endverbatim
106*>
107*> \param[in,out] X
108*> \verbatim
109*>          X is DOUBLE PRECISION array, dimension at least
110*>           ( 1 + ( n - 1 )*abs( INCX ) ).
111*>           Before entry, the incremented array X must contain the n
112*>           element right-hand side vector b. On exit, X is overwritten
113*>           with the solution vector x.
114*> \endverbatim
115*>
116*> \param[in] INCX
117*> \verbatim
118*>          INCX is INTEGER
119*>           On entry, INCX specifies the increment for the elements of
120*>           X. INCX must not be zero.
121*>
122*>  Level 2 Blas routine.
123*>
124*>  -- Written on 22-October-1986.
125*>     Jack Dongarra, Argonne National Lab.
126*>     Jeremy Du Croz, Nag Central Office.
127*>     Sven Hammarling, Nag Central Office.
128*>     Richard Hanson, Sandia National Labs.
129*> \endverbatim
130*
131*  Authors:
132*  ========
133*
134*> \author Univ. of Tennessee
135*> \author Univ. of California Berkeley
136*> \author Univ. of Colorado Denver
137*> \author NAG Ltd.
138*
139*> \ingroup double_blas_level1
140*
141*  =====================================================================
142      SUBROUTINE DTRSV(UPLO,TRANS,DIAG,N,A,LDA,X,INCX)
143*
144*  -- Reference BLAS level1 routine --
145*  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
146*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
147*
148*     .. Scalar Arguments ..
149      INTEGER INCX,LDA,N
150      CHARACTER DIAG,TRANS,UPLO
151*     ..
152*     .. Array Arguments ..
153      DOUBLE PRECISION A(LDA,*),X(*)
154*     ..
155*
156*  =====================================================================
157*
158*     .. Parameters ..
159      DOUBLE PRECISION ZERO
160      PARAMETER (ZERO=0.0D+0)
161*     ..
162*     .. Local Scalars ..
163      DOUBLE PRECISION TEMP
164      INTEGER I,INFO,IX,J,JX,KX
165      LOGICAL NOUNIT
166*     ..
167*     .. External Functions ..
168      LOGICAL LSAME
169      EXTERNAL LSAME
170*     ..
171*     .. External Subroutines ..
172      EXTERNAL XERBLA
173*     ..
174*     .. Intrinsic Functions ..
175      INTRINSIC MAX
176*     ..
177*
178*     Test the input parameters.
179*
180      INFO = 0
181      IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
182          INFO = 1
183      ELSE IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND.
184     +         .NOT.LSAME(TRANS,'C')) THEN
185          INFO = 2
186      ELSE IF (.NOT.LSAME(DIAG,'U') .AND. .NOT.LSAME(DIAG,'N')) THEN
187          INFO = 3
188      ELSE IF (N.LT.0) THEN
189          INFO = 4
190      ELSE IF (LDA.LT.MAX(1,N)) THEN
191          INFO = 6
192      ELSE IF (INCX.EQ.0) THEN
193          INFO = 8
194      END IF
195      IF (INFO.NE.0) THEN
196          CALL XERBLA('DTRSV ',INFO)
197          RETURN
198      END IF
199*
200*     Quick return if possible.
201*
202      IF (N.EQ.0) RETURN
203*
204      NOUNIT = LSAME(DIAG,'N')
205*
206*     Set up the start point in X if the increment is not unity. This
207*     will be  ( N - 1 )*INCX  too small for descending loops.
208*
209      IF (INCX.LE.0) THEN
210          KX = 1 - (N-1)*INCX
211      ELSE IF (INCX.NE.1) THEN
212          KX = 1
213      END IF
214*
215*     Start the operations. In this version the elements of A are
216*     accessed sequentially with one pass through A.
217*
218      IF (LSAME(TRANS,'N')) THEN
219*
220*        Form  x := inv( A )*x.
221*
222          IF (LSAME(UPLO,'U')) THEN
223              IF (INCX.EQ.1) THEN
224                  DO 20 J = N,1,-1
225                      IF (X(J).NE.ZERO) THEN
226                          IF (NOUNIT) X(J) = X(J)/A(J,J)
227                          TEMP = X(J)
228                          DO 10 I = J - 1,1,-1
229                              X(I) = X(I) - TEMP*A(I,J)
230   10                     CONTINUE
231                      END IF
232   20             CONTINUE
233              ELSE
234                  JX = KX + (N-1)*INCX
235                  DO 40 J = N,1,-1
236                      IF (X(JX).NE.ZERO) THEN
237                          IF (NOUNIT) X(JX) = X(JX)/A(J,J)
238                          TEMP = X(JX)
239                          IX = JX
240                          DO 30 I = J - 1,1,-1
241                              IX = IX - INCX
242                              X(IX) = X(IX) - TEMP*A(I,J)
243   30                     CONTINUE
244                      END IF
245                      JX = JX - INCX
246   40             CONTINUE
247              END IF
248          ELSE
249              IF (INCX.EQ.1) THEN
250                  DO 60 J = 1,N
251                      IF (X(J).NE.ZERO) THEN
252                          IF (NOUNIT) X(J) = X(J)/A(J,J)
253                          TEMP = X(J)
254                          DO 50 I = J + 1,N
255                              X(I) = X(I) - TEMP*A(I,J)
256   50                     CONTINUE
257                      END IF
258   60             CONTINUE
259              ELSE
260                  JX = KX
261                  DO 80 J = 1,N
262                      IF (X(JX).NE.ZERO) THEN
263                          IF (NOUNIT) X(JX) = X(JX)/A(J,J)
264                          TEMP = X(JX)
265                          IX = JX
266                          DO 70 I = J + 1,N
267                              IX = IX + INCX
268                              X(IX) = X(IX) - TEMP*A(I,J)
269   70                     CONTINUE
270                      END IF
271                      JX = JX + INCX
272   80             CONTINUE
273              END IF
274          END IF
275      ELSE
276*
277*        Form  x := inv( A**T )*x.
278*
279          IF (LSAME(UPLO,'U')) THEN
280              IF (INCX.EQ.1) THEN
281                  DO 100 J = 1,N
282                      TEMP = X(J)
283                      DO 90 I = 1,J - 1
284                          TEMP = TEMP - A(I,J)*X(I)
285   90                 CONTINUE
286                      IF (NOUNIT) TEMP = TEMP/A(J,J)
287                      X(J) = TEMP
288  100             CONTINUE
289              ELSE
290                  JX = KX
291                  DO 120 J = 1,N
292                      TEMP = X(JX)
293                      IX = KX
294                      DO 110 I = 1,J - 1
295                          TEMP = TEMP - A(I,J)*X(IX)
296                          IX = IX + INCX
297  110                 CONTINUE
298                      IF (NOUNIT) TEMP = TEMP/A(J,J)
299                      X(JX) = TEMP
300                      JX = JX + INCX
301  120             CONTINUE
302              END IF
303          ELSE
304              IF (INCX.EQ.1) THEN
305                  DO 140 J = N,1,-1
306                      TEMP = X(J)
307                      DO 130 I = N,J + 1,-1
308                          TEMP = TEMP - A(I,J)*X(I)
309  130                 CONTINUE
310                      IF (NOUNIT) TEMP = TEMP/A(J,J)
311                      X(J) = TEMP
312  140             CONTINUE
313              ELSE
314                  KX = KX + (N-1)*INCX
315                  JX = KX
316                  DO 160 J = N,1,-1
317                      TEMP = X(JX)
318                      IX = KX
319                      DO 150 I = N,J + 1,-1
320                          TEMP = TEMP - A(I,J)*X(IX)
321                          IX = IX - INCX
322  150                 CONTINUE
323                      IF (NOUNIT) TEMP = TEMP/A(J,J)
324                      X(JX) = TEMP
325                      JX = JX - INCX
326  160             CONTINUE
327              END IF
328          END IF
329      END IF
330*
331      RETURN
332*
333*     End of DTRSV
334*
335      END
336