1: *> \brief \b ZPTTS2 solves a tridiagonal system of the form AX=B using the L D LH factorization computed by spttrf.
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download ZPTTS2 + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zptts2.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zptts2.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zptts2.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZPTTS2( IUPLO, N, NRHS, D, E, B, LDB )
22: *
23: * .. Scalar Arguments ..
24: * INTEGER IUPLO, LDB, N, NRHS
25: * ..
26: * .. Array Arguments ..
27: * DOUBLE PRECISION D( * )
28: * COMPLEX*16 B( LDB, * ), E( * )
29: * ..
30: *
31: *
32: *> \par Purpose:
33: * =============
34: *>
35: *> \verbatim
36: *>
37: *> ZPTTS2 solves a tridiagonal system of the form
38: *> A * X = B
39: *> using the factorization A = U**H *D*U or A = L*D*L**H computed by ZPTTRF.
40: *> D is a diagonal matrix specified in the vector D, U (or L) is a unit
41: *> bidiagonal matrix whose superdiagonal (subdiagonal) is specified in
42: *> the vector E, and X and B are N by NRHS matrices.
43: *> \endverbatim
44: *
45: * Arguments:
46: * ==========
47: *
48: *> \param[in] IUPLO
49: *> \verbatim
50: *> IUPLO is INTEGER
51: *> Specifies the form of the factorization and whether the
52: *> vector E is the superdiagonal of the upper bidiagonal factor
53: *> U or the subdiagonal of the lower bidiagonal factor L.
54: *> = 1: A = U**H *D*U, E is the superdiagonal of U
55: *> = 0: A = L*D*L**H, E is the subdiagonal of L
56: *> \endverbatim
57: *>
58: *> \param[in] N
59: *> \verbatim
60: *> N is INTEGER
61: *> The order of the tridiagonal matrix A. N >= 0.
62: *> \endverbatim
63: *>
64: *> \param[in] NRHS
65: *> \verbatim
66: *> NRHS is INTEGER
67: *> The number of right hand sides, i.e., the number of columns
68: *> of the matrix B. NRHS >= 0.
69: *> \endverbatim
70: *>
71: *> \param[in] D
72: *> \verbatim
73: *> D is DOUBLE PRECISION array, dimension (N)
74: *> The n diagonal elements of the diagonal matrix D from the
75: *> factorization A = U**H *D*U or A = L*D*L**H.
76: *> \endverbatim
77: *>
78: *> \param[in] E
79: *> \verbatim
80: *> E is COMPLEX*16 array, dimension (N-1)
81: *> If IUPLO = 1, the (n-1) superdiagonal elements of the unit
82: *> bidiagonal factor U from the factorization A = U**H*D*U.
83: *> If IUPLO = 0, the (n-1) subdiagonal elements of the unit
84: *> bidiagonal factor L from the factorization A = L*D*L**H.
85: *> \endverbatim
86: *>
87: *> \param[in,out] B
88: *> \verbatim
89: *> B is COMPLEX*16 array, dimension (LDB,NRHS)
90: *> On entry, the right hand side vectors B for the system of
91: *> linear equations.
92: *> On exit, the solution vectors, X.
93: *> \endverbatim
94: *>
95: *> \param[in] LDB
96: *> \verbatim
97: *> LDB is INTEGER
98: *> The leading dimension of the array B. LDB >= max(1,N).
99: *> \endverbatim
100: *
101: * Authors:
102: * ========
103: *
104: *> \author Univ. of Tennessee
105: *> \author Univ. of California Berkeley
106: *> \author Univ. of Colorado Denver
107: *> \author NAG Ltd.
108: *
109: *> \ingroup complex16PTcomputational
110: *
111: * =====================================================================
112: SUBROUTINE ZPTTS2( IUPLO, N, NRHS, D, E, B, LDB )
113: *
114: * -- LAPACK computational routine --
115: * -- LAPACK is a software package provided by Univ. of Tennessee, --
116: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
117: *
118: * .. Scalar Arguments ..
119: INTEGER IUPLO, LDB, N, NRHS
120: * ..
121: * .. Array Arguments ..
122: DOUBLE PRECISION D( * )
123: COMPLEX*16 B( LDB, * ), E( * )
124: * ..
125: *
126: * =====================================================================
127: *
128: * .. Local Scalars ..
129: INTEGER I, J
130: * ..
131: * .. External Subroutines ..
132: EXTERNAL ZDSCAL
133: * ..
134: * .. Intrinsic Functions ..
135: INTRINSIC DCONJG
136: * ..
137: * .. Executable Statements ..
138: *
139: * Quick return if possible
140: *
141: IF( N.LE.1 ) THEN
142: IF( N.EQ.1 )
143: $ CALL ZDSCAL( NRHS, 1.D0 / D( 1 ), B, LDB )
144: RETURN
145: END IF
146: *
147: IF( IUPLO.EQ.1 ) THEN
148: *
149: * Solve A * X = B using the factorization A = U**H *D*U,
150: * overwriting each right hand side vector with its solution.
151: *
152: IF( NRHS.LE.2 ) THEN
153: J = 1
154: 10 CONTINUE
155: *
156: * Solve U**H * x = b.
157: *
158: DO 20 I = 2, N
159: B( I, J ) = B( I, J ) - B( I-1, J )*DCONJG( E( I-1 ) )
160: 20 CONTINUE
161: *
162: * Solve D * U * x = b.
163: *
164: DO 30 I = 1, N
165: B( I, J ) = B( I, J ) / D( I )
166: 30 CONTINUE
167: DO 40 I = N - 1, 1, -1
168: B( I, J ) = B( I, J ) - B( I+1, J )*E( I )
169: 40 CONTINUE
170: IF( J.LT.NRHS ) THEN
171: J = J + 1
172: GO TO 10
173: END IF
174: ELSE
175: DO 70 J = 1, NRHS
176: *
177: * Solve U**H * x = b.
178: *
179: DO 50 I = 2, N
180: B( I, J ) = B( I, J ) - B( I-1, J )*DCONJG( E( I-1 ) )
181: 50 CONTINUE
182: *
183: * Solve D * U * x = b.
184: *
185: B( N, J ) = B( N, J ) / D( N )
186: DO 60 I = N - 1, 1, -1
187: B( I, J ) = B( I, J ) / D( I ) - B( I+1, J )*E( I )
188: 60 CONTINUE
189: 70 CONTINUE
190: END IF
191: ELSE
192: *
193: * Solve A * X = B using the factorization A = L*D*L**H,
194: * overwriting each right hand side vector with its solution.
195: *
196: IF( NRHS.LE.2 ) THEN
197: J = 1
198: 80 CONTINUE
199: *
200: * Solve L * x = b.
201: *
202: DO 90 I = 2, N
203: B( I, J ) = B( I, J ) - B( I-1, J )*E( I-1 )
204: 90 CONTINUE
205: *
206: * Solve D * L**H * x = b.
207: *
208: DO 100 I = 1, N
209: B( I, J ) = B( I, J ) / D( I )
210: 100 CONTINUE
211: DO 110 I = N - 1, 1, -1
212: B( I, J ) = B( I, J ) - B( I+1, J )*DCONJG( E( I ) )
213: 110 CONTINUE
214: IF( J.LT.NRHS ) THEN
215: J = J + 1
216: GO TO 80
217: END IF
218: ELSE
219: DO 140 J = 1, NRHS
220: *
221: * Solve L * x = b.
222: *
223: DO 120 I = 2, N
224: B( I, J ) = B( I, J ) - B( I-1, J )*E( I-1 )
225: 120 CONTINUE
226: *
227: * Solve D * L**H * x = b.
228: *
229: B( N, J ) = B( N, J ) / D( N )
230: DO 130 I = N - 1, 1, -1
231: B( I, J ) = B( I, J ) / D( I ) -
232: $ B( I+1, J )*DCONJG( E( I ) )
233: 130 CONTINUE
234: 140 CONTINUE
235: END IF
236: END IF
237: *
238: RETURN
239: *
240: * End of ZPTTS2
241: *
242: END
CVSweb interface <joel.bertrand@systella.fr>