upgrade AST

1 files changed, 0 insertions, 218 deletions

diff --git a/ast/cminpack/qrsolv.c b/ast/cminpack/qrsolv.c
deleted file mode 100644
index 6ab9e98..0000000
--- a/ast/cminpack/qrsolv.c
+++ /dev/null
@@ -1,218 +0,0 @@
-#include "cminpack.h"
-#include <math.h>
-#include "cminpackP.h"
-
-__cminpack_attr__
-void __cminpack_func__(qrsolv)(int n, real *r, int ldr, 
-	const int *ipvt, const real *diag, const real *qtb, real *x, 
-	real *sdiag, real *wa)
-{
-    /* Initialized data */
-
-#define p5 .5
-#define p25 .25
-
-    /* Local variables */
-    int i, j, k, l;
-    real cos, sin, sum, temp;
-    int nsing;
-    real qtbpj;
-
-/*     ********** */
-
-/*     subroutine qrsolv */
-
-/*     given an m by n matrix a, an n by n diagonal matrix d, */
-/*     and an m-vector b, the problem is to determine an x which */
-/*     solves the system */
-
-/*           a*x = b ,     d*x = 0 , */
-
-/*     in the least squares sense. */
-
-/*     this subroutine completes the solution of the problem */
-/*     if it is provided with the necessary information from the */
-/*     qr factorization, with column pivoting, of a. that is, if */
-/*     a*p = q*r, where p is a permutation matrix, q has orthogonal */
-/*     columns, and r is an upper triangular matrix with diagonal */
-/*     elements of nonincreasing magnitude, then qrsolv expects */
-/*     the full upper triangle of r, the permutation matrix p, */
-/*     and the first n components of (q transpose)*b. the system */
-/*     a*x = b, d*x = 0, is then equivalent to */
-
-/*                  t       t */
-/*           r*z = q *b ,  p *d*p*z = 0 , */
-
-/*     where x = p*z. if this system does not have full rank, */
-/*     then a least squares solution is obtained. on output qrsolv */
-/*     also provides an upper triangular matrix s such that */
-
-/*            t   t               t */
-/*           p *(a *a + d*d)*p = s *s . */
-
-/*     s is computed within qrsolv and may be of separate interest. */
-
-/*     the subroutine statement is */
-
-/*       subroutine qrsolv(n,r,ldr,ipvt,diag,qtb,x,sdiag,wa) */
-
-/*     where */
-
-/*       n is a positive integer input variable set to the order of r. */
-
-/*       r is an n by n array. on input the full upper triangle */
-/*         must contain the full upper triangle of the matrix r. */
-/*         on output the full upper triangle is unaltered, and the */
-/*         strict lower triangle contains the strict upper triangle */
-/*         (transposed) of the upper triangular matrix s. */
-
-/*       ldr is a positive integer input variable not less than n */
-/*         which specifies the leading dimension of the array r. */
-
-/*       ipvt is an integer input array of length n which defines the */
-/*         permutation matrix p such that a*p = q*r. column j of p */
-/*         is column ipvt(j) of the identity matrix. */
-
-/*       diag is an input array of length n which must contain the */
-/*         diagonal elements of the matrix d. */
-
-/*       qtb is an input array of length n which must contain the first */
-/*         n elements of the vector (q transpose)*b. */
-
-/*       x is an output array of length n which contains the least */
-/*         squares solution of the system a*x = b, d*x = 0. */
-
-/*       sdiag is an output array of length n which contains the */
-/*         diagonal elements of the upper triangular matrix s. */
-
-/*       wa is a work array of length n. */
-
-/*     subprograms called */
-
-/*       fortran-supplied ... dabs,dsqrt */
-
-/*     argonne national laboratory. minpack project. march 1980. */
-/*     burton s. garbow, kenneth e. hillstrom, jorge j. more */
-
-/*     ********** */
-
-/*     copy r and (q transpose)*b to preserve input and initialize s. */
-/*     in particular, save the diagonal elements of r in x. */
-
-    for (j = 0; j < n; ++j) {
-	for (i = j; i < n; ++i) {
-	    r[i + j * ldr] = r[j + i * ldr];
-	}
-	x[j] = r[j + j * ldr];
-	wa[j] = qtb[j];
-    }
-
-/*     eliminate the diagonal matrix d using a givens rotation. */
-
-    for (j = 0; j < n; ++j) {
-
-/*        prepare the row of d to be eliminated, locating the */
-/*        diagonal element using p from the qr factorization. */
-
-	l = ipvt[j]-1;
-	if (diag[l] != 0.) {
-            for (k = j; k < n; ++k) {
-                sdiag[k] = 0.;
-            }
-            sdiag[j] = diag[l];
-
-/*        the transformations to eliminate the row of d */
-/*        modify only a single element of (q transpose)*b */
-/*        beyond the first n, which is initially zero. */
-
-            qtbpj = 0.;
-            for (k = j; k < n; ++k) {
-
-/*           determine a givens rotation which eliminates the */
-/*           appropriate element in the current row of d. */
-
-                if (sdiag[k] != 0.) {
-#                 ifdef USE_LAPACK
-                    dlartg_( &r[k + k * ldr], &sdiag[k], &cos, &sin, &temp );
-#                 else /* !USE_LAPACK */
-                    if (fabs(r[k + k * ldr]) < fabs(sdiag[k])) {
-                        real cotan;
-                        cotan = r[k + k * ldr] / sdiag[k];
-                        sin = p5 / sqrt(p25 + p25 * (cotan * cotan));
-                        cos = sin * cotan;
-                    } else {
-                        real tan;
-                        tan = sdiag[k] / r[k + k * ldr];
-                        cos = p5 / sqrt(p25 + p25 * (tan * tan));
-                        sin = cos * tan;
-                    }
-
-/*           compute the modified diagonal element of r and */
-/*           the modified element of ((q transpose)*b,0). */
-
-#                 endif /* !USE_LAPACK */
-                    temp = cos * wa[k] + sin * qtbpj;
-                    qtbpj = -sin * wa[k] + cos * qtbpj;
-                    wa[k] = temp;
-
-/*           accumulate the tranformation in the row of s. */
-#                 ifdef USE_CBLAS
-                    cblas_drot( n-k, &r[k + k * ldr], 1, &sdiag[k], 1, cos, sin );
-#                 else /* !USE_CBLAS */
-                    r[k + k * ldr] = cos * r[k + k * ldr] + sin * sdiag[k];
-                    if (n > k+1) {
-                        for (i = k+1; i < n; ++i) {
-                            temp = cos * r[i + k * ldr] + sin * sdiag[i];
-                            sdiag[i] = -sin * r[i + k * ldr] + cos * sdiag[i];
-                            r[i + k * ldr] = temp;
-                        }
-                    }
-#                 endif /* !USE_CBLAS */
-                }
-            }
-        }
-
-/*        store the diagonal element of s and restore */
-/*        the corresponding diagonal element of r. */
-
-	sdiag[j] = r[j + j * ldr];
-	r[j + j * ldr] = x[j];
-    }
-
-/*     solve the triangular system for z. if the system is */
-/*     singular, then obtain a least squares solution. */
-
-    nsing = n;
-    for (j = 0; j < n; ++j) {
-	if (sdiag[j] == 0. && nsing == n) {
-	    nsing = j;
-	}
-	if (nsing < n) {
-	    wa[j] = 0.;
-	}
-    }
-    if (nsing >= 1) {
-        for (k = 1; k <= nsing; ++k) {
-            j = nsing - k;
-            sum = 0.;
-            if (nsing > j+1) {
-                for (i = j+1; i < nsing; ++i) {
-                    sum += r[i + j * ldr] * wa[i];
-                }
-            }
-            wa[j] = (wa[j] - sum) / sdiag[j];
-        }
-    }
-
-/*     permute the components of z back to components of x. */
-
-    for (j = 0; j < n; ++j) {
-	l = ipvt[j]-1;
-	x[l] = wa[j];
-    }
-    return;
-
-/*     last card of subroutine qrsolv. */
-
-} /* qrsolv_ */
-