Getting started - wlevine/nmatrix GitHub Wiki
Observation
Both this tutorial and How to create an NMatrix should be merged, probably into Tentative NMatrix Tutorial. We need to centralize NMatrix documentation, so it gets easier for other people to know where to get help.
NMatrix is a linear algebra library for Ruby, written mostly in C. It can be used with or without SciRuby, but is part of the SciRuby project. NMatrix was inspired by and based heavily upon NArray, by Masahiro Tanaka.
Installation
See Installation for details.
Generally, you'd only need gem install nmatrix. If you want to install from source:
git clone git://github.com/SciRuby/nmatrix.git
cd nmatrix
rake install
Features
In addition to dense matrices, NMatrix supports two types of sparse storage (list-of-list and Yale). We support the following data classes:
- Byte (unsigned 8-bit integer)
- Integer (8, 16, 32, and 64-bit)
- Floating point (32 and 64-bit)
- Complex (64 and 128-bit)
- Rational (32, 64, and 128-bit)
- Ruby objects
NMatrix makes use of and exposes the CBLAS API included in ATLAS. However, we have also written BLAS templates for NMatrix-supported data types that are not included in ATLAS or LAPACK (e.g., rationals).
Differences from other matrix libraries
In contrast to ruby's Matrix class, and scipy.matrix, the * operator is used for element-wise multiplication and the #dot method for matrix-matrix or matrix-vector multiplication.
Brief Tutorial
Creating matrices
Currently, you create NMatrix objects with the NMatrix#new method. The first parameter is the dimensions of the matrix, and the second parameter contains the initial values. For example, the following creates a 3x4 matrix filled with zeros:
>> n = NMatrix.new( [3, 4], 0)
0 0 0 0
0 0 0 0
0 0 0 0
=> nil
If the first parameter is an integer, NMatrix.new returns a square matrix.
>> n = NMatrix.new( 4, 0) # Same as NMatrix.new( [4,4], 0 )
0 0 0 0
0 0 0 0
0 0 0 0
0 0 0 0
=> nil
Matrices need not have only two dimensions:
# Create a 3D matrix (4 x 4 x 4)
b = NMatrix.new([4,4,4], 0)
In addition, the second parameter to NMatrix.new can be an array which specifies all the values in the matrix:
>> B = NMatrix.new([2,5], [1,2,3,4,5,6,7,8,9,0])
1 2 3 4 5
6 7 8 9 0
If the array of values has fewer elements than the matrix requires, NMatrix.new will repeat the array as many times as required:
>> q = NMatrix.new([3,6], [1,2,3])
1 2 3 1 2 3
1 2 3 1 2 3
1 2 3 1 2 3
Simpler constructor for matrices
Thanks to Daniel Carrera, we have an in-line constructor that makes creating matrices a very simple and readable operation:
a = N[ [1,2,3,4] ] # => [ [1.0 2.0 3.0 4.0] ]
a = N[ [1,2,3,4], dtype: :int32 ] # => [ [1 2 3 4] ]
a = N[ [1,2,3], [3,4,5] ] # => 1.0 2.0 3.0
# 3.0 4.0 5.0
Note that N[1,2,3,4] will create a "dimensionless" (technically 1-dimensional) vector. We advise creating a vector with the same dimensionality as the matrices you're working with. For two dimensions:
N[ [1,2,3,4] ] # row vector
N[ [1],[2],[3],[4] ] # column vector
N[ [1,2,3,4] ].transpose # column vector
One-dimensional matrices
One-dimensional matrices aren't really one-dimensional. They are n-dimensional, but n-1 of their dimensions are length 1.
Here's how you create vectors with different orientations:
>> x = NMatrix.new([4], [1,2,3,4]) # no orientation (dimensionless/1-dimensional)
>> x = NMatrix.new([4,1], [1,2,3,4]) # column vector in two dimensions (4 rows, 1 column)
>> x = NMatrix.new([1,4], [1,2,3,4]) # row vector in two dimensions (1 row, 4 columns)
>> x = NMatrix.new([4,1,1], [1,2,3,4]) # column vector in three dimensions (4 rows, 1 column, 1 layer)
Matrix functions and operations
NMatrix supports the standard arithmetic operations on matrices:
# a = NMatrix.new([3,3], 2) # 3x3 matrix containing all 2's
a + a # Element-wise / Matrix addition.
a - a # Element-wise / Matrix subtraction.
a ** 3 # Element-wise exponentiation.
a / 2 # Element-wise division.
a * a # Element-wise multiplication.
a.dot a # Matrix dot-product multiplication.
Shortcuts
There are several shortcuts for creating NMatrix or NVector objects. They are as follows:
NMatrix.ones [3, 3] # A 3x3 matrix of ones.
NVector.ones 3 # A 3-vector of ones.
NMatrix.zeros [4, 4] (or zeroes [4, 4]) # Creates a matrix of zeros with dimensions as parameters.
NVector.zeros 4 # Creates a vector of zeros
NMatrix.identity [5, 5] (or eye [5, 5]) # A 5x5 identity matrix. Only works with NMatrix.
NMatrix.seq [3, 4] # A 3x4 matrix with values from 0 to 11.
NVector.seq 8 # An 8-vector with a sequence of values from 0 to 7.
random [6,6] # A 6x6 matrix of random values between 0 and 1.
m.row(0) # Return the 1st column of the matrix m. Only works with NMatrix.
m.column(2) # Return the 3rd column of the matrix m. Only works with NMatrix.
NVector.indgen(8) # To make IDL users happy. Same as `NVector.seq(8, :int32)`.
NVector.findgen(8) # Same as `NVector.seq(8, :float32)`.
NVector.bindgen(8) # Same as `NVector.seq(8, :byte)`.
NVector.cindgen(8) # Same as `NVector.seq(8, :complex64)`.
NVector.linspace(0, pi, 100) # A vector of 100 values from 0 to pi, inclusive. Only works with NVector.
Solving linear systems
For now, the best way is probably to use NMatrix::LAPACK::clapack_gesv, which is the Ruby-exposed API for ATLAS's CLAPACK implementation of SGESV/DGESV/CGESV/ZGESV.
FIXME: Needs documentation (and a spec)
Matrix Decomposition
Here is LU decomposition (using CLAPACK):
n = NMatrix.new(3, [4,9,2,3,5,7,8,1,6], dtype: :float64)
lu = n.factorize_lu # the result are the Lower and Upper components stored in a single matrix
There are multiple methods for singular value decomposition, all of which use LAPACK instead of CLAPACK -- so they should work even if you haven't installed your own ATLAS (e.g., on Mac OS X):
u, sigma, v_transpose = n.gesvd # LAPACK's GESVD function
u, sigma, v_transpose = n.gesdd # Divide-and-conquer strategy
You can also do eigenvalue decomposition using NMatrix::LAPACK::lapack_geev. See the specs for an example.
Some of these functions still need NMatrix instance methods. Please consider contributing one!
Advanced features
NMatrix supports several data types. NMatrix will try to guess the dtype based on the first initial value you provide (e.g., NMatrix.new(4, [0.0, 1]) will be :float32), but you can choose to provide a dtype in addition to or in lieu of initial values:
>> n = NMatrix.new(4, [0,1], dtype: :rational128)
0/1 1/1 0/1 1/1
0/1 1/1 0/1 1/1
0/1 1/1 0/1 1/1
0/1 1/1 0/1 1/1
>> m = NMatrix.new(4, dtype: :int64) # Caution! No initialization of contents. Example output:
8, 140642513449760, 0, 0
0, 0, 0, 0
0, 0, 0, 140642505889952
140642512627712, -4294967290, 1, 140734747434816
>> m = NMatrix.new(4, dtype: :rational128) # not recommended (but allowed)
>> m = NMatrix.new(4) # 4x4 Ruby :object matrix containing only nil
The storage type (stype) can also be specified, prior to the dimension argument. However, with sparse storage formats, initial values don’t make sense, and these matrices will contain zeros by default.
# empty list-of-lists-of-lists 4x3x4 matrix
n = NMatrix.new([4,3,4], dtype: :int64, stype: :list)
# Ruby objects in a 'Yale' sparse matrix
m = NMatrix.new([5,4], stype: :yale)
# A byte matrix containing a gradient
o = NMatrix.new(5, [0,1,2,3,4], dtype: :byte)
The matrix m created above is a Yale-format sparse matrix, or more specifically, “new Yale,” which differs from “old Yale” in that the diagonal is stored separately from the non-diagonal elements. Thus, diagonals can be accessed and set in constant time.
Currently, all storage is row-based.
You can also convert between any of these three stypes using cast:
n = NMatrix.new(4, dtype: :int64, stype: :list)
n[0,0] = 5
n[0,3] = -2
dense = n.cast(:dense) # you may also supply other information to cast, e.g.,
dense_double = n.cast(stype: :dense, dtype: :float64)
Let's say you have a dense matrix consisting of only 1s with a few 0s. There are two ways to cast this to sparse matrices. You can have the sparse storage store only the non-zero values, or only the non-one values, by specifying a default value for cast.
list = dense.cast(stype: :list, default: 1) # 1-based matrix
list = dense.cast(:list) # default is 0-based
The same works for Yale matrices.