|
5164
|
1 SUBROUTINE DPTTS2( N, NRHS, D, E, B, LDB ) |
|
|
2 * |
|
|
3 * -- LAPACK routine (version 3.0) -- |
|
|
4 * Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., |
|
|
5 * Courant Institute, Argonne National Lab, and Rice University |
|
|
6 * June 30, 1999 |
|
|
7 * |
|
|
8 * .. Scalar Arguments .. |
|
|
9 INTEGER LDB, N, NRHS |
|
|
10 * .. |
|
|
11 * .. Array Arguments .. |
|
|
12 DOUBLE PRECISION B( LDB, * ), D( * ), E( * ) |
|
|
13 * .. |
|
|
14 * |
|
|
15 * Purpose |
|
|
16 * ======= |
|
|
17 * |
|
|
18 * DPTTS2 solves a tridiagonal system of the form |
|
|
19 * A * X = B |
|
|
20 * using the L*D*L' factorization of A computed by DPTTRF. D is a |
|
|
21 * diagonal matrix specified in the vector D, L is a unit bidiagonal |
|
|
22 * matrix whose subdiagonal is specified in the vector E, and X and B |
|
|
23 * are N by NRHS matrices. |
|
|
24 * |
|
|
25 * Arguments |
|
|
26 * ========= |
|
|
27 * |
|
|
28 * N (input) INTEGER |
|
|
29 * The order of the tridiagonal matrix A. N >= 0. |
|
|
30 * |
|
|
31 * NRHS (input) INTEGER |
|
|
32 * The number of right hand sides, i.e., the number of columns |
|
|
33 * of the matrix B. NRHS >= 0. |
|
|
34 * |
|
|
35 * D (input) DOUBLE PRECISION array, dimension (N) |
|
|
36 * The n diagonal elements of the diagonal matrix D from the |
|
|
37 * L*D*L' factorization of A. |
|
|
38 * |
|
|
39 * E (input) DOUBLE PRECISION array, dimension (N-1) |
|
|
40 * The (n-1) subdiagonal elements of the unit bidiagonal factor |
|
|
41 * L from the L*D*L' factorization of A. E can also be regarded |
|
|
42 * as the superdiagonal of the unit bidiagonal factor U from the |
|
|
43 * factorization A = U'*D*U. |
|
|
44 * |
|
|
45 * B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS) |
|
|
46 * On entry, the right hand side vectors B for the system of |
|
|
47 * linear equations. |
|
|
48 * On exit, the solution vectors, X. |
|
|
49 * |
|
|
50 * LDB (input) INTEGER |
|
|
51 * The leading dimension of the array B. LDB >= max(1,N). |
|
|
52 * |
|
|
53 * ===================================================================== |
|
|
54 * |
|
|
55 * .. Local Scalars .. |
|
|
56 INTEGER I, J |
|
|
57 * .. |
|
|
58 * .. External Subroutines .. |
|
|
59 EXTERNAL DSCAL |
|
|
60 * .. |
|
|
61 * .. Executable Statements .. |
|
|
62 * |
|
|
63 * Quick return if possible |
|
|
64 * |
|
|
65 IF( N.LE.1 ) THEN |
|
|
66 IF( N.EQ.1 ) |
|
|
67 $ CALL DSCAL( NRHS, 1.D0 / D( 1 ), B, LDB ) |
|
|
68 RETURN |
|
|
69 END IF |
|
|
70 * |
|
|
71 * Solve A * X = B using the factorization A = L*D*L', |
|
|
72 * overwriting each right hand side vector with its solution. |
|
|
73 * |
|
|
74 DO 30 J = 1, NRHS |
|
|
75 * |
|
|
76 * Solve L * x = b. |
|
|
77 * |
|
|
78 DO 10 I = 2, N |
|
|
79 B( I, J ) = B( I, J ) - B( I-1, J )*E( I-1 ) |
|
|
80 10 CONTINUE |
|
|
81 * |
|
|
82 * Solve D * L' * x = b. |
|
|
83 * |
|
|
84 B( N, J ) = B( N, J ) / D( N ) |
|
|
85 DO 20 I = N - 1, 1, -1 |
|
|
86 B( I, J ) = B( I, J ) / D( I ) - B( I+1, J )*E( I ) |
|
|
87 20 CONTINUE |
|
|
88 30 CONTINUE |
|
|
89 * |
|
|
90 RETURN |
|
|
91 * |
|
|
92 * End of DPTTS2 |
|
|
93 * |
|
|
94 END |
