OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
double
LAPACK » General tridiagonal Matrix » Linear Solve

Functions

subroutine dgtsv (n, nrhs, dl, d, du, b, ldb, info)
  DGTSV computes the solution to system of linear equations A * X = B for GT matrices
subroutine dgtsvx (fact, trans, n, nrhs, dl, d, du, dlf, df, duf, du2, ipiv, b, ldb, x, ldx, rcond, ferr, berr, work, iwork, info)
  DGTSVX computes the solution to system of linear equations A * X = B for GT matrices

Detailed Description

This is the group of double solve driver functions for GT matrices

Function Documentation

◆ dgtsv()

subroutine dgtsv ( integer n,
integer nrhs,
double precision, dimension( * ) dl,
double precision, dimension( * ) d,
double precision, dimension( * ) du,
double precision, dimension( ldb, * ) b,
integer ldb,
integer info )

DGTSV computes the solution to system of linear equations A * X = B for GT matrices

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

Purpose:
!>
!> DGTSV  solves the equation
!>
!>    A*X = B,
!>
!> where A is an n by n tridiagonal matrix, by Gaussian elimination with
!> partial pivoting.
!>
!> Note that the equation  A**T*X = B  may be solved by interchanging the
!> order of the arguments DU and DL.
!>
Parameters
[in]N
!>          N is INTEGER
!>          The order of the matrix A.  N >= 0.
!>
[in]NRHS
!>          NRHS is INTEGER
!>          The number of right hand sides, i.e., the number of columns
!>          of the matrix B.  NRHS >= 0.
!>
[in,out]DL
!>          DL is DOUBLE PRECISION array, dimension (N-1)
!>          On entry, DL must contain the (n-1) sub-diagonal elements of
!>          A.
!>
!>          On exit, DL is overwritten by the (n-2) elements of the
!>          second super-diagonal of the upper triangular matrix U from
!>          the LU factorization of A, in DL(1), ..., DL(n-2).
!>
[in,out]D
!>          D is DOUBLE PRECISION array, dimension (N)
!>          On entry, D must contain the diagonal elements of A.
!>
!>          On exit, D is overwritten by the n diagonal elements of U.
!>
[in,out]DU
!>          DU is DOUBLE PRECISION array, dimension (N-1)
!>          On entry, DU must contain the (n-1) super-diagonal elements
!>          of A.
!>
!>          On exit, DU is overwritten by the (n-1) elements of the first
!>          super-diagonal of U.
!>
[in,out]B
!>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
!>          On entry, the N by NRHS matrix of right hand side matrix B.
!>          On exit, if INFO = 0, the N by NRHS solution matrix X.
!>
[in]LDB
!>          LDB is INTEGER
!>          The leading dimension of the array B.  LDB >= max(1,N).
!>
[out]INFO
!>          INFO is INTEGER
!>          = 0: successful exit
!>          < 0: if INFO = -i, the i-th argument had an illegal value
!>          > 0: if INFO = i, U(i,i) is exactly zero, and the solution
!>               has not been computed.  The factorization has not been
!>               completed unless i = N.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 126 of file dgtsv.f.

