MODLIB (CCP4: Library)

NAME

modlib - Subroutine Library for handling mathematical operations

DESCRIPTION

The library modlib contains subroutines for a number of mathematical operations, mainly operations with vectors and matrices. This document is the first stage in an attempt to document and rationalise the mathematical subroutines found in the CCP4 libraries and programs.

The subroutines have been arranged as much as possible into groups with similar functions or origins.

Overview of the subroutines in modlib

Vector operations
Matrix-vector operations
Matrix operations
Harwell routines
Miscellaneous routines

Information on BLAS and LAPACK routines

Introduction
BLAS
LAPACK
BLAS and LAPACK in CCP4
Installation issues for non-IEEE compliant machines

List of subroutines in modlib

Overview of the MODLIB routines

The routines are divided into four sections below: vector operations, matrix operations, Harwell routines, and miscellaneous (i.e. those which don't fit in anywhere else).

a. Vector Operations

Type Name Function Call using Arguments
s/r cross compute vector product A = B x C CROSS(A,B,C) REAL: A(3), B(3), C(3)
fun dot dot product of two vectors A.B DOT(A,B) REAL: A(3), B(3)
s/r icross cross product (integer version) A = B x C ICROSS(A,B,C) INTEGER: A(3), B(3), C(3)
fun idot dot product (integer version) A.B IDOT(A,B) INTEGER: A(3), B(3)
s/r unit Vector V reduced to unit vector UNIT(V) REAL V(3)
s/r scalev Scale vector B with scalar X and put result in A SCALEV(A,X,B) REAL: A(3), B(3), X
s/r vsum Add two vectors A = B + C VSUM(A,B,C) REAL: A(3),B(3),C(3)
s/r vdif Subtract two vectors A = B - C VDIF(A,B,C) REAL: A(3),B(3),C(3)
s/r vset Copy a vector from B to A VSET(A,B) REAL A(3),B(3)

b. Matrix-Vector Operations

Type Name Function Call using Arguments
s/r matvec Post-multiply a 3x3 matrix by a vector V = AB MATVEC(V,A,B) REAL: A(3,3), B(3), V(3)
s/r imatvec Post-multiply a 3x3 matrix by a vector V=AB (integer version) IMATVEC(V,A,B) INTEGER A(3,3),B(3),V(3)
s/r matvc4 Matrix x Vector, A(3) = R(4x4) . B(3) (see notes) MATVC4(A,R,B) REAL A(3), R(4,4), B(3)
s/r transfrm Transform vector X(3) by quine matrix MAT(4,4), return transformed vector in X (see notes) TRANSFRM(X,MAT) REAL X(3),MAT(4,4)

Notes on Matrix-Vector operations:

The subroutines MATVC4 and TRANSFRM perform the same operation, and differ only in that in the first case the result vector is returned separately from the input vector, whilst in the second the result overwrites the input vector.

The transformation of the input 3-vector by the 4x4 matrix is effected by treating the input vector as a 4-vector, with 1 in the last position, and performing a standard matrix times vector operation.

4x4 matrices are commonly used in CCP4 programs to represent symmetry operations, which may have both rotational and translational components (screw operations). For more information on this, see the documentation for SYMLIB, in particular the section on How Symmetry Operations are Stored and Applied.

c. Matrix Operations

Type Name Function Call using Arguments
s/r minv3 Invert a general 3x3 matrix and return determinant in D
A=(B)^-1 MINV3(A,B,D) REAL: A(3,3), B(3,3), D
s/r iminv3 Invert a general 3x3 matrix and return determinant in D (integer version) A = (B)^-1 IMINV3(A,B,D) INTEGER: A(3,3),B(3,3), D
s/r minvn Invert a matrix A of order N (overwritten on output) and return determinant D MINVN(A,N,D,L,M) INTEGER: N
REAL: A(N,N), D, L(N), M(N)
L,M are work-vectors
s/r matmul Multiply two 3x3 matrices A = BC MATMUL(A,B,C) REAL: A(3,3),B(3,3),C(3,3)
s/r matmulnm Multiply NxM MXN matrices A = BC MATMULNM(N,M,A,B,C) INTEGER: N, M
REAL: A(N,N), B(N,M), C(M,N)
s/r matmulgen Generalised matrix multiplication subroutine
Multiplies a NbxMbc matrix (B) by a MbcXNc (C) matrix, so that A = BC MATMULGEN(Nb,Mbc,Nc,A,B,C) INTEGER: Nb, Mbc, Nc
REAL: A(Nb,Nc),B(Nb,Mbc),C(Mbc,Nc)
s/r transp Transpose a 3x3 matrix A = B^T TRANSP(A,B) REAL: A(3,3),B(3,3)
s/r gmprd general matrix product R(N,L) = A(N,M) * B(M,L) GMPRD(A,B,R,N,M,L) REAL A(N*M),B(M*L),R(N*L)
INTEGER L,M,N
s/r matmul4 Multiply two 4x4 matrices, A=B*C MATMUL4(A,B,C) REAL A(4,4),B(4,4),C(4,4)
s/r matmln Multiply two nxn matrices, a = b . c MATMLN(N,A,B,C) INTEGER N
REAL A(N,N),B(N,N),C(N,N)
s/r detmat Calculate determinant DET of 3x3 matrix RMAT DETMAT(RMAT,DET) REAL DET,RMAT(3,3)
s/r ml3mat Multiply together three matrices of any size D = A.B.C ML3MAT(P,A,Q,B,R,C,S,D) REAL A(P,Q),B(Q,R),C(R,S),D(P,S)
INTEGER P,Q,R,S
s/r inv44 4x4 matrix inversion, returns IA as inverse of A INV44(A,AI) REAL A(4,4),AI(4,4)
s/r matmli Multiply two 3x3 matrices A = BC (integer version) MATMLI(A,B,C) INTEGER A(3,3),B(3,3),C(3,3)
s/r matmultrans A=B*C^T for 3x3 matrices MATMULTrans(A,B,C) REAL A(3,3),B(3,3),C(3,3)
s/r eigen_rs_asc eigenvalues and eigenvectors of real symmetric matrix EIGEN_RS_ASC(A, R, N, MV) REAL A(*), R(*)
INTEGER N, MV

d. Harwell Subroutines

The Harwell Subroutine Library was started at Harwell by the precursor of the Numerical Analysis Group at Rutherford Appleton Laboratory in 1963, to provide high quality numerical software to the scientists and engineers working there. Since then it has developed into a portable fully documented Fortran library. See http://www.cse.clrc.ac.uk/Activity/HSL for more details.

Some of the Harwell routines are already been distributed as part of CCP4, in modlib.f (ea06c, ea08c, ea09c, fa01as, fa01bs, fa01cs, fa01ds, fm02ad, mc04b)

The modlib header has the following comment on the conditions of use of the routines:

The routines ea06c, ea08c, ea09c, fa01as, fa01bs, fa01cs, fa01ds,
fm02ad, mc04b, (and possibly others) are from the
Harwell Subroutine library.  The conditions on their external use,
reproduced from the Harwell manual are:
* due acknowledgement is made of the use of subroutines in any
research publications resulting from their use.
* the subroutines may be modified for use in research applications
by external users.  The nature of such modifiactions should be
indicated in writing for information to the liaison officer.  At
no time however, shall the subroutines or modifications thereof
become the property of the external user.
The liaison officer for the library's external affairs is listed
as: Mr. S. Marlow, Building 8.9, Harwell Laboratory, Didcot, 
Oxon OX11 0RA, UK.

Information on some of the Harwell routines can also be found at http://www.netlib.org/harwell/index.html.

Type Name Function
s/r ea06c Given a real MxM symmetric matrix A = {aij} this routine finds all its eigenvalues (lambda)i i=1,2,.....,m and eigenvectors xj i=1,2,...,m. i.e. finds the non-trivial solutions of Ax=(lambda)x
s/r ea08c This uses QR iteration to find the all the eigenvalues and eigenvectors of a real symmetric tri-diagonal matrix
s/r ea09c Purpose is unspecified (no comments in header)
fun fa01as Purpose is unspecified (no comments in header)
s/r fa01bs Purpose is unspecified (no comments in header)
s/r fa01cs Purpose is unspecified (no comments in header)
s/r fa01ds Purpose is unspecified (no comments in header)
fun fm02ad Compute the inner product of two double precision real vectors accumulating the result double precision, when the elements of each vector are stored at some fixed displacement from neighbouring elements.
s/r mc04b Transforms a real symmetric matrix A={a(i,j)}, i, j=1..IA into a tri-diagonal matrix having the same eigenvalues as A using Householder's method.

e. Miscellaneous

Type Name Function Call using Arguments
s/r ranmar Universal random number generator; returns a random vector RVEC(LEN) RANMAR(RVEC,LEN) INTEGER: LEN
REAL: RVEC(LEN)
ent rmarin initialisation for ranmar RMARIN(IJ,KL) INTEGER: IJ, KL (seeds for ranmar)
s/r zipin Fast binary read on unit ID into real array BUF of length N ZIPIN(ID,N,BUF) INTEGER: ID, N
REAL: BUF(N)
s/r zipout Fast binary write to unit ID of real array BUF length N ZIPOUT(ID,N,BUF) INTEGER: ID, N
REAL: BUF(N)

2. Information on BLAS and LAPACK routines

(i) Introduction

From version 4.2, installation of the CCP4 suite is configured by default to include the BLAS and LAPACK libraries, which are used by some of the programs (e.g. SCALA, REFMAC5, BEAST). The advantage of using these external libraries is their robustness and efficiency of implementation, which should give greater reliability and higher speed.

(ii) BLAS = Basic Linear Algebra Sub-programs

BLAS contains routines for vector-vector, vector-matrix and matrix-matrix operations. Pre-compiled optimised versions of BLAS are distributed and should be used in preference to the versions available from NetLib, since they should be considerably faster.

The NetLib BLAS FAQ is at http://www.netlib.org/blas/faq.html which includes a list of links to vendor BLAS libraries for various systems, and to the BLAS Quick Reference Guide.

For systems where optimised ("vendor") BLAS libraries are unavailable, investigation of the Automatically Tuned Linear Algebra Software project (ATLAS) is recommended for improved BLAS performance- see http://math-atlas.sourceforge.net/ for more information. ATLAS can be also be used on systems which have a vendor BLAS library provided.

(iii) LAPACK = Linear Algebra PACKage

LAPACK is built on the BLAS routines and provides routines for solving systems of simultaneous equations, least-squares solutions of linear systems of equations, eigenvalue problems and singular value problems (plus others).

The NetLib LAPACK FAQ is at http://www.netlib.org/lapack/faq.html

There is also a html version of the LAPACK User Guide.

(iv) BLAS and LAPACK in CCP4

Installation of the CCP4 suite by default includes the BLAS and LAPACK libraries. This section describes how this done, and may be of interest if you have problems. LAPACK configuration can be disabled at build time by running the main CCP4 configure with the --disable-lapack option, but this is not advised as some programs will fail to build as a result.

Many systems carry their own version of the BLAS and/or LAPACK libraries, and configure will try to use these preferentially. Some examples for common systems are given in the table:

System	LAPACK and BLAS	BLAS only	See also:
OSF1 (DEC/Compaq)	[none known]	http://www5.compaq.com/math/documentation/index.html
IRIX (SGI)	libscs	libblas	http://www.sgi.com/software/scsl.html
HPUX	HP mlib	libblas	http://www.hp.com/rsn/mlib/mlibhome.html
SUNOS	Sun Performance Library	[none known]	http://docs.sun.com/htmlcoll/coll.118.3/iso-8859-1/PERFLIBUG/plug_bookTOC.html

If no system LAPACK library is found then configure will set up to build the `reference' version from Fortran source code originally taken from the Netlib archive. Source code for both BLAS and LAPACK v3.0 is included in the $CCP4/lib/lapack directory. Some comments:

The Netlib LAPACK routines will normally only be built if configure cannot find a system LAPACK library.
The NetLib BLAS routines will only be built if configure cannot detect a native (i.e. ``vendor'') BLAS library on the system. Vendor BLAS are used preferentially over the NetLib versions because they are optimised for the specific platform.
Only the single (real), double precision and complex versions of the routines are built, although the code is also supplied for complex*16 versions. To build these additional routines you will need to edit the blaslib and lapacklib targets in $CCP4/lib/lapack/Makefile.in, to add complex16, before re-running configure.
If you are installing on a non-IEEE compliant machine then there is a small but very important change which must be made to prevent run-time crashes; see below. Configure should normally make the change automatically, but if you encounter problems then please let us know.

Note that building of the Netlib routines can also be forced by specifying the --with-netlib-lapack when running configure.

Since the configuration and build has only been tested on a limited number of platforms, CCP4 welcomes any feedback or information about the LAPACK libraries under the Suite. If you have queries or comments please e-mail then to ccp4@ccp4.ac.uk

(v) Installation issues for non-IEEE compliant machines

By default the NetLib version of the LAPACK routine ILAENV (distributed with CCP4) assumes an IEEE-compliant machine, which implies certain behaviour for NaN and "infinity arithmetic". On non-IEEE compliant machines or O/S the LAPACK library will compile successfully but may crash at run-time when ILAENV is called.

The CCP4 distribution contains two versions of the ILAENV routine, in $CCP4/lib/lapack/src/. ilaenv-dist.f is the standard Netlib version of the routine, whilst ilaenv-non-ieee.f is the version modified for a non-IEEE compliant machine.

The CCP4 configure will copy the appropriate version of the routine to ilaenv.f on installation.

3. List of modlib routines:

SUBROUTINE CROSS(A,B,C)

Compute vector product A = B x C

     .. Array Arguments ..
      REAL             A(3),B(3),C(3)

REAL FUNCTION DOT(A,B)

dot product of two vectors

     .. Array Arguments ..
      REAL              A(3),B(3)

SUBROUTINE EA06C(A,VALUE,VECTOR,M,IA,IV,W)

C
C**  18/03/70 LAST LIBRARY UPDATE
C
C       ( Calls EA08C(W,W(M1),VALUE,VECTOR,M,IV,W(M+M1))
C       and MC04B(A,W,W(M1),M,IA,W(M+M1)) )
C
C       Given a real MxM symmetric matrix A = {aij} this routine
C       finds all its eigenvalues (lambda)i i=1,2,.....,m  and
C       eigenvectors xj i=1,2,...,m.  i.e. finds the non-trivial
C       solutions of Ax=(lambda)x
C
C       The matrix is reduced to tri-diagonal form by applying
C       Householder transformations. The eigenvalue problem for
C       the reduced problem is then solved using the QR algorithm
C       by calling EA08C.
C
C  Argument list
C  -------------
C
C   IA    (I) (integer) should be set to the first dimension
C                       of the array A, i.e. if the allocation
C                       for the array A was specified by
C                       DIMENSION A(100,50)
C                       then IA would be set to 100
C
C   M     (I) (integer) should be set to the order m of the matrix
C
C   IV    (I) (integer) should be set to the first dimension
C                       of the 2-dimensional array VECTOR
C
C   VECTOR(IV,M) (O) (real) 2-dimensional array, with first dimension IV,
C                           containing the eigenvectors. The components
C                           of the eigenvector vector(i) corresponding
C                           to the eigenvalue (lambda)i (in VALUE(I))
C                           are placed in VECTOR(J,I) J=1,2,...,M.
C                           The eigenvectors are normalized so that
C                           xT(i)x(i)=1 i=1,2,...,m.
C
C   VALUE(M) (O) (real)  array in which the routine puts
C                        the eigenvalues (lambda)i, i=1,2,...,m.
C                        These are not necessarily in any order.
C
C   W    (I) (real(*)) working array used by the routine for
C                      work space. The dimension must be set
C                      to at least 5*M.

SUBROUTINE EA08C(A,B,VALUE,VEC,M,IV,W)

C       (Calls EA09C(A,B,W(M+1),M,W))
C
C       This uses QR iteration to find the all the eigenvalues and
C       eigenvectors of the real symmetric tri-diagonal matrix
C       whose diagonal elements are A(i), i=1,M and off-diagonal
C       elements are B(i),i=2,M. The eigenvalues will have unit
C       length. The array W is used for workspace and must have
C       dimension at least 2*M. We treat VEC as if it had
C       dimensions (IV,M).
C
C       First EA09, which uses the QR algorithm, is used to find
C       the eigenvalues; using these as shifts the QR algorithm is
C       again applied but now using the plane rotations to generate
C       the eigenvectors. Finally the eigenvalues are refined
C       by taking Rayleigh quotients of the vectors.
C
C    Argument list
C    -------------
C
C       A(M)    (I) (real)    Diagonal elements
C
C       B(M)    (I) (real)    Off-diagonal elements
C
C       IV      (I)  (integer)  should be set to the first dimen-
C                               sion of the 2-dimensional array VEC
C
C       M       (I) (integer) should be set to the order m of the
C                             matrix
C
C       VALUE(M) (O) (real)   Eigenvalues
C
C       VEC     (O) (real)    Eigenvectors. The dimensions
C                             should be set to (IV,M).
C
C       W(*)    (I) (real)    Working array.The dimension must be
C                             set to at least 2*M.

SUBROUTINE EA09C(A,B,VALUE,M,OFF)

C    18/03/70 LAST LIBRARY UPDATE

[No info]

FUNCTION FA01AS(I)

[No info]

SUBROUTINE FA01BS(MAX,NRAND)

[No info]

SUBROUTINE FA01CS(IL,IR)

[No info]

SUBROUTINE FA01CS(IL,IR)

[No info]

SUBROUTINE FA01DS(IL,IR)

[No info]

DOUBLE PRECISION FUNCTION FM02AD(N,A,IA,B,IB)

C     Compute the inner product of two double precision real
C     vectors accumulating the result double precision, when the
C     elements of each vector are stored at some fixed displacement
C     from neighbouring elements. Given vectors A={a(j)},
C     B={b(j)} of length N, evaluates w=a(j)b(j) summed over
C     j=1..N. Can be used to evaluate inner products involving
C     rows of multi-dimensional arrays.
C     It can be used as an alternative to the assembler version,
C     but note that it is likely to be significantly slower in execution.
C
C     Argument list
C     -------------
C
C       N  (I) (integer) The length of the vectors (if N <= 0 FM02AD = 0)
C
C       A  (I) (double precision) The first vector
C
C       IA (I) (integer) Subscript displacement between elements of A
C
C       B  (I) (double precision) The second vector
C
C       IB (I) (integer) Subscript displacement between elements of B
C
C       FM02AD  the result

SUBROUTINE ICROSS(A,B,C)

C    Cross product (integer version)
C
C    A = B x C
C
C     .. Array Arguments ..
      INTEGER           A(3),B(3),C(3)

INTEGER FUNCTION IDOT(A,B)

C      Dot product (integer version)
C
C      IDOT = A . B
C
C     .. Array Arguments ..
      INTEGER               A(3),B(3)

SUBROUTINE IMINV3(A,B,D)

C     Invert a general 3x3 matrix and return determinant in D
C     (integer version)
C
C     A = (B)-1
C
C     .. Scalar Arguments ..
      INTEGER          D
C     ..
C     .. Array Arguments ..
      INTEGER          A(3,3),B(3,3)

SUBROUTINE MATMUL(A,B,C)

C      Multiply two 3x3 matrices
C
C      A = BC
C
C     .. Array Arguments ..
      REAL              A(3,3),B(3,3),C(3,3)

SUBROUTINE MATMULNM(N,M,A,B,C)

C      Multiply  NxM  MXN matrices
C
C      A = BC
C
C     .. Array Arguments ..
      REAL              A(N,N),B(N,M),C(M,N)

SUBROUTINE MATVEC(V,A,B)

C      Post-multiply a 3x3 matrix by a vector
C
C      V = AB
C
C     .. Array Arguments ..
      REAL              A(3,3),B(3),V(3)

SUBROUTINE MATMULGEN(Nb,Mbc,Nc,A,B,C)

C      Generalised matrix multiplication subroutine
C      Multiplies a NbxMbc matrix (B) by a  MbcXNc
C      (C) matrix, so that
C
C      A = BC
C
      IMPLICIT NONE
C     ..
C     .. Scalar Arguments ..
      INTEGER           Nb,Mbc,Nc
C     ..
C     .. Array Arguments ..
      REAL              A(Nb,Nc),B(Nb,Mbc),C(Mbc,Nc)

SUBROUTINE MC04B(A,ALPHA,BETA,M,IA,Q)

C     Transforms a real symmetric matrix A={a(i,j)}, i, j=1..IA
C     into a tri-diagonal matrix having the same eigenvalues as A
C     using Householder's method.
C
C     .. Scalar Arguments ..
      INTEGER          IA,M
C     ..
C     .. Array Arguments ..
      REAL             A(IA,1),ALPHA(1),BETA(1),Q(1)

SUBROUTINE MINVN(A,N,D,L,M)

C---- Purpose
C     =======
C
C           invert a matrix
C
C---- Usage
C     ======
C
C           CALL MINVN(A,N,D,L,M)
C
C---- Description of parameters
C     =========================
C
C    A - input matrix, destroyed in computation and replaced by
C        resultant inverse.
C
C    N - order of matrix A
C
C    D - resultant determinant
C
C    L - work vector of length n
C
C    M - work vector of length n
C
C---- Remarks
C     =======
C
C     Matrix a must be a general matrix
C
C---- Subroutines and function subprograms required
C     =============================================
C
C           NONE
C
C---- Method
C     ======
C
C     The standard gauss-jordan method is used. the determinant
C     is also calculated. a determinant of zero indicates that
C     the matrix is singular.
C
C
C---- Note
C     =====
C
C     If a double precision version of this routine is desired, the
C     c in column 1 should be removed from the double precision
C     statement which follows.
C
C     double precision a,d,biga,hold
C
C        the c must also be removed from double precision statements
C        appearing in other routines used in conjunction with this
C        routine.
C
C        The double precision version of this subroutine must also
C        contain double precision fortran functions.  abs in statement
C        10 must be changed to dabs.
C
ccc     REAL*8 D

SUBROUTINE MINV3(A,B,D)

C     Invert a general 3x3 matrix and return determinant in D
C
C     A = (B)-1
C
C     .. Scalar Arguments ..
      REAL            D
C     ..
C     .. Array Arguments ..
      REAL            A(3,3),B(3,3)

SUBROUTINE RANMAR(RVEC,LEN)

C     Universal random number generator proposed by Marsaglia and Zaman
C     in report FSU-SCRI-87-50
C     slightly modified by F. James, 1988 to generate a vector
C     of pseudorandom numbers RVEC of length LEN
C     and making the COMMON block include everything needed to
C     specify completely the state of the generator.
C     Transcribed from CERN report DD/88/22.
C     Rather inelegant messing about added by D. Love, Jan. 1989 to
C     make sure initialisation always occurs.
C     *** James says that this is the preferred generator.
C     Gives bit-identical results on all machines with at least
C     24-bit mantissas in the flotaing point representation (i.e.
C     all common 32-bit computers. Fairly fast, satisfies very
C     stringent tests, has very long period and makes it very
C     simple to generate independly disjoint sequences.
C     See also RANECU.
C     The state of the generator may be saved/restored using the
C     whole contents of /RASET1/.
C     Call RANMAR to get a vector, RMARIN to initialise. 
C     
C  Argument list
C  -------------
C
C     VREC (O)                 (REAL)   Random Vector
C
C     LEN  (I)              (INTEGER)   Length of random vector
C
C
C  For ENTRY point RMARIN
C  ----------------------
C
C     Initialisation for RANMAR.  The input values should
C     be in the ranges: 0<=ij<=31328, 0<=kl<=30081
C     This shows the correspondence between the simplified input seeds
C     IJ, KL and the original Marsaglia-Zaman seeds i,j,k,l
C     To get standard values in Marsaglia-Zaman paper,
C     (I=12, J=34, K=56, L=78) put IJ=1802, KL=9373
C
C     IJ   (I)              (INTEGER)   Seed for random number generator
C
C     KL   (I)              (INTEGER)   Seed for randon number generator

SUBROUTINE SCALEV(A,X,B)

C     Scale vector B with scalar X and put result in A
C
C     .. Scalar Arguments ..
      REAL              X
C     ..
C     .. Array Arguments ..
      REAL              A(3),B(3)

SUBROUTINE TRANSP(A,B)

C---- Transpose a 3x3 matrix
C
C     A = BT
C
C     .. Array Arguments ..
      REAL              A(3,3),B(3,3)

SUBROUTINE UNIT(V)

C    Vector V reduced to unit vector
C
C     .. Array Arguments ..
      REAL            V(3)

SUBROUTINE VDIF(A,B,C)

C      Subtract two vectors
C
C      A = B - C
C
C     .. Array Arguments ..
      REAL            A(3),B(3),C(3)

SUBROUTINE VSET(A,B)

C      Copy a vector from B to A
C
C     .. Array Arguments ..
      REAL            A(3),B(3)

SUBROUTINE VSUM(A,B,C)

C     Add two vectors
C
C     A = B + C
C
C     .. Array Arguments ..
      REAL            A(3),B(3),C(3)

SUBROUTINE ZIPIN(ID,N,BUF)

C     Fast binary read on unit ID into real array BUF of length N
C
C     .. Scalar Arguments ..
      INTEGER          ID,N
C     ..
C     .. Array Arguments ..
      REAL             BUF(N)

SUBROUTINE ZIPOUT(ID,N,BUF)

C     Fast binary write to unit ID of real array BUF length N
C
C     .. Scalar Arguments ..
      INTEGER           ID,N
C     ..
C     .. Array Arguments ..
      REAL              BUF(N)

SUBROUTINE GMPRD(A,B,R,N,M,L)

C---- Ssp general matrix product
C
C   R(N,L) = A(N,M) * B(M,L)
C
C     .. Scalar Arguments ..
      INTEGER L,M,N
C     ..
C     .. Array Arguments ..
      REAL A(N*M),B(M*L),R(N*L)

SUBROUTINE MATMUL4(A,B,C)

C     Multiply two 4x4 matrices, A=B*C
C
C     .. Array Arguments ..
      REAL A(4,4),B(4,4),C(4,4)

subroutine matmln(n,a,b,c)

C Multiply two nxn matrices
C  a = b . c
C
      integer n
      real a(n,n),b(n,n),c(n,n)

SUBROUTINE DETMAT(RMAT,DET)

C---- Calculate determinant DET of 3x3 matrix RMAT
C
C
C     .. Scalar Arguments ..
      REAL DET
C     ..
C     .. Array Arguments ..
      REAL RMAT(3,3)

SUBROUTINE MATVC4(A,R,B)

C Matrix x Vector
C       A(3)  =  R(4x4) . B(3)
C
      REAL A(3), R(4,4), B(3)

SUBROUTINE ML3MAT(P,A,Q,B,R,C,S,D)

C 26-Nov-1988       J. W. Pflugrath           Cold Spring Harbor Laboratory
C    Edited to conform to Fortran 77.  Renamed from Multiply_3_matrices to
C    ML3MAT
C
C ==============================================================================
C
C! to multiply three matrices
C
      SUBROUTINE ML3MAT
C          ! input: 1st side of 1st matrix
     1   (P
C          ! input: first matrix
     2   ,A
C          ! input: 2nd side of 1st matrix & 1st side of 2nd matrix
     3   ,Q
C          ! input: second matrix
     4   ,B
C          ! input: 2nd side of 2nd matrix & 1st side of 3rd matrix
     5   ,R
C          ! input: third matrix
     6   ,C
C          ! input: 2nd side of 3rd matrix
     7   ,S
C          ! output: product matrix
     8   ,D)
C
CEE Multiplies three real matrices of any dimensions.  It is not optimised
C   for very large matrices.
C Multiply_3_matrices
C*** this routine is inefficient!
C Multiply_3_matrices
Created: 15-NOV-1985 D.J.Thomas,   MRC Laboratory of Molecular Biology,
C                                  Hills Road, Cambridge, CB2 2QH, England

SUBROUTINE INV44(A,AI)

C SUBROUTINE TO INVERT 4*4 MATRICES FOR CONVERSION BETWEEN
C FRACTIONAL AND ORTHOGONAL AXES
C PARAMETERS
C
C           A (I)   4*4 MATRIX TO BE INVERTED
C          AI (O)   INVERSE MATRIX

SUBROUTINE MATMLI(A,B,C)

C Integer matrix multiply
      INTEGER A(3,3),B(3,3),C(3,3)

SUBROUTINE MATMULTrans(A,B,C)

C  A=B*C(transpose) for 3x3 matrices
C
C     ..Array arguments
      REAL A(3,3),B(3,3),C(3,3)

SUBROUTINE IMATVEC(V,A,B)

C---- Post-multiply a 3x3 matrix by a vector
C
C       V=AB
C
C     .. Array Arguments ..
      INTEGER A(3,3),B(3),V(3)

SUBROUTINE TRANSFRM(X,MAT)

C  Transform vector X(3) by quine matrix MAT(4,4)
C  Return transformed vector in X.
C
C     ..Array arguments..
      REAL X(3),MAT(4,4)

SUBROUTINE EIGEN_RS_ASC(A, R, N, MV)


C---- SUBROUTINE TO COMPUTE EIGENVALUES & EIGENVECTORS OF A REAL
C---- SYMMETRIC MATRIX, FROM IBM SSP MANUAL (SEE P165).
C---- DESCRIPTION OF PARAMETERS -
C---- A - ORIGINAL MATRIX STORED COLUMNWISE AS UPPER TRIANGLE ONLY,
C---- I.E. "STORAGE MODE" = 1.  EIGENVALUES ARE WRITTEN INTO DIAGONAL
C---- ELEMENTS OF A  I.E.  A(1)  A(3)  A(6)  FOR A 3*3 MATRIX.
C---- R - RESULTANT MATRIX OF EIGENVECTORS STORED COLUMNWISE IN SAME
C---- ORDER AS EIGENVALUES.
C---- N - ORDER OF MATRICES A & R.
C---- MV = 0 TO COMPUTE EIGENVALUES & EIGENVECTORS.
C
      REAL        A(*), R(*)
      INTEGER     N, MV

MODLIB (CCP4: Library)

NAME

DESCRIPTION

CONTENTS

Overview of the MODLIB routines

a. Vector Operations

b. Matrix-Vector Operations

c. Matrix Operations

d. Harwell Subroutines

e. Miscellaneous

2. Information on BLAS and LAPACK routines

(i) Introduction

(ii) BLAS = Basic Linear Algebra Sub-programs

(iii) LAPACK = Linear Algebra PACKage

(iv) BLAS and LAPACK in CCP4

(v) Installation issues for non-IEEE compliant machines

3. List of modlib routines:

SUBROUTINE CROSS(A,B,C)

REAL FUNCTION DOT(A,B)

SUBROUTINE EA06C(A,VALUE,VECTOR,M,IA,IV,W)

SUBROUTINE EA08C(A,B,VALUE,VEC,M,IV,W)

SUBROUTINE EA09C(A,B,VALUE,M,OFF)

FUNCTION FA01AS(I)

SUBROUTINE FA01BS(MAX,NRAND)

SUBROUTINE FA01CS(IL,IR)

SUBROUTINE FA01CS(IL,IR)

SUBROUTINE FA01DS(IL,IR)

DOUBLE PRECISION FUNCTION FM02AD(N,A,IA,B,IB)

SUBROUTINE ICROSS(A,B,C)

INTEGER FUNCTION IDOT(A,B)

SUBROUTINE IMINV3(A,B,D)

SUBROUTINE MATMUL(A,B,C)

SUBROUTINE MATMULNM(N,M,A,B,C)

SUBROUTINE MATVEC(V,A,B)

SUBROUTINE MATMULGEN(Nb,Mbc,Nc,A,B,C)

SUBROUTINE MC04B(A,ALPHA,BETA,M,IA,Q)

SUBROUTINE MINVN(A,N,D,L,M)

SUBROUTINE MINV3(A,B,D)

SUBROUTINE RANMAR(RVEC,LEN)

SUBROUTINE SCALEV(A,X,B)

SUBROUTINE TRANSP(A,B)

SUBROUTINE UNIT(V)

SUBROUTINE VDIF(A,B,C)

SUBROUTINE VSET(A,B)

SUBROUTINE VSUM(A,B,C)

SUBROUTINE ZIPIN(ID,N,BUF)

SUBROUTINE ZIPOUT(ID,N,BUF)

SUBROUTINE GMPRD(A,B,R,N,M,L)

SUBROUTINE MATMUL4(A,B,C)

subroutine matmln(n,a,b,c)

SUBROUTINE DETMAT(RMAT,DET)

SUBROUTINE MATVC4(A,R,B)

SUBROUTINE ML3MAT(P,A,Q,B,R,C,S,D)

SUBROUTINE INV44(A,AI)

SUBROUTINE MATMLI(A,B,C)

SUBROUTINE MATMULTrans(A,B,C)

SUBROUTINE IMATVEC(V,A,B)

SUBROUTINE TRANSFRM(X,MAT)

SUBROUTINE EIGEN_RS_ASC(A, R, N, MV)