LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
ctpsv.f
Go to the documentation of this file.
1 *> \brief \b CTPSV
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 CTPSV(UPLO,TRANS,DIAG,N,AP,X,INCX)
12 *
13 * .. Scalar Arguments ..
14 * INTEGER INCX,N
15 * CHARACTER DIAG,TRANS,UPLO
16 * ..
17 * .. Array Arguments ..
18 * COMPLEX AP(*),X(*)
19 * ..
20 *
21 *
22 *> \par Purpose:
23 * =============
24 *>
25 *> \verbatim
26 *>
27 *> CTPSV solves one of the systems of equations
28 *>
29 *> A*x = b, or A**T*x = b, or A**H*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, supplied in packed form.
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**H*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] AP
85 *> \verbatim
86 *> AP is COMPLEX array of DIMENSION at least
87 *> ( ( n*( n + 1 ) )/2 ).
88 *> Before entry with UPLO = 'U' or 'u', the array AP must
89 *> contain the upper triangular matrix packed sequentially,
90 *> column by column, so that AP( 1 ) contains a( 1, 1 ),
91 *> AP( 2 ) and AP( 3 ) contain a( 1, 2 ) and a( 2, 2 )
92 *> respectively, and so on.
93 *> Before entry with UPLO = 'L' or 'l', the array AP must
94 *> contain the lower triangular matrix packed sequentially,
95 *> column by column, so that AP( 1 ) contains a( 1, 1 ),
96 *> AP( 2 ) and AP( 3 ) contain a( 2, 1 ) and a( 3, 1 )
97 *> respectively, and so on.
98 *> Note that when DIAG = 'U' or 'u', the diagonal elements of
99 *> A are not referenced, but are assumed to be unity.
100 *> \endverbatim
101 *>
102 *> \param[in,out] X
103 *> \verbatim
104 *> X is COMPLEX array of dimension at least
105 *> ( 1 + ( n - 1 )*abs( INCX ) ).
106 *> Before entry, the incremented array X must contain the n
107 *> element right-hand side vector b. On exit, X is overwritten
108 *> with the solution vector x.
109 *> \endverbatim
110 *>
111 *> \param[in] INCX
112 *> \verbatim
113 *> INCX is INTEGER
114 *> On entry, INCX specifies the increment for the elements of
115 *> X. INCX must not be zero.
116 *> \endverbatim
117 *
118 * Authors:
119 * ========
120 *
121 *> \author Univ. of Tennessee
122 *> \author Univ. of California Berkeley
123 *> \author Univ. of Colorado Denver
124 *> \author NAG Ltd.
125 *
126 *> \date November 2011
127 *
128 *> \ingroup complex_blas_level2
129 *
130 *> \par Further Details:
131 * =====================
132 *>
133 *> \verbatim
134 *>
135 *> Level 2 Blas routine.
136 *>
137 *> -- Written on 22-October-1986.
138 *> Jack Dongarra, Argonne National Lab.
139 *> Jeremy Du Croz, Nag Central Office.
140 *> Sven Hammarling, Nag Central Office.
141 *> Richard Hanson, Sandia National Labs.
142 *> \endverbatim
143 *>
144 * =====================================================================
145  SUBROUTINE ctpsv(UPLO,TRANS,DIAG,N,AP,X,INCX)
146 *
147 * -- Reference BLAS level2 routine (version 3.4.0) --
148 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
149 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
150 * November 2011
151 *
152 * .. Scalar Arguments ..
153  INTEGER incx,n
154  CHARACTER diag,trans,uplo
155 * ..
156 * .. Array Arguments ..
157  COMPLEX ap(*),x(*)
158 * ..
159 *
160 * =====================================================================
161 *
162 * .. Parameters ..
163  COMPLEX zero
164  parameter(zero= (0.0e+0,0.0e+0))
165 * ..
166 * .. Local Scalars ..
167  COMPLEX temp
168  INTEGER i,info,ix,j,jx,k,kk,kx
169  LOGICAL noconj,nounit
170 * ..
171 * .. External Functions ..
172  LOGICAL lsame
173  EXTERNAL lsame
174 * ..
175 * .. External Subroutines ..
176  EXTERNAL xerbla
177 * ..
178 * .. Intrinsic Functions ..
179  INTRINSIC conjg
180 * ..
181 *
182 * Test the input parameters.
183 *
184  info = 0
185  IF (.NOT.lsame(uplo,'U') .AND. .NOT.lsame(uplo,'L')) THEN
186  info = 1
187  ELSE IF (.NOT.lsame(trans,'N') .AND. .NOT.lsame(trans,'T') .AND.
188  + .NOT.lsame(trans,'C')) THEN
189  info = 2
190  ELSE IF (.NOT.lsame(diag,'U') .AND. .NOT.lsame(diag,'N')) THEN
191  info = 3
192  ELSE IF (n.LT.0) THEN
193  info = 4
194  ELSE IF (incx.EQ.0) THEN
195  info = 7
196  END IF
197  IF (info.NE.0) THEN
198  CALL xerbla('CTPSV ',info)
199  return
200  END IF
201 *
202 * Quick return if possible.
203 *
204  IF (n.EQ.0) return
205 *
206  noconj = lsame(trans,'T')
207  nounit = lsame(diag,'N')
208 *
209 * Set up the start point in X if the increment is not unity. This
210 * will be ( N - 1 )*INCX too small for descending loops.
211 *
212  IF (incx.LE.0) THEN
213  kx = 1 - (n-1)*incx
214  ELSE IF (incx.NE.1) THEN
215  kx = 1
216  END IF
217 *
218 * Start the operations. In this version the elements of AP are
219 * accessed sequentially with one pass through AP.
220 *
221  IF (lsame(trans,'N')) THEN
222 *
223 * Form x := inv( A )*x.
224 *
225  IF (lsame(uplo,'U')) THEN
226  kk = (n* (n+1))/2
227  IF (incx.EQ.1) THEN
228  DO 20 j = n,1,-1
229  IF (x(j).NE.zero) THEN
230  IF (nounit) x(j) = x(j)/ap(kk)
231  temp = x(j)
232  k = kk - 1
233  DO 10 i = j - 1,1,-1
234  x(i) = x(i) - temp*ap(k)
235  k = k - 1
236  10 continue
237  END IF
238  kk = kk - j
239  20 continue
240  ELSE
241  jx = kx + (n-1)*incx
242  DO 40 j = n,1,-1
243  IF (x(jx).NE.zero) THEN
244  IF (nounit) x(jx) = x(jx)/ap(kk)
245  temp = x(jx)
246  ix = jx
247  DO 30 k = kk - 1,kk - j + 1,-1
248  ix = ix - incx
249  x(ix) = x(ix) - temp*ap(k)
250  30 continue
251  END IF
252  jx = jx - incx
253  kk = kk - j
254  40 continue
255  END IF
256  ELSE
257  kk = 1
258  IF (incx.EQ.1) THEN
259  DO 60 j = 1,n
260  IF (x(j).NE.zero) THEN
261  IF (nounit) x(j) = x(j)/ap(kk)
262  temp = x(j)
263  k = kk + 1
264  DO 50 i = j + 1,n
265  x(i) = x(i) - temp*ap(k)
266  k = k + 1
267  50 continue
268  END IF
269  kk = kk + (n-j+1)
270  60 continue
271  ELSE
272  jx = kx
273  DO 80 j = 1,n
274  IF (x(jx).NE.zero) THEN
275  IF (nounit) x(jx) = x(jx)/ap(kk)
276  temp = x(jx)
277  ix = jx
278  DO 70 k = kk + 1,kk + n - j
279  ix = ix + incx
280  x(ix) = x(ix) - temp*ap(k)
281  70 continue
282  END IF
283  jx = jx + incx
284  kk = kk + (n-j+1)
285  80 continue
286  END IF
287  END IF
288  ELSE
289 *
290 * Form x := inv( A**T )*x or x := inv( A**H )*x.
291 *
292  IF (lsame(uplo,'U')) THEN
293  kk = 1
294  IF (incx.EQ.1) THEN
295  DO 110 j = 1,n
296  temp = x(j)
297  k = kk
298  IF (noconj) THEN
299  DO 90 i = 1,j - 1
300  temp = temp - ap(k)*x(i)
301  k = k + 1
302  90 continue
303  IF (nounit) temp = temp/ap(kk+j-1)
304  ELSE
305  DO 100 i = 1,j - 1
306  temp = temp - conjg(ap(k))*x(i)
307  k = k + 1
308  100 continue
309  IF (nounit) temp = temp/conjg(ap(kk+j-1))
310  END IF
311  x(j) = temp
312  kk = kk + j
313  110 continue
314  ELSE
315  jx = kx
316  DO 140 j = 1,n
317  temp = x(jx)
318  ix = kx
319  IF (noconj) THEN
320  DO 120 k = kk,kk + j - 2
321  temp = temp - ap(k)*x(ix)
322  ix = ix + incx
323  120 continue
324  IF (nounit) temp = temp/ap(kk+j-1)
325  ELSE
326  DO 130 k = kk,kk + j - 2
327  temp = temp - conjg(ap(k))*x(ix)
328  ix = ix + incx
329  130 continue
330  IF (nounit) temp = temp/conjg(ap(kk+j-1))
331  END IF
332  x(jx) = temp
333  jx = jx + incx
334  kk = kk + j
335  140 continue
336  END IF
337  ELSE
338  kk = (n* (n+1))/2
339  IF (incx.EQ.1) THEN
340  DO 170 j = n,1,-1
341  temp = x(j)
342  k = kk
343  IF (noconj) THEN
344  DO 150 i = n,j + 1,-1
345  temp = temp - ap(k)*x(i)
346  k = k - 1
347  150 continue
348  IF (nounit) temp = temp/ap(kk-n+j)
349  ELSE
350  DO 160 i = n,j + 1,-1
351  temp = temp - conjg(ap(k))*x(i)
352  k = k - 1
353  160 continue
354  IF (nounit) temp = temp/conjg(ap(kk-n+j))
355  END IF
356  x(j) = temp
357  kk = kk - (n-j+1)
358  170 continue
359  ELSE
360  kx = kx + (n-1)*incx
361  jx = kx
362  DO 200 j = n,1,-1
363  temp = x(jx)
364  ix = kx
365  IF (noconj) THEN
366  DO 180 k = kk,kk - (n- (j+1)),-1
367  temp = temp - ap(k)*x(ix)
368  ix = ix - incx
369  180 continue
370  IF (nounit) temp = temp/ap(kk-n+j)
371  ELSE
372  DO 190 k = kk,kk - (n- (j+1)),-1
373  temp = temp - conjg(ap(k))*x(ix)
374  ix = ix - incx
375  190 continue
376  IF (nounit) temp = temp/conjg(ap(kk-n+j))
377  END IF
378  x(jx) = temp
379  jx = jx - incx
380  kk = kk - (n-j+1)
381  200 continue
382  END IF
383  END IF
384  END IF
385 *
386  return
387 *
388 * End of CTPSV .
389 *
390  END