127*
128* -- LAPACK driver routine --
129* -- LAPACK is a software package provided by Univ. of Tennessee, --
130* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
131*
132* .. Scalar Arguments ..
133 INTEGER INFO, LDB, N, NRHS
134* ..
135* .. Array Arguments ..
136 DOUBLE PRECISION B( LDB, * ), D( * ), DL( * ), DU( * )
137* ..
138*
139* =====================================================================
140*
141* .. Parameters ..
142 DOUBLE PRECISION ZERO
143 parameter( zero = 0.0d+0 )
144* ..
145* .. Local Scalars ..
146 INTEGER I, J
147 DOUBLE PRECISION FACT, TEMP
148* ..
149* .. Intrinsic Functions ..
150 INTRINSIC abs, max
151* ..
152* .. External Subroutines ..
153 EXTERNAL xerbla
154* ..
155* .. Executable Statements ..
156*
157 info = 0
158 IF( n.LT.0 ) THEN
159 info = -1
160 ELSE IF( nrhs.LT.0 ) THEN
161 info = -2
162 ELSE IF( ldb.LT.max( 1, n ) ) THEN
163 info = -7
164 END IF
165 IF( info.NE.0 ) THEN
166 CALL xerbla( 'DGTSV ', -info )
167 RETURN
168 END IF
169*
170 IF( n.EQ.0 )
171 $ RETURN
172*
173 IF( nrhs.EQ.1 ) THEN
174 DO 10 i = 1, n - 2
175 IF( abs( d( i ) ).GE.abs( dl( i ) ) ) THEN
176*
177* No row interchange required
178*
179 IF( d( i ).NE.zero ) THEN
180 fact = dl( i ) / d( i )
181 d( i+1 ) = d( i+1 ) - fact*du( i )
182 b( i+1, 1 ) = b( i+1, 1 ) - fact*b( i, 1 )
183 ELSE
184 info = i
185 RETURN
186 END IF
187 dl( i ) = zero
188 ELSE
189*
190* Interchange rows I and I+1
191*
192 fact = d( i ) / dl( i )
193 d( i ) = dl( i )
194 temp = d( i+1 )
195 d( i+1 ) = du( i ) - fact*temp
196 dl( i ) = du( i+1 )
197 du( i+1 ) = -fact*dl( i )
198 du( i ) = temp
199 temp = b( i, 1 )
200 b( i, 1 ) = b( i+1, 1 )
201 b( i+1, 1 ) = temp - fact*b( i+1, 1 )
202 END IF
203 10 CONTINUE
204 IF( n.GT.1 ) THEN
205 i = n - 1
206 IF( abs( d( i ) ).GE.abs( dl( i ) ) ) THEN
207 IF( d( i ).NE.zero ) THEN
208 fact = dl( i ) / d( i )
209 d( i+1 ) = d( i+1 ) - fact*du( i )
210 b( i+1, 1 ) = b( i+1, 1 ) - fact*b( i, 1 )
211 ELSE
212 info = i
213 RETURN
214 END IF
215 ELSE
216 fact = d( i ) / dl( i )
217 d( i ) = dl( i )
218 temp = d( i+1 )
219 d( i+1 ) = du( i ) - fact*temp
220 du( i ) = temp
221 temp = b( i, 1 )
222 b( i, 1 ) = b( i+1, 1 )
223 b( i+1, 1 ) = temp - fact*b( i+1, 1 )
224 END IF
225 END IF
226 IF( d( n ).EQ.zero ) THEN
227 info = n
228 RETURN
229 END IF
230 ELSE
231 DO 40 i = 1, n - 2
232 IF( abs( d( i ) ).GE.abs( dl( i ) ) ) THEN
233*
234* No row interchange required
235*
236 IF( d( i ).NE.zero ) THEN
237 fact = dl( i ) / d( i )
238 d( i+1 ) = d( i+1 ) - fact*du( i )
239 DO 20 j = 1, nrhs
240 b( i+1, j ) = b( i+1, j ) - fact*b( i, j )
241 20 CONTINUE
242 ELSE
243 info = i
244 RETURN
245 END IF
246 dl( i ) = zero
247 ELSE
248*
249* Interchange rows I and I+1
250*
251 fact = d( i ) / dl( i )
252 d( i ) = dl( i )
253 temp = d( i+1 )
254 d( i+1 ) = du( i ) - fact*temp
255 dl( i ) = du( i+1 )
256 du( i+1 ) = -fact*dl( i )
257 du( i ) = temp
258 DO 30 j = 1, nrhs
259 temp = b( i, j )
260 b( i, j ) = b( i+1, j )
261 b( i+1, j ) = temp - fact*b( i+1, j )
262 30 CONTINUE
263 END IF
264 40 CONTINUE
265 IF( n.GT.1 ) THEN
266 i = n - 1
267 IF( abs( d( i ) ).GE.abs( dl( i ) ) ) THEN
268 IF( d( i ).NE.zero ) THEN
269 fact = dl( i ) / d( i )
270 d( i+1 ) = d( i+1 ) - fact*du( i )
271 DO 50 j = 1, nrhs
272 b( i+1, j ) = b( i+1, j ) - fact*b( i, j )
273 50 CONTINUE
274 ELSE
275 info = i
276 RETURN
277 END IF
278 ELSE
279 fact = d( i ) / dl( i )
280 d( i ) = dl( i )
281 temp = d( i+1 )
282 d( i+1 ) = du( i ) - fact*temp
283 du( i ) = temp
284 DO 60 j = 1, nrhs
285 temp = b( i, j )
286 b( i, j ) = b( i+1, j )
287 b( i+1, j ) = temp - fact*b( i+1, j )
288 60 CONTINUE
289 END IF
290 END IF
291 IF( d( n ).EQ.zero ) THEN
292 info = n
293 RETURN
294 END IF
295 END IF
296*
297* Back solve with the matrix U from the factorization.
298*
299 IF( nrhs.LE.2 ) THEN
300 j = 1
301 70 CONTINUE
302 b( n, j ) = b( n, j ) / d( n )
303 IF( n.GT.1 )
304 $ b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) / d( n-1 )
305 DO 80 i = n - 2, 1, -1
306 b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-dl( i )*
307 $ b( i+2, j ) ) / d( i )
308 80 CONTINUE
309 IF( j.LT.nrhs ) THEN
310 j = j + 1
311 GO TO 70
312 END IF
313 ELSE
314 DO 100 j = 1, nrhs
315 b( n, j ) = b( n, j ) / d( n )
316 IF( n.GT.1 )
317 $ b( n-1, j ) = ( b( n-1, j )-du( n-1 )*b( n, j ) ) /
318 $ d( n-1 )
319 DO 90 i = n - 2, 1, -1
320 b( i, j ) = ( b( i, j )-du( i )*b( i+1, j )-dl( i )*
321 $ b( i+2, j ) ) / d( i )
322 90 CONTINUE
323 100 CONTINUE
324 END IF
325*
326 RETURN
327*
328* End of DGTSV
329*
xerbla
subroutine xerbla(srname, info)
XERBLA
Definition xerbla.f:60
max
#define max(a, b)
Definition macros.h:21

