Two of the great names in numerical methods – James Wilkinson
of the National Physical Laboratory, U.K., and C. Reinsch of the
Mathematical Institute of TH, Munchen, Germany – combined to assemble
a masterful suite of numerical procedures in
Handbook for Automatic Computation, Volume II,
Linear Algebra.
The algorithms cover linear equations, least squares, linear
programming, and eigenvalue and eigenvector problems.
Note: in some of the following procedures, the original comment
documentation and description of the parameters – which were
missing – have been restored. As well, the original procedure
names have been restored. The routines have not been retested.
Most of these routines have active links; the remainder
will be added as soon as they have been converted.
 Cholesky symmetric decomposition of a
positive definite matrix.
 Cholesky symmetric decomposition of a
positive definite matrix, using compact storage.
(Otherwise, it's the same as the previous item.)
 LDLt symmetric decomposition of a positive
definite matrix.
 Iterative refinement of the solution
of a positive definite system.
 Inversion of positive definite matrices by
the GaussJordan method. Two procedures are provided:
gjdef1 takes an n × n matrix as input. gjdef2
takes the lower triangle stored in compact form as an
n(n+1)/2 vector.
 Symmetric decomposition of positive definite
band matrices.
 The Conjugate gradient method
for solving a system of linear equations that have
symmetric positive definite coefficients.
 The solution of symmetric and unsymmetric
band equations & the calculation of eigenvectors of
band matrices.
 The solution of a real
system of linear equations by the Crout method.
 The solution of a complex
system of linear equations by the Crout method.
 Linear least squares solutions by
Householder transformations. This algorithm can also be
used for systems of linear equations of n equations in
n unknowns, and for matrix inversion.
 Elimination with weighted row
combinations.
Four procedures are provided for solving least squares problems,
and for solving n equations in n unknowns,
and for matrix inversion. Choose the one best suited to
your application:
ORTHOLIN1 is for solving for X a least
squares problem AX = B, where A is an n × m
matrix, X is an m × k matrix, and B is an
n × k matrix.
Orthlin2 is used when B is a single column.
ORTHO1 is used for matrix inversion.
ORTHO2 is for matrix inversion, and economises on storage.
 Singular value decomposition and
least squares solutions.
 Simplex method based on triangular
decompositions.
 The Jacobi method for a real
symmetric matrix converts it to tridiagonal form.
 Householder's tridiaagonalization
of a symmetric matrix.
Several algorithms are provided, some requiring the
matrix in compressed form, stored as a vector.

Eigenvalues and eigenvectors

IBM's Scientific Subroutine Library.
Most if not all of these
procedures are available. See
a list of their titles.
For the procedures, contact the undersigned.
STATISTICAL ALGORITHMS
These algorithms were originally published in
Journal of Royal Statistical Society (Series C):
Applied Statistics.
Alan Miller converted them to Fortran 90, and I have converted
those to PL/I.
 M. J. R. Healy's algorithm,
AS 6: Given a symmetric matrix a
of order n as lower triangle,
calculates an upper triangle, u, such that
uprime * u = a.
a must be positive semidefinite.
 P. R. Freeman's algorithm, AS 7:
Forms as a lower triangle, a generalised inverse
of the positive semidefinite symmetric matrix a of
order n, stored as lower triangle.
Tests both AS 6 and AS 7.
 AS 27: Calculate the upper tail area under
Student's tdistribution,
by G. A. R. Taylor.
 AS 60: Calculate the
eigenvalues and eigenvectors of a real symmetric matrix,
by D. N. Spanks & A. D. Dodd.
 AS 63: Log of the beta function (includes
log of the gamma function, by
K. L. Majumder & G. P. BhattaCharjee.
 AS 66:
Evaluates the tail area of the standardised normal curve
from x to infinity or from minus infinity to x,
by I. D. Hill.
 AS 91:
Evaluates the percentage points of the chisquared
probability distribution function, by
D. J. Best & D. E. Roberts.
 AS 110:
LpNORM Fit of straight line by extension of Schlossmacher,
particularly for 1 <= p <= 2, by
V. A. Sposito, W. J. Kennedy, and J. E. Gentle.
 J. Barnard's algorithm, AS 126:
Computes the probability of the normal range given t, the
upper limit of integration, and n, the sample size.
 R. D. Armstrong and M. K. Tung's algorithm,
AS 132:
Fit Y = ALPHA + BETA.X + error
 AS 135:
MinMax estimates for a linear multiple regression problem,
R. D. Armstrong and D. S. Tung.
 AS 136:
Divides M points in Ndimensional space into K clusters so that
the within cluster the sum of squares is minimized, by
J. A. Hartigan & M. A. Wong.
 R. E. Lund's algorithm, AS 152:
Cumulative hypergeometric probabilities
(Replaces AS 59 and AS 152, and incorporates AS R86.)
 AS 154:
Algorithm for exact maximum likelihood estimation of
autoregressivemoving average models by means of Kalman filtering,
by G. Gardner, A. C. Harvey, & G. D. A. Phillips.
 AS 155:
Distribution of a linear combination of noncentral chisquared
random variables, by R. B. Davies.
 AS 157:
The RunsUp and RunsDown Test.
 AS 177:
Exact calculation of Normal Scores.
 AS 181:
Calculates the ShapiroWilk W test and its significance level.
 AS 190:
Evaluates the probability from 0 to q for a studentized range
having v degrees of freedom and r samples.
 AS 192:
Computes approximate significance points of a Pearson curve with given
first four moments, or first three moments and left or right boundary.
 AS 205:
Enumerates all R*C contingency tables with given row totals N(I)
and column totals M(J) and calculates the hypergeometric
probability of each table.
 AS 207:
Fitting a generalized loglinear model
to fully or partially classified frequencies.
 AS 217:
Does the dip calculation for an ordered vector X using the
greatest convex minorant and the least concave majorant, skipping
through the data using the change points of these distributions.
 AS 227:
Generates all possible Nbit binary codes, and applies a users
procedure for each code generated.
 AS 239: Incomplete Gamma Integral.
 AS 241:
Produces the normal deviate Z corresponding to a given lower
tail area of P;
Z is accurate to about 1 part in 10**7.
 AS 245: Logarithm of the gamma function.
 AS 260:
Computes the C.D.F. for the distribution of the square of the
multiple correlation coefficient with parameters X, IP, N, and RHO2.
X is the sample value of R**2.
IP is the number of predictors, including 1 for the constant if one is
being fitted.
N is the number of cases.
RHO2 is the population value of the squared multiple correlation
coefficient (often set = 0).
 AS 261:
Computes the quantile of the distribution of the square of the
sample multiple correlation coefficient for given number of
random variables M, sample size SIZE, square of the population
multiple correlation coefficient RHO2, and lower tail area P.
 AS 275:
Computes the noncentral chisquare distribution function with positive
real degrees of freedom f and nonnegative noncentrality parameter theta.
 AS 282:
High breakdown regression and multivariate estimation,
by D. M. Hawkins.
Calculates the least median of squares regression, minimum volume
ellipsoid, and associated statistics.
 AS 285:
Finds the probability that a normally distributed random
Nvector with mean 0 and covariance COVAR falls
in area enclosed by the external userdefined function F,
S. L. Lohr.
 A. J. Miller's AS 290:
Generate a rectangular 2D grid of variance ratios from which to plot
confidence regions for two parameters using Halperin's method.
If the model is linear in all parameters other than the two selected,
the confidence regions are exact; otherwise they are approximate and the
user should test the sensitivity of the confidence regions to variation
in the other parameters.
 AS 295: Alan Miller's & Nam Nguyen's
Heuristic algorithm to pick N rows of X out of NCAND to maximize
the determinant of X'X, using the Fedorov exchange algorithm.
(This is a modified version of the published algorithm.)
 AS 298:
Hybrid minimization routine using simulated annealing,
by S. P. Brooks.
 AS 304:
Fisher's nonparametric randomization test for two small
independent random samples, by L. E. Richards and J. Byrd.
 AS 310:
Computes the cumulative distribution function of a
noncentral beta random variable, by R. Chattamvelli and R. Shanmugam.
 AS 319
Function minimization without derivatives.
R. A. Vowels, 25 February 2006. Updated 7 March 2006, 29 March 2006,
8 April 2006, 15 April 2006, 14 October 2006, 21 December 2006,
29 March 2007, 29 May 2007, 5 June 2007, 29 October 2007, 1 December 2007,
25 April 2008, 23 June 2008, 28 December 2008, 15 March 2009, 14 July 2009,
28 October 2009, 1 December 2009, 3 December 2010, 27 December 2010.