**Numeric Lua** is a numerical package for the Lua programming language. It includes support
for complex numbers, multidimensional matrices, random number generation, and
special functions. Most of the routines are simple wrappers for stable and
well known libraries from Netlib.

For instance, the special function package and the complex number package use routines from Netlib's SLATEC and DCDFLIB, the random number package uses Takuji Nishimura and Makoto Matsumoto's Mersenne Twister generator as the "engine" (uniform deviates) and Netlib's RANLIB for the remaining deviates, and the matrix package draws most of its numeric intensive routines from Netlib's ubiquitous BLAS and LAPACK packages.

Numeric Lua tries to maintain Lua's minimalist approach by providing bare-bone wrappers to the numerical routines. The user can use the outputs for further computations and is then fully responsible for the results (for a simple example, check usage below). Other Lua features are also available, such as OO simulation through metamethods and functional facilities. A basic API is provided in order to promote extensibility.

Numeric Lua is licensed under the same license as Lua -- the MIT license -- and so can be freely used for academic and commercial purposes.

Numeric Lua can be obtained from Luaforge.

The win32 binaries are provided for a quick and ready evaluation (and use); they contain a standalone version and the libraries. It should be noted that the binaries are generic, that is, they do not use any optimized version of BLAS/LAPACK.

If you got the source files, proceed to the next section for instructions.

**Note:** This is a very beta release, with
minimal testing. Any bugs, patches, comments and suggestions are welcome.

The installation process should be simple and easy on any flavor of Unix
environment (that includes MinGW and
Cygwin on Windows): in a nutshell, a plain
`make`

and a further `make install`

should do it. The
process can be divided in two: basic library compilation and package
compilation and installation.

After downloading the source `numlua-X.Y.tar.gz`

, where X.Y is
the current version, issue:

$ tar -zxf numlua-X.Y.tar.gz $ cd numlua-X.Y

and edit the `TARGET`

variable in Makefile to your platform
(currently only "linux" and "win32" are available as
options).

The first step is to compile the basic libraries in the `lib`

folder. These include `dcdflib`

, `fnlib`

(a Netlib's
SLATEC subset) and `ranlib`

(Mersenne Twister RNG + some routines
from Netlib's RANLIB). A light version of BLAS and LAPACK, containing only
double and complex double routines, is also available for completeness' sake;
note, however, that an optimized version of these libraries is preferred,
like ATLAS or Intel's MKL.

To compile the basic libraries, run:

$ cd lib; make

from `numlua-X.Y`

(root) dir, and wait a few minutes. :)

Some level of customization can be achieved by editing
`numlua-X.Y/lib/config`

. If you want to use an optimized version
of BLAS/LAPACK, you should also comment out `libblas.a`

and
`liblapack.a`

from `LIBS`

in
`numlua-X.Y/lib/Makefile`

, or issue `make noblas`

, to
save some time.

After compiling the basic libraries, change dir to `src`

:

$ cd ../src

Edit `Makefile`

to your needs/taste: set `TARGET`

as
in `numlua-X.Y/Makefile`

, set your BLAS library ("blas"
is set as default, in case you want to use the provided BLAS/LAPACK libraries
in `numlua-X.Y/lib`

), set installation paths, compilation flags and
programs. You are now ready to compile and install:

$ make && make install

A standalone version of Numeric Lua is also available in
`numlua-X.Y/etc`

. The drill is the same: cd
`numlua-X.Y/etc`

, edit `Makefile`

, ```
make &&
make install
```

.

There is no testing routine yet, unfortunately. However, you can try the following commands for a quick "test drive":

