This is to announce version 4.1 of LinAlg, Linear Algebra and
Optimization classlib. The library:
- defines Matrix,Vector, subMatrices, and LAStreams over real domain;
- contains efficient and fool-proof implementations of level 1
and 2 BLAS (element-wise operations + various multiplications),
transposition, determinant evaluation and matrix inverse.
- There are operations on a single row/col/diagonal of a matrix.
- The package sports new concepts of Matrix views, various
Matrix streams, and a "new style" of returning matrices (via
LazyMatrix) and filling them out.
- The package also implements Brent's univariate optimization
and root finding, Hooke-Jeeves multidimensional optimization of
functors, and Aitken-Lagrange interpolation. The new version indeed
has the old body (interface) of the previous incarnation, but a
changed spirit.
LinAlg 4.1 has been redesigned following a grander vision of
doing Numerical Math in C++. You may find the vision rather
unconventional. It has lofty goals though: efficiency, safety, _and_
expressiveness. I think LinAlg 4.1 comes close.
IMHO, there ought to be a clean separation between storage
classes that own data, and view classes, which provide some particular
access to the data. For example, Matrix is a purely storage class. It
has rather few methods, pertaining mostly to allocating, disposing of,
resizing of a 2D grid of elements, answering a variety of related
queries, etc. All other LinAlg classes provide different views of/into
the grids: in 2D, in 1D (that is, a view of all matrix elements
linearly arranged in a column major order), a view of a single matrix
row, column, or a diagonal, a view of a certain matrix blocks,
etc. The views are designed to be as lightweight as possible: most of
the view functionality could be inlined. The views turned out to be so
flexible that many former Matrix methods can now be implemented
outside of the Matrix class. These methods can do their job without a
privileged access to Matrix' internals while not sacrificing
performance.
Another overarching idea is that of Streams. It has become my
conviction that many Numerical Math algorithms don't actually require
a direct access to matrix elements. Surprisingly many algorithms in
essence merely sweep along a matrix or some parts of it. This _serial_
access is faster than the random one: one does not need to bother
computing index expressions. The whole issue of index boundary
checking simply goes away: thus safety is _built-in_. To prove that it
is indeed possible to use the serial access in non-trivial
applications, I rewrote the computations of a determinant and of a
matrix inverse with only a sequential access to a matrix, _without
ever_ using a familiar 'A(i,j)'.
That's why the package offers a whole variety of streams,
striding the whole matrix or separate matrix rows, columns, the
diagonal in arbitrary steps. Streams can be user-paced or
library-paced. They offer a tremendous flexibility in accessing matrix
elements. One can specify exactly what he wants:
A *= B; to_every(A) *= of_every(B);
both statements do their own thing, as efficiently as
possible. Construction of streams does not require any heap storage,
and all operations on streams are in-lined (which permits a whole
range of compiler optimizations). It appears that it is possible
indeed to achieve safety, efficiency and expressiveness _all at the
same time_.
Lazy matrices also turned out more general than I originally
thought. You can write now
A = unit(B); Matrix D = transposed(B); D = inverse(B);
and even C = A * B; Matrix hht = haar * haar_t;
In all the cases, the result (of multiplication, transposition,
inverse) will be computed right into its intended destination, no
temporary matrix storage is ever allocated/used.
The README file talks at great length about the new concepts,
and gives many examples of their use. Let me briefly show just very
few:
// Checking up a difference of a matrix row and
// a given vector
assert( of_every(ConstMatrixRow(m,i)).max_abs(of_every(vr)) == 2.0 );
// multiplying the i-th matrix column by a vector,
// elementwise
to_every(MatrixColumn(m,i)) *= of_every(vc);
// Computing a max abs diag value of a matrix m
double max_diag = of_every(ConstMatrixDiag(m)).max_abs();
// Add 1 to the first matrix row:
to_every(MatrixRow(m,1)) += 1;
// An operation during Hooke-Jeeves minimization, hjmin.cc
void update_in_direction(Vector& from, Vector& to)
{
for(LAStreamOut tos = to, froms = from; !tos.eof(); )
{
register const double t = tos.peek();
tos.get() += (t - froms.peek());
froms.get() = t;
}
}
All matrices and vectors in these examples are accessed strictly
sequentially; no heap storage is ever allocated; all the operations on
"streams" are inlined.
Despite so many changes, the package interface remained almost
the same, as evidenced by very few changes I had to make to the old
validation code.
FFT and SVD parts aren't verified yet as I'm planning on some
significant extensions and optimizations.
LinAlg 4.1 also includes a new version of Brent univariate
optimizers, which used to be available separately. The code was
re-written to make it more in the genuine C++ style and spirit. The
most significant change is the way a function to minimize/investigate
- an argument to zeroin() and fminbr() - is specified. It is now a
functor, rather than an ordinary function pointer. This functor may be
a full-blown "clause", and even a, well, "lambda-expression". The
README file talks about this in detail.
The same source code base compiles on UNIX (gcc 2.7.2) and on
a Mac (PowerMac 8500/132, Codewarrior (12) Pro). It also compiles on
Linux, however, you have to take care of the div() "bug" in libc (or
compile without '-freg-struct-return' gcc flag).
As usual, the complete verification code and its output are
included (on HP-UX B.10.10 9000/770).
Downloads:
http://pobox.com/~oleg/ftp/packages/LinAlg41.tar.gz
The new README file, which includes detailed discussions of the new
concepts, examples, and a change log, can also be viewed separately:
http://pobox.com/~oleg/ftp/LinAlg.README.txt
http://pobox.com/~oleg/ftp/packages.html
The interface can be browsed separately
http://pobox.com/~oleg/ftp/packages/LinAlg.h
I will really appreciate any suggestions, comments, questions and
trouble reports.
Oleg
This archive was generated by hypermail 2b29 : Wed Feb 20 2002 - 03:20:06 EST