LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dtptri.f
Go to the documentation of this file.
1*> \brief \b DTPTRI
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DTPTRI + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtptri.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtptri.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtptri.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DTPTRI( UPLO, DIAG, N, AP, INFO )
22*
23* .. Scalar Arguments ..
24* CHARACTER DIAG, UPLO
25* INTEGER INFO, N
26* ..
27* .. Array Arguments ..
28* DOUBLE PRECISION AP( * )
29* ..
30*
31*
32*> \par Purpose:
33* =============
34*>
35*> \verbatim
36*>
37*> DTPTRI computes the inverse of a real upper or lower triangular
38*> matrix A stored in packed format.
39*> \endverbatim
40*
41* Arguments:
42* ==========
43*
44*> \param[in] UPLO
45*> \verbatim
46*> UPLO is CHARACTER*1
47*> = 'U': A is upper triangular;
48*> = 'L': A is lower triangular.
49*> \endverbatim
50*>
51*> \param[in] DIAG
52*> \verbatim
53*> DIAG is CHARACTER*1
54*> = 'N': A is non-unit triangular;
55*> = 'U': A is unit triangular.
56*> \endverbatim
57*>
58*> \param[in] N
59*> \verbatim
60*> N is INTEGER
61*> The order of the matrix A. N >= 0.
62*> \endverbatim
63*>
64*> \param[in,out] AP
65*> \verbatim
66*> AP is DOUBLE PRECISION array, dimension (N*(N+1)/2)
67*> On entry, the upper or lower triangular matrix A, stored
68*> columnwise in a linear array. The j-th column of A is stored
69*> in the array AP as follows:
70*> if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
71*> if UPLO = 'L', AP(i + (j-1)*((2*n-j)/2) = A(i,j) for j<=i<=n.
72*> See below for further details.
73*> On exit, the (triangular) inverse of the original matrix, in
74*> the same packed storage format.
75*> \endverbatim
76*>
77*> \param[out] INFO
78*> \verbatim
79*> INFO is INTEGER
80*> = 0: successful exit
81*> < 0: if INFO = -i, the i-th argument had an illegal value
82*> > 0: if INFO = i, A(i,i) is exactly zero. The triangular
83*> matrix is singular and its inverse can not be computed.
84*> \endverbatim
85*
86* Authors:
87* ========
88*
89*> \author Univ. of Tennessee
90*> \author Univ. of California Berkeley
91*> \author Univ. of Colorado Denver
92*> \author NAG Ltd.
93*
94*> \ingroup tptri
95*
96*> \par Further Details:
97* =====================
98*>
99*> \verbatim
100*>
101*> A triangular matrix A can be transferred to packed storage using one
102*> of the following program segments:
103*>
104*> UPLO = 'U': UPLO = 'L':
105*>
106*> JC = 1 JC = 1
107*> DO 2 J = 1, N DO 2 J = 1, N
108*> DO 1 I = 1, J DO 1 I = J, N
109*> AP(JC+I-1) = A(I,J) AP(JC+I-J) = A(I,J)
110*> 1 CONTINUE 1 CONTINUE
111*> JC = JC + J JC = JC + N - J + 1
112*> 2 CONTINUE 2 CONTINUE
113*> \endverbatim
114*>
115* =====================================================================
116 SUBROUTINE dtptri( UPLO, DIAG, N, AP, INFO )
117*
118* -- LAPACK computational routine --
119* -- LAPACK is a software package provided by Univ. of Tennessee, --
120* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
121*
122* .. Scalar Arguments ..
123 CHARACTER DIAG, UPLO
124 INTEGER INFO, N
125* ..
126* .. Array Arguments ..
127 DOUBLE PRECISION AP( * )
128* ..
129*
130* =====================================================================
131*
132* .. Parameters ..
133 DOUBLE PRECISION ONE, ZERO
134 parameter( one = 1.0d+0, zero = 0.0d+0 )
135* ..
136* .. Local Scalars ..
137 LOGICAL NOUNIT, UPPER
138 INTEGER J, JC, JCLAST, JJ
139 DOUBLE PRECISION AJJ
140* ..
141* .. External Functions ..
142 LOGICAL LSAME
143 EXTERNAL lsame
144* ..
145* .. External Subroutines ..
146 EXTERNAL dscal, dtpmv, xerbla
147* ..
148* .. Executable Statements ..
149*
150* Test the input parameters.
151*
152 info = 0
153 upper = lsame( uplo, 'U' )
154 nounit = lsame( diag, 'N' )
155 IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
156 info = -1
157 ELSE IF( .NOT.nounit .AND. .NOT.lsame( diag, 'U' ) ) THEN
158 info = -2
159 ELSE IF( n.LT.0 ) THEN
160 info = -3
161 END IF
162 IF( info.NE.0 ) THEN
163 CALL xerbla( 'DTPTRI', -info )
164 RETURN
165 END IF
166*
167* Check for singularity if non-unit.
168*
169 IF( nounit ) THEN
170 IF( upper ) THEN
171 jj = 0
172 DO 10 info = 1, n
173 jj = jj + info
174 IF( ap( jj ).EQ.zero )
175 $ RETURN
176 10 CONTINUE
177 ELSE
178 jj = 1
179 DO 20 info = 1, n
180 IF( ap( jj ).EQ.zero )
181 $ RETURN
182 jj = jj + n - info + 1
183 20 CONTINUE
184 END IF
185 info = 0
186 END IF
187*
188 IF( upper ) THEN
189*
190* Compute inverse of upper triangular matrix.
191*
192 jc = 1
193 DO 30 j = 1, n
194 IF( nounit ) THEN
195 ap( jc+j-1 ) = one / ap( jc+j-1 )
196 ajj = -ap( jc+j-1 )
197 ELSE
198 ajj = -one
199 END IF
200*
201* Compute elements 1:j-1 of j-th column.
202*
203 CALL dtpmv( 'Upper', 'No transpose', diag, j-1, ap,
204 $ ap( jc ), 1 )
205 CALL dscal( j-1, ajj, ap( jc ), 1 )
206 jc = jc + j
207 30 CONTINUE
208*
209 ELSE
210*
211* Compute inverse of lower triangular matrix.
212*
213 jc = n*( n+1 ) / 2
214 DO 40 j = n, 1, -1
215 IF( nounit ) THEN
216 ap( jc ) = one / ap( jc )
217 ajj = -ap( jc )
218 ELSE
219 ajj = -one
220 END IF
221 IF( j.LT.n ) THEN
222*
223* Compute elements j+1:n of j-th column.
224*
225 CALL dtpmv( 'Lower', 'No transpose', diag, n-j,
226 $ ap( jclast ), ap( jc+1 ), 1 )
227 CALL dscal( n-j, ajj, ap( jc+1 ), 1 )
228 END IF
229 jclast = jc
230 jc = jc - n + j - 2
231 40 CONTINUE
232 END IF
233*
234 RETURN
235*
236* End of DTPTRI
237*
238 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dscal(n, da, dx, incx)
DSCAL
Definition dscal.f:79
subroutine dtpmv(uplo, trans, diag, n, ap, x, incx)
DTPMV
Definition dtpmv.f:142
subroutine dtptri(uplo, diag, n, ap, info)
DTPTRI
Definition dtptri.f:117