$ lua Lua 5.1.1 Copyright (C) 1994-2006 Lua.org, PUC-Rio > -- let's generate some random numbers (fast!) > require "numlua.rng" > r = rng(os.time()) > for i=1,5 do print(r:unif()) end -- numbers can vary, of course :) 0.47017181466799 0.63861816644203 0.55608703091275 0.68107319541741 0.24100433581043 > > -- now, let's add some math to spice up > require "numlua.spfun" > print(spfun.log1p(1e-16), math.log(1 + 1e-16)) -- keep precision 1e-16 0 > print(spfun.isinf(1/0), spfun.isnan(0/0)) true true > -- normal quantiles, mean = 0, var = 1 > print(spfun.qnorm(0.025, 0, 1), spfun.pnorm(spfun.qnorm(0.025))) -1.9599639845401 0.025 > > -- complex numbers > require "numlua.complex" > print(complex(1,2) + complex.i) 1 + 3i > print(complex.exp(complex.sqrt(-2))) -- cos(sqrt(2)) + sin(sqrt(2))*i 0.15594369476537 + 0.98776594599274i > > -- matrices > require "numlua.matrix" > x = matrix(2,3,4) -- 2 x 3 x 4 matrix > for k, i, j in x:entries() do x[k][i][j] = k * (i + j) end > matrix.foreach(x, print) -- list x ... (contents are listed) > -- let's try the type system > x = matrix.zeros(1000) -- 1000 x 1000 matrix > x:apply(function(i, j) return i + j end) -- x[i][j] = i+j > print(x:gettype()) -- x is general G > -- let's profile transpose(x) * x > local t = os.clock(); print((#x) * x); print(os.clock() - t) matrix: 0xb586001c 1.32 > x:settype() -- guess x's type > print(x:gettype()) -- x is symmetric S > local t = os.clock(); print((#x) * x); print(os.clock() - t) 1.28 > -- a bit faster... let's compute the inner product directly > local t = os.clock(); print(x:inner()); print(os.clock() - t) 0.66

In this first example it is shown that a matrix can be used in a "standard" way, similar to a regular Lua table:

-- returns a symm posdef matrix with integer entries from Pascal's triangle function pascal(n) assert(type(n) == "number", "number expected, got " .. type(n)) assert(n > 0, "order must be positive") local m = matrix(n,n) -- first row for i = 1, n do m[1][i] = 1 end -- remaining rows for i = 2, n do m[i][1] = 1 for j = 2, n do m[i][j] = m[i][j-1] + m[i-1][j] end end return m end

Now, to demonstrate some functional facilities:

-- prints the entries of a matrix, right aligned function pretty(m) assert(matrix.check(m), "matrix expected, got " .. type(m)) assert(matrix.getn(m) <= 2) -- maximum entry length local ml = matrix.reduce(m, function(x, y) local s = string.len(y); if x < s then x = s end return x end, 0) if matrix.getn(m) == 1 then -- vector? matrix.foreache(m, function(e) -- for each entry print(string.rep(" ", 2 + ml - string.len(e)) .. e) end) else -- array for i = 1, m:size(1) do -- print each row print(matrix.reduce(m[i], function(x, y) return x .. string.rep(" ", 2 + ml - string.len(y)) .. y end, "")) end end end

Finally, some more specialized routines:

-- create random number generator local r = rng(os.time()) -- generate n multivariate normal deviates X_i ~ N(mean, var) function mvnorm(mean, var, n) -- do some boring testing on the args ... -- compute square root of var using Cholesky decomposition local u, info = matrix.chol(var) -- #u * u = var assert(info == 0, "matrix is not positive definite") -- generate y where y:col(i) ~ N(0_m, I_m) local m = mean:size() local y = matrix(m, n):foreache(function(e) return r:norm(0, 1) end) -- compute #u * y + mean ~ N(mean, var) y = #u * y for i, j, e in y:entries() do y[i][j] = e + mean[i] end return y end

-- solve system of linear equations A*X = B function lsolve(A, B) assert(matrix.check(A), "matrix expected, got " .. type(A)) assert(matrix.check(B), "matrix expected, got " .. type(B)) -- try to solve by elimination first local x, info = matrix.linsolve(A, B) if info < matrix.EPS then -- A is nearly singular? print("WARNING: matrix is singular to machine precision: rcond = " .. info) -- try by least squares then x, info = matrix.lss(A, B) end return x end

These routines can be found at
`lua/matrix.lua`

.

Unfortunately, there is no documentation available (yet!).