◆ dgtsvx()

subroutine dgtsvx ( character fact,
character trans,
integer n,
integer nrhs,
double precision, dimension( * ) dl,
double precision, dimension( * ) d,
double precision, dimension( * ) du,
double precision, dimension( * ) dlf,
double precision, dimension( * ) df,
double precision, dimension( * ) duf,
double precision, dimension( * ) du2,
integer, dimension( * ) ipiv,
double precision, dimension( ldb, * ) b,
integer ldb,
double precision, dimension( ldx, * ) x,
integer ldx,
double precision rcond,
double precision, dimension( * ) ferr,
double precision, dimension( * ) berr,
double precision, dimension( * ) work,
integer, dimension( * ) iwork,
integer info )

DGTSVX computes the solution to system of linear equations A * X = B for GT matrices

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

Purpose:
!>
!> DGTSVX uses the LU factorization to compute the solution to a real
!> system of linear equations A * X = B or A**T * X = B,
!> where A is a tridiagonal matrix of order N and X and B are N-by-NRHS
!> matrices.
!>
!> Error bounds on the solution and a condition estimate are also
!> provided.
!>
Description:
!>
!> The following steps are performed:
!>
!> 1. If FACT = 'N', the LU decomposition is used to factor the matrix A
!>    as A = L * U, where L is a product of permutation and unit lower
!>    bidiagonal matrices and U is upper triangular with nonzeros in
!>    only the main diagonal and first two superdiagonals.
!>
!> 2. If some U(i,i)=0, so that U is exactly singular, then the routine
!>    returns with INFO = i. Otherwise, the factored form of A is used
!>    to estimate the condition number of the matrix A.  If the
!>    reciprocal of the condition number is less than machine precision,
!>    INFO = N+1 is returned as a warning, but the routine still goes on
!>    to solve for X and compute error bounds as described below.
!>
!> 3. The system of equations is solved for X using the factored form
!>    of A.
!>
!> 4. Iterative refinement is applied to improve the computed solution
!>    matrix and calculate error bounds and backward error estimates
!>    for it.
!>
Parameters
[in]FACT
!>          FACT is CHARACTER*1
!>          Specifies whether or not the factored form of A has been
!>          supplied on entry.
!>          = 'F':  DLF, DF, DUF, DU2, and IPIV contain the factored
!>                  form of A; DL, D, DU, DLF, DF, DUF, DU2 and IPIV
!>                  will not be modified.
!>          = 'N':  The matrix will be copied to DLF, DF, and DUF
!>                  and factored.
!>
[in]TRANS
!>          TRANS is CHARACTER*1
!>          Specifies the form of the system of equations:
!>          = 'N':  A * X = B     (No transpose)
!>          = 'T':  A**T * X = B  (Transpose)
!>          = 'C':  A**H * X = B  (Conjugate transpose = Transpose)
!>
[in]N
!>          N is INTEGER
!>          The order of the matrix A.  N >= 0.
!>
[in]NRHS
!>          NRHS is INTEGER
!>          The number of right hand sides, i.e., the number of columns
!>          of the matrix B.  NRHS >= 0.
!>
[in]DL
!>          DL is DOUBLE PRECISION array, dimension (N-1)
!>          The (n-1) subdiagonal elements of A.
!>
[in]D
!>          D is DOUBLE PRECISION array, dimension (N)
!>          The n diagonal elements of A.
!>
[in]DU
!>          DU is DOUBLE PRECISION array, dimension (N-1)
!>          The (n-1) superdiagonal elements of A.
!>
[in,out]DLF
!>          DLF is DOUBLE PRECISION array, dimension (N-1)
!>          If FACT = 'F', then DLF is an input argument and on entry
!>          contains the (n-1) multipliers that define the matrix L from
!>          the LU factorization of A as computed by DGTTRF.
!>
!>          If FACT = 'N', then DLF is an output argument and on exit
!>          contains the (n-1) multipliers that define the matrix L from
!>          the LU factorization of A.
!>
[in,out]DF
!>          DF is DOUBLE PRECISION array, dimension (N)
!>          If FACT = 'F', then DF is an input argument and on entry
!>          contains the n diagonal elements of the upper triangular
!>          matrix U from the LU factorization of A.
!>
!>          If FACT = 'N', then DF is an output argument and on exit
!>          contains the n diagonal elements of the upper triangular
!>          matrix U from the LU factorization of A.
!>
[in,out]DUF
!>          DUF is DOUBLE PRECISION array, dimension (N-1)
!>          If FACT = 'F', then DUF is an input argument and on entry
!>          contains the (n-1) elements of the first superdiagonal of U.
!>
!>          If FACT = 'N', then DUF is an output argument and on exit
!>          contains the (n-1) elements of the first superdiagonal of U.
!>
[in,out]DU2
!>          DU2 is DOUBLE PRECISION array, dimension (N-2)
!>          If FACT = 'F', then DU2 is an input argument and on entry
!>          contains the (n-2) elements of the second superdiagonal of
!>          U.
!>
!>          If FACT = 'N', then DU2 is an output argument and on exit
!>          contains the (n-2) elements of the second superdiagonal of
!>          U.
!>
[in,out]IPIV
!>          IPIV is INTEGER array, dimension (N)
!>          If FACT = 'F', then IPIV is an input argument and on entry
!>          contains the pivot indices from the LU factorization of A as
!>          computed by DGTTRF.
!>
!>          If FACT = 'N', then IPIV is an output argument and on exit
!>          contains the pivot indices from the LU factorization of A;
!>          row i of the matrix was interchanged with row IPIV(i).
!>          IPIV(i) will always be either i or i+1; IPIV(i) = i indicates
!>          a row interchange was not required.
!>
[in]B
!>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
!>          The N-by-NRHS right hand side matrix B.
!>
[in]LDB
!>          LDB is INTEGER
!>          The leading dimension of the array B.  LDB >= max(1,N).
!>
[out]X
!>          X is DOUBLE PRECISION array, dimension (LDX,NRHS)
!>          If INFO = 0 or INFO = N+1, the N-by-NRHS solution matrix X.
!>
[in]LDX
!>          LDX is INTEGER
!>          The leading dimension of the array X.  LDX >= max(1,N).
!>
[out]RCOND
!>          RCOND is DOUBLE PRECISION
!>          The estimate of the reciprocal condition number of the matrix
!>          A.  If RCOND is less than the machine precision (in
!>          particular, if RCOND = 0), the matrix is singular to working
!>          precision.  This condition is indicated by a return code of
!>          INFO > 0.
!>
[out]FERR
!>          FERR is DOUBLE PRECISION array, dimension (NRHS)
!>          The estimated forward error bound for each solution vector
!>          X(j) (the j-th column of the solution matrix X).
!>          If XTRUE is the true solution corresponding to X(j), FERR(j)
!>          is an estimated upper bound for the magnitude of the largest
!>          element in (X(j) - XTRUE) divided by the magnitude of the
!>          largest element in X(j).  The estimate is as reliable as
!>          the estimate for RCOND, and is almost always a slight
!>          overestimate of the true error.
!>
[out]BERR
!>          BERR is DOUBLE PRECISION array, dimension (NRHS)
!>          The componentwise relative backward error of each solution
!>          vector X(j) (i.e., the smallest relative change in
!>          any element of A or B that makes X(j) an exact solution).
!>
[out]WORK
!>          WORK is DOUBLE PRECISION array, dimension (3*N)
!>
[out]IWORK
!>          IWORK is INTEGER array, dimension (N)
!>
[out]INFO
!>          INFO is INTEGER
!>          = 0:  successful exit
!>          < 0:  if INFO = -i, the i-th argument had an illegal value
!>          > 0:  if INFO = i, and i is
!>                <= N:  U(i,i) is exactly zero.  The factorization
!>                       has not been completed unless i = N, but the
!>                       factor U is exactly singular, so the solution
!>                       and error bounds could not be computed.
!>                       RCOND = 0 is returned.
!>                = N+1: U is nonsingular, but RCOND is less than machine
!>                       precision, meaning that the matrix is singular
!>                       to working precision.  Nevertheless, the
!>                       solution and error bounds are computed because
!>                       there are a number of situations where the
!>                       computed solution can be more accurate than the
!>                       value of RCOND would suggest.
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 290 of file dgtsvx.f.

293*
294* -- LAPACK driver routine --
295* -- LAPACK is a software package provided by Univ. of Tennessee, --
296* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
297*
298* .. Scalar Arguments ..
299 CHARACTER FACT, TRANS
300 INTEGER INFO, LDB, LDX, N, NRHS
301 DOUBLE PRECISION RCOND
302* ..
303* .. Array Arguments ..
304 INTEGER IPIV( * ), IWORK( * )
305 DOUBLE PRECISION B( LDB, * ), BERR( * ), D( * ), DF( * ),
306 $ DL( * ), DLF( * ), DU( * ), DU2( * ), DUF( * ),
307 $ FERR( * ), WORK( * ), X( LDX, * )
308* ..
309*
310* =====================================================================
311*
312* .. Parameters ..
313 DOUBLE PRECISION ZERO
314 parameter( zero = 0.0d+0 )
315* ..
316* .. Local Scalars ..
317 LOGICAL NOFACT, NOTRAN
318 CHARACTER NORM
319 DOUBLE PRECISION ANORM
320* ..
321* .. External Functions ..
322 LOGICAL LSAME
323 DOUBLE PRECISION DLAMCH, DLANGT
324 EXTERNAL lsame, dlamch, dlangt
325* ..
326* .. External Subroutines ..
327 EXTERNAL dcopy, dgtcon, dgtrfs, dgttrf, dgttrs, dlacpy,
328 $ xerbla
329* ..
330* .. Intrinsic Functions ..
331 INTRINSIC max
332* ..
333* .. Executable Statements ..
334*
335 info = 0
336 nofact = lsame( fact, 'N' )
337 notran = lsame( trans, 'N' )
338 IF( .NOT.nofact .AND. .NOT.lsame( fact, 'F' ) ) THEN
339 info = -1
340 ELSE IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
341 $ lsame( trans, 'C' ) ) THEN
342 info = -2
343 ELSE IF( n.LT.0 ) THEN
344 info = -3
345 ELSE IF( nrhs.LT.0 ) THEN
346 info = -4
347 ELSE IF( ldb.LT.max( 1, n ) ) THEN
348 info = -14
349 ELSE IF( ldx.LT.max( 1, n ) ) THEN
350 info = -16
351 END IF
352 IF( info.NE.0 ) THEN
353 CALL xerbla( 'DGTSVX', -info )
354 RETURN
355 END IF
356*
357 IF( nofact ) THEN
358*
359* Compute the LU factorization of A.
360*
361 CALL dcopy( n, d, 1, df, 1 )
362 IF( n.GT.1 ) THEN
363 CALL dcopy( n-1, dl, 1, dlf, 1 )
364 CALL dcopy( n-1, du, 1, duf, 1 )
365 END IF
366 CALL dgttrf( n, dlf, df, duf, du2, ipiv, info )
367*
368* Return if INFO is non-zero.
369*
370 IF( info.GT.0 )THEN
371 rcond = zero
372 RETURN
373 END IF
374 END IF
375*
376* Compute the norm of the matrix A.
377*
378 IF( notran ) THEN
379 norm = '1'
380 ELSE
381 norm = 'I'
382 END IF
383 anorm = dlangt( norm, n, dl, d, du )
384*
385* Compute the reciprocal of the condition number of A.
386*
387 CALL dgtcon( norm, n, dlf, df, duf, du2, ipiv, anorm, rcond, work,
388 $ iwork, info )
389*
390* Compute the solution vectors X.
391*
392 CALL dlacpy( 'Full', n, nrhs, b, ldb, x, ldx )
393 CALL dgttrs( trans, n, nrhs, dlf, df, duf, du2, ipiv, x, ldx,
394 $ info )
395*
396* Use iterative refinement to improve the computed solutions and
397* compute error bounds and backward error estimates for them.
398*
399 CALL dgtrfs( trans, n, nrhs, dl, d, du, dlf, df, duf, du2, ipiv,
400 $ b, ldb, x, ldx, ferr, berr, work, iwork, info )
401*
402* Set INFO = N+1 if the matrix is singular to working precision.
403*
404 IF( rcond.LT.dlamch( 'Epsilon' ) )
405 $ info = n + 1
406*
407 RETURN
408*
409* End of DGTSVX
410*
norm
norm(diag(diag(diag(inv(mat))) -id.SOL), 2) % destroy mumps instance id.JOB
dlacpy
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:53
dgtcon
subroutine dgtcon(norm, n, dl, d, du, du2, ipiv, anorm, rcond, work, iwork, info)
DGTCON
Definition dgtcon.f:146
dgttrs
subroutine dgttrs(trans, n, nrhs, dl, d, du, du2, ipiv, b, ldb, info)
DGTTRS
Definition dgttrs.f:138
dgtrfs
subroutine dgtrfs(trans, n, nrhs, dl, d, du, dlf, df, duf, du2, ipiv, b, ldb, x, ldx, ferr, berr, work, iwork, info)
DGTRFS
Definition dgtrfs.f:209
dgttrf
subroutine dgttrf(n, dl, d, du, du2, ipiv, info)
DGTTRF
Definition dgttrf.f:124
dlangt
double precision function dlangt(norm, n, dl, d, du)
DLANGT returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition dlangt.f:106
dcopy
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
dlamch
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69