X-Git-Url: https://vcs.maemo.org/git/?a=blobdiff_plain;f=3rdparty%2Flapack%2Fslagtf.c;fp=3rdparty%2Flapack%2Fslagtf.c;h=d24f691831655ccef67e0ce840ed5dfca497894f;hb=e4c14cdbdf2fe805e79cd96ded236f57e7b89060;hp=0000000000000000000000000000000000000000;hpb=454138ff8a20f6edb9b65a910101403d8b520643;p=opencv diff --git a/3rdparty/lapack/slagtf.c b/3rdparty/lapack/slagtf.c new file mode 100644 index 0000000..d24f691 --- /dev/null +++ b/3rdparty/lapack/slagtf.c @@ -0,0 +1,210 @@ +#include "clapack.h" + +/* Subroutine */ int slagtf_(integer *n, real *a, real *lambda, real *b, real + *c__, real *tol, real *d__, integer *in, integer *info) +{ + /* System generated locals */ + integer i__1; + real r__1, r__2; + + /* Local variables */ + integer k; + real tl, eps, piv1, piv2, temp, mult, scale1, scale2; + extern doublereal slamch_(char *); + extern /* Subroutine */ int xerbla_(char *, integer *); + + +/* -- LAPACK routine (version 3.1) -- */ +/* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ +/* November 2006 */ + +/* .. Scalar Arguments .. */ +/* .. */ +/* .. Array Arguments .. */ +/* .. */ + +/* Purpose */ +/* ======= */ + +/* SLAGTF factorizes the matrix (T - lambda*I), where T is an n by n */ +/* tridiagonal matrix and lambda is a scalar, as */ + +/* T - lambda*I = PLU, */ + +/* where P is a permutation matrix, L is a unit lower tridiagonal matrix */ +/* with at most one non-zero sub-diagonal elements per column and U is */ +/* an upper triangular matrix with at most two non-zero super-diagonal */ +/* elements per column. */ + +/* The factorization is obtained by Gaussian elimination with partial */ +/* pivoting and implicit row scaling. */ + +/* The parameter LAMBDA is included in the routine so that SLAGTF may */ +/* be used, in conjunction with SLAGTS, to obtain eigenvectors of T by */ +/* inverse iteration. */ + +/* Arguments */ +/* ========= */ + +/* N (input) INTEGER */ +/* The order of the matrix T. */ + +/* A (input/output) REAL array, dimension (N) */ +/* On entry, A must contain the diagonal elements of T. */ + +/* On exit, A is overwritten by the n diagonal elements of the */ +/* upper triangular matrix U of the factorization of T. */ + +/* LAMBDA (input) REAL */ +/* On entry, the scalar lambda. */ + +/* B (input/output) REAL array, dimension (N-1) */ +/* On entry, B must contain the (n-1) super-diagonal elements of */ +/* T. */ + +/* On exit, B is overwritten by the (n-1) super-diagonal */ +/* elements of the matrix U of the factorization of T. */ + +/* C (input/output) REAL array, dimension (N-1) */ +/* On entry, C must contain the (n-1) sub-diagonal elements of */ +/* T. */ + +/* On exit, C is overwritten by the (n-1) sub-diagonal elements */ +/* of the matrix L of the factorization of T. */ + +/* TOL (input) REAL */ +/* On entry, a relative tolerance used to indicate whether or */ +/* not the matrix (T - lambda*I) is nearly singular. TOL should */ +/* normally be chose as approximately the largest relative error */ +/* in the elements of T. For example, if the elements of T are */ +/* correct to about 4 significant figures, then TOL should be */ +/* set to about 5*10**(-4). If TOL is supplied as less than eps, */ +/* where eps is the relative machine precision, then the value */ +/* eps is used in place of TOL. */ + +/* D (output) REAL array, dimension (N-2) */ +/* On exit, D is overwritten by the (n-2) second super-diagonal */ +/* elements of the matrix U of the factorization of T. */ + +/* IN (output) INTEGER array, dimension (N) */ +/* On exit, IN contains details of the permutation matrix P. If */ +/* an interchange occurred at the kth step of the elimination, */ +/* then IN(k) = 1, otherwise IN(k) = 0. The element IN(n) */ +/* returns the smallest positive integer j such that */ + +/* abs( u(j,j) ).le. norm( (T - lambda*I)(j) )*TOL, */ + +/* where norm( A(j) ) denotes the sum of the absolute values of */ +/* the jth row of the matrix A. If no such j exists then IN(n) */ +/* is returned as zero. If IN(n) is returned as positive, then a */ +/* diagonal element of U is small, indicating that */ +/* (T - lambda*I) is singular or nearly singular, */ + +/* INFO (output) INTEGER */ +/* = 0 : successful exit */ +/* .lt. 0: if INFO = -k, the kth argument had an illegal value */ + +/* ===================================================================== */ + +/* .. Parameters .. */ +/* .. */ +/* .. Local Scalars .. */ +/* .. */ +/* .. Intrinsic Functions .. */ +/* .. */ +/* .. External Functions .. */ +/* .. */ +/* .. External Subroutines .. */ +/* .. */ +/* .. Executable Statements .. */ + + /* Parameter adjustments */ + --in; + --d__; + --c__; + --b; + --a; + + /* Function Body */ + *info = 0; + if (*n < 0) { + *info = -1; + i__1 = -(*info); + xerbla_("SLAGTF", &i__1); + return 0; + } + + if (*n == 0) { + return 0; + } + + a[1] -= *lambda; + in[*n] = 0; + if (*n == 1) { + if (a[1] == 0.f) { + in[1] = 1; + } + return 0; + } + + eps = slamch_("Epsilon"); + + tl = dmax(*tol,eps); + scale1 = dabs(a[1]) + dabs(b[1]); + i__1 = *n - 1; + for (k = 1; k <= i__1; ++k) { + a[k + 1] -= *lambda; + scale2 = (r__1 = c__[k], dabs(r__1)) + (r__2 = a[k + 1], dabs(r__2)); + if (k < *n - 1) { + scale2 += (r__1 = b[k + 1], dabs(r__1)); + } + if (a[k] == 0.f) { + piv1 = 0.f; + } else { + piv1 = (r__1 = a[k], dabs(r__1)) / scale1; + } + if (c__[k] == 0.f) { + in[k] = 0; + piv2 = 0.f; + scale1 = scale2; + if (k < *n - 1) { + d__[k] = 0.f; + } + } else { + piv2 = (r__1 = c__[k], dabs(r__1)) / scale2; + if (piv2 <= piv1) { + in[k] = 0; + scale1 = scale2; + c__[k] /= a[k]; + a[k + 1] -= c__[k] * b[k]; + if (k < *n - 1) { + d__[k] = 0.f; + } + } else { + in[k] = 1; + mult = a[k] / c__[k]; + a[k] = c__[k]; + temp = a[k + 1]; + a[k + 1] = b[k] - mult * temp; + if (k < *n - 1) { + d__[k] = b[k + 1]; + b[k + 1] = -mult * d__[k]; + } + b[k] = temp; + c__[k] = mult; + } + } + if (dmax(piv1,piv2) <= tl && in[*n] == 0) { + in[*n] = k; + } +/* L10: */ + } + if ((r__1 = a[*n], dabs(r__1)) <= scale1 * tl && in[*n] == 0) { + in[*n] = *n; + } + + return 0; + +/* End of SLAGTF */ + +} /* slagtf_ */