LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ dtptri()

subroutine dtptri ( character  UPLO,
character  DIAG,
integer  N,
double precision, dimension( * )  AP,
integer  INFO 
)

DTPTRI

Download DTPTRI + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DTPTRI computes the inverse of a real upper or lower triangular
 matrix A stored in packed format.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  A is upper triangular;
          = 'L':  A is lower triangular.
[in]DIAG
          DIAG is CHARACTER*1
          = 'N':  A is non-unit triangular;
          = 'U':  A is unit triangular.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]AP
          AP is DOUBLE PRECISION array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangular matrix A, stored
          columnwise in a linear array.  The j-th column of A is stored
          in the array AP as follows:
          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
          if UPLO = 'L', AP(i + (j-1)*((2*n-j)/2) = A(i,j) for j<=i<=n.
          See below for further details.
          On exit, the (triangular) inverse of the original matrix, in
          the same packed storage format.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, A(i,i) is exactly zero.  The triangular
                matrix is singular and its inverse can not be computed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
December 2016
Further Details:
  A triangular matrix A can be transferred to packed storage using one
  of the following program segments:

  UPLO = 'U':                      UPLO = 'L':

        JC = 1                           JC = 1
        DO 2 J = 1, N                    DO 2 J = 1, N
           DO 1 I = 1, J                    DO 1 I = J, N
              AP(JC+I-1) = A(I,J)              AP(JC+I-J) = A(I,J)
      1    CONTINUE                    1    CONTINUE
           JC = JC + J                      JC = JC + N - J + 1
      2 CONTINUE                       2 CONTINUE

Definition at line 119 of file dtptri.f.

119 *
120 * -- LAPACK computational routine (version 3.7.0) --
121 * -- LAPACK is a software package provided by Univ. of Tennessee, --
122 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
123 * December 2016
124 *
125 * .. Scalar Arguments ..
126  CHARACTER diag, uplo
127  INTEGER info, n
128 * ..
129 * .. Array Arguments ..
130  DOUBLE PRECISION ap( * )
131 * ..
132 *
133 * =====================================================================
134 *
135 * .. Parameters ..
136  DOUBLE PRECISION one, zero
137  parameter( one = 1.0d+0, zero = 0.0d+0 )
138 * ..
139 * .. Local Scalars ..
140  LOGICAL nounit, upper
141  INTEGER j, jc, jclast, jj
142  DOUBLE PRECISION ajj
143 * ..
144 * .. External Functions ..
145  LOGICAL lsame
146  EXTERNAL lsame
147 * ..
148 * .. External Subroutines ..
149  EXTERNAL dscal, dtpmv, xerbla
150 * ..
151 * .. Executable Statements ..
152 *
153 * Test the input parameters.
154 *
155  info = 0
156  upper = lsame( uplo, 'U' )
157  nounit = lsame( diag, 'N' )
158  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
159  info = -1
160  ELSE IF( .NOT.nounit .AND. .NOT.lsame( diag, 'U' ) ) THEN
161  info = -2
162  ELSE IF( n.LT.0 ) THEN
163  info = -3
164  END IF
165  IF( info.NE.0 ) THEN
166  CALL xerbla( 'DTPTRI', -info )
167  RETURN
168  END IF
169 *
170 * Check for singularity if non-unit.
171 *
172  IF( nounit ) THEN
173  IF( upper ) THEN
174  jj = 0
175  DO 10 info = 1, n
176  jj = jj + info
177  IF( ap( jj ).EQ.zero )
178  $ RETURN
179  10 CONTINUE
180  ELSE
181  jj = 1
182  DO 20 info = 1, n
183  IF( ap( jj ).EQ.zero )
184  $ RETURN
185  jj = jj + n - info + 1
186  20 CONTINUE
187  END IF
188  info = 0
189  END IF
190 *
191  IF( upper ) THEN
192 *
193 * Compute inverse of upper triangular matrix.
194 *
195  jc = 1
196  DO 30 j = 1, n
197  IF( nounit ) THEN
198  ap( jc+j-1 ) = one / ap( jc+j-1 )
199  ajj = -ap( jc+j-1 )
200  ELSE
201  ajj = -one
202  END IF
203 *
204 * Compute elements 1:j-1 of j-th column.
205 *
206  CALL dtpmv( 'Upper', 'No transpose', diag, j-1, ap,
207  $ ap( jc ), 1 )
208  CALL dscal( j-1, ajj, ap( jc ), 1 )
209  jc = jc + j
210  30 CONTINUE
211 *
212  ELSE
213 *
214 * Compute inverse of lower triangular matrix.
215 *
216  jc = n*( n+1 ) / 2
217  DO 40 j = n, 1, -1
218  IF( nounit ) THEN
219  ap( jc ) = one / ap( jc )
220  ajj = -ap( jc )
221  ELSE
222  ajj = -one
223  END IF
224  IF( j.LT.n ) THEN
225 *
226 * Compute elements j+1:n of j-th column.
227 *
228  CALL dtpmv( 'Lower', 'No transpose', diag, n-j,
229  $ ap( jclast ), ap( jc+1 ), 1 )
230  CALL dscal( n-j, ajj, ap( jc+1 ), 1 )
231  END IF
232  jclast = jc
233  jc = jc - n + j - 2
234  40 CONTINUE
235  END IF
236 *
237  RETURN
238 *
239 * End of DTPTRI
240 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:81
subroutine dtpmv(UPLO, TRANS, DIAG, N, AP, X, INCX)
DTPMV
Definition: dtpmv.f:144
Here is the call graph for this function:
Here is the caller graph for this function: