02 Matrix

Introduction

Package matrix could be loaded via the standalone binary, or in Lua with require("aprilann.matrix").

A matrix is a multidimensional data container, by default float data. It is similar to the concept of tensor, as defined in libraries like Torch. This notion of tensor is not to be confused with tensors in physics and engineering, known as tensor fields.

The data could be stored following row_major or col_major orders, but from the outside there is no difference for the user. However, CUBLAS implementation of fast mathematical operations are only allowed in col_major order. Additionally a matrix could be transposed or not, being the transposition a symbolical state, sharing the same data reference as the non-transposed matrix.

Matrices in col_major are used as input and output of ANN components. WARNING!!! because of the col_major nature of ANNs, the data generated by row_major matrices is interpreted as transposed, which is very important when your data matrix has more than one dimension. This issue will be discussed properly on ANN related wiki pages.

From Lua, a matrix is declared using one of the two available constructors, one for row_major and the other for col_major:

> -- row major constructor
> m1 = matrix(3,3) -- this is a 3x3 matrix of floats
> -- It is also possible to receive a table with data (in row-major order)
> m2 = matrix(3,3, {1, 2, 3, 4, 5, 6, 7, 8, 9})
> -- it is also possible to use the following equivalent constructor
> m1 = matrix.row_major(3,3, {1, 2, 3, 4, 5, 6, 7, 8, 9})
> -- and this is the col-major constructor
> m3 = matrix.col_major(3,3, {1, 2, 3, 4, 5, 6, 7, 8, 9})
> print(m2)
1 2 3
4 5 6
7 8 9
# Matrix of size [3,3] in row_major [0x23d5b20 data= 0x23d5e90]
> print(m3)
1 2 3
4 5 6
7 8 9
# Matrix of size [3,3] in col_major [0x23d62c0 data= 0x23d6630]

Observe that print function shows the same for m2 and m3, but internally the data is in different order. The pretty print of a matrix shows the data and a commented line with the size of the matrix and two memory pointers values, the first is the pointer to the C++ object related with the given matrix, and the second is a pointer to the C++ data object where values are stored.

The matrix and its data is separated to allow the declaration of sub-matrices:

> m4 = m2:slice({2,1},{1,3})
> print(m4)
4 5 6
# Matrix of size [1,3] in row_major [0x2218a90 data= 0x23d5e90]

In this case, the matrix m4 is a slice which begins at position {2,1} of matrix m2, and has sizes {1,3} in each dimension. Note that the matrix pointer is different, but the data pointer is the same as for m2 (any change to m4 will be reflected at m2.)

Besides, it is possible to do a sub-matrix cloning the data (deep copy) if you add to slice method a new boolean argument with true value:

> m5 = m2:slice({2,1},{1,3}, true)
> print(m5)
4 5 6
# Matrix of size [1,3] in row_major [0x220c6f0 data= 0x2250d60]

A shortcut for slicing is to use the operator m(...), which has a meaning similar to Matlab or Octave. It is overloaded and it is parsed into a slice method call, so, it is slower than using directly the slice method, but it is easier to understand. This operator receives a variable number of arguments, as many as dimensions has the caller matrix object. In every dimension position, three possible values could be given:

• A number, indicating an exact dimension position: m = m2(2,3)

• A table, with start and end positions in the dimension: m = m2(2, {1,3})

• A string, with start and end positions in the dimension: m = m2(2, '1:3')

The following examples are equivalent:

> = m2:slice({2,1},{1,3})
 4           5           6         
# Matrix of size [1,3] in row_major [0xf44470 data= 0xdeae90]
> = m2(2, {1,3})
 4           5           6         
# Matrix of size [1,3] in row_major [0xf3c8d0 data= 0xdeae90]
> = m2(2, '1:3')
 4           5           6         
# Matrix of size [1,3] in row_major [0xe32dc0 data= 0xdeae90]

It is possible to use the string shortcut to indicate a whole dimension, or only start/end positions:

> = m2(':', '2:3')
 2           3         
 5           6         
 8           9         
# Matrix of size [3,2] in row_major [0xef1260 data= 0xdeae90]
> = m2('2:', '2:3')
 5           6         
 8           9         
# Matrix of size [2,2] in row_major [0xe08f50 data= 0xdeae90]
> = m2(':2', '2:3')
 2           3         
 5           6         
# Matrix of size [2,2] in row_major [0xf04290 data= 0xdeae90]

Hints

NOTICE: Almost all matrix methods returns the caller matrix (when it is possible), allowing to chain transformation sequences.

NOTICE: In Lua the arrays start at 1 instead of 0, so, in the matrix methods the dimensions start at 1.

MMapped matrix

It is possible to force the declaration of matrix memory as a mmapped anonymous file. First you need the package require("aprilann.mathcore").

> mathcore.set_mmap_allocation(true)
> -- the following matrix will be allocated as mmapped memory
> m = matrix(2,2):linear()
> print(m)
> mathcore.set_mmap_allocation(false)

Another way is to serialize a matrix in MMap format (see serialization section).

Basic matrix methods

[table | number] = m:dim([number])

It returns the size of matrix dimensions. Without arguments, it returns a Lua table with the sizes. If an argument is given, it returns the size of the given dimension.

> a = matrix(4,3,2)
> = a:dim()
table: 0x23a4780
> = table.concat(a:dim(), " ")
4 3 2
> = a:dim(1), a:dim(2), a:dim(3)
4	3	2

number = m:get(p1, p2, ...)

This method returns the value of a given matrix position.

> a = matrix(3,4,{1,2,3,4, 5,6,7,8, 10,11,12,13})
> = a:get(1,1)
1
> = a:get(2,3)
7

matrix = m:set(p1, p2, ..., value)

This method sets the value of a matrix position, and returns the caller matrix, allowing a sequence of sets.

> a = matrix(3,4,{1,2,3,4, 5,6,7,8, 10,11,12,13})
> a:set(2,3, 10000)
> a:set(2,4, 500):set(4,1, 200)
> = a
1 2 3 4
5 6 10000 500
200 11 12 13
# Matrix of size [3,4] in row_major [0x27093d0 data= 0x2709960]

matrix = m:clone([major_order])

It allows to clone matrices (deep copy). Besides, this method allows to change major order of data (row_major or col_major). If the caller matrix was transposed, the resulting clone will contain the data in a transposed shape, but transposed flag will be false.

> a = matrix(2,3,{1,2,3, 4,5,6}) -- row major matrix
> b = a:clone()                  -- clone (or deep copy) of a
> c = a:clone("col_major")       -- clone of a in col major order

matrix = m:copy_from_table(table)

This method copies the data in the given table into the caller matrix, traversing the matrix in row_major order, as in matrix constructor. The table must fit in matrix size. The caller matrix is returned.

> a = matrix(2,3)
> a:copy_from_table({1,2,3, 4,5,6})

matrix = m:map(m1, m2, ..., function)

Maps the matrix values by a given list of matrices and a Lua map function. The Lua function will be called for every possible matrix position. The Lua function receives the caller matrix value at the given position, the value of the second matrix, the value of the third matrix, and so on. The Lua function returns nil, or only one value which will be assigned to the caller matrix in-place. All the matrices must have the same dimension sizes. The number of given matrices could be >= 0.

> m  = matrix(2,2):linear()
> m2 = matrix(2,2):linear(10)
> m3 = matrix(2,2):linear(100)
> = m
 0           1         
 2           3         
# Matrix of size [2,2] in row_major [0x1f12050 data= 0x1f0f6a0]
> = m2
 10          11        
 12          13        
# Matrix of size [2,2] in row_major [0x1f11cc0 data= 0x1f12110]
> = m3
 100         101       
 102         103       
# Matrix of size [2,2] in row_major [0x1f12740 data= 0x1f11e00]
> m:map(m2,m3,function(x,y,z) return x+y+z end)
> = m
 110         113       
 116         119       
# Matrix of size [2,2] in row_major [0x1f12050 data= 0x1f0f6a0]

matrix = m:rewrap(size1, size2, ...)

This method only works if the data is contiguous in memory. The caller matrix is reinterpreted as if it was of another number of dimensions and sizes. A different matrix instance is returned, but the data pointer is shared.

> a = matrix(2,3,{1,2,3, 4,5,6})
> = a
1 2 3
4 5 6
# Matrix of size [2,3] in row_major [0x2700850 data= 0x2700900]
> b = a:rewrap(3,2)
> = b
1 2
3 4
5 6
# Matrix of size [3,2] in row_major [0x2701360 data= 0x2700900]

matrix = m:select(dimension, index [, matrix])

This methods returns a matrix with one less dimension, resulting of select at the caller matrix the indicated dimension at the given index. The resulting matrix references the internal data of original matrix. If given, the third argument must be a matrix which will be used to store the result of the select call, and must fit the expected dimensions. In this last case, the computation effort is dismissed to constant.

> m = matrix.col_major(4,3):zeros()
> = m
0 0 0
0 0 0
0 0 0
0 0 0
# Matrix of size [4,3] in col_major [0x23dcab0 data= 0x23727e0]
> = m:select(2,2):fill(9)
9 9 9 9
# Matrix of size [4] in col_major [0x23dd330 data= 0x23727e0]
> = m:select(1,3):fill(4)
4 4 4
# Matrix of size [3] in col_major [0x23dd790 data= 0x23727e0]
> = m
0 9 0
0 9 0
4 4 4
0 9 0
# Matrix of size [4,3] in col_major [0x23dcab0 data= 0x23727e0]

NOTE that the third argument matrix must be created by a previous call to select over the same dimension (but not the same index). As example, the following design pattern gives the same variable as third argument result and as left-side of the expression, allocating the memory in the first loop iteration, and reusing it in the following:

> m = matrix(4,5):linear()
> for i=1,m:dim(2) do
    result = m:select(2,i,result)
  end

matrix = m:transpose()

This method returns a matrix which is a transposition of the caller object. Note that both, the caller and the transposition, reference the same data.

> m = matrix(3,4):linear()
> = m
 0           1           2           3         
 4           5           6           7         
 8           9           10          11        
# Matrix of size [3,4] in row_major [0x2777140 data= 0x27799b0]
> = m:transpose()
 0           4           8         
 1           5           9         
 2           6           10        
 3           7           11        
# Matrix of size [4,3] in row_major [0x274e620 data= 0x27799b0]

matrix = m:slice(position, size [, clone])

This methods produces a sub-matrix of the caller matrix. By default, the returned sub-matrix shares the data pointer with the caller, but it is also possible to do a deep copy sub-matrix. The syntax is:

m:slice(pos_table, size_table, clone=false)

being pos_table a Lua table with the position of first element (starting at 1, not 0), and size_table a Lua table with the size of each dimension. The last argument, clone, is an optional boolean (by default false) indicating if the resulting matrix will be a clone or not.

> a = matrix(3,4,{1,2,3,4, 5,6,7,8, 10,11,12,13}) -- row major matrix
> = a
1 2 3 4
5 6 7 8
10 11 12 13
# Matrix of size [3,4] in row_major [0x2706530 data= 0x2706b00]
> b = a:slice({2,1},{2,2}) -- slice at position (2,1) with size 2x2
> = b
5 6
10 11
# Matrix of size [2,2] in row_major [0x2707cd0 data= 0x2706b00]
> -- same slice as before but making a clone (deep copy)
> b = a:slice({2,1},{2,2}, true)
> = b
5 6
10 11
# Matrix of size [2,2] in row_major [0x2708a20 data= 0x2708ad0]

matrix = m(s1, s2, ...)

A shortcut for slicing, using the __call metamethod. It is similar to Matlab or Octave slice operators. The call is parsed and converted into a slice method call, so, it is slower than using directly the slice method, but it is easier to understand. This operator receives a variable number of arguments, as many as dimensions has the caller matrix object. In every dimension position, three possible values could be given:

• A number, indicating an exact dimension position: m = m2(2,3)

• A table, with start and end positions in the dimension: m = m2(2, {1,3})

• A string, with start and end positions in the dimension: m = m2(2, '1:3')

The following examples are equivalent:

> = m2:slice({2,1},{1,3})
 4           5           6         
# Matrix of size [1,3] in row_major [0xf44470 data= 0xdeae90]
> = m2(2, {1,3})
 4           5           6         
# Matrix of size [1,3] in row_major [0xf3c8d0 data= 0xdeae90]
> = m2(2, '1:3')
 4           5           6         
# Matrix of size [1,3] in row_major [0xe32dc0 data= 0xdeae90]

It is possible to use the string shortcut to indicate a whole dimension, or only start/end positions:

> = m2(':', '2:3')
 2           3         
 5           6         
 8           9         
# Matrix of size [3,2] in row_major [0xef1260 data= 0xdeae90]
> = m2('2:', '2:3')
 5           6         
 8           9         
# Matrix of size [2,2] in row_major [0xe08f50 data= 0xdeae90]
> = m2(':2', '2:3')
 2           3         
 5           6         
# Matrix of size [2,2] in row_major [0xf04290 data= 0xdeae90]

matrix = matrix.join(dimension, m1, m2, ...)

This function joins the given matrices by the given dimension. All the dimensions of the matrices must be the same, except the given dimension, which could differ. Warning, this method duplicates the memory needed, because all the matrices are copied to the destination matrix.

> m1 = matrix(10,2):linear()
> m2 = matrix(10,3):linear()
> outm = matrix.join(2, m1, m2)
> = outm
 0           1           0           1           2         
 2           3           3           4           5         
 4           5           6           7           8         
 6           7           9           10          11        
 8           9           12          13          14        
 10          11          15          16          17        
 12          13          18          19          20        
 14          15          21          22          23        
 16          17          24          25          26        
 18          19          27          28          29        
# Matrix of size [10,5] in row_major [0x1f9c100 data= 0x1f9c1c0]

matrix = m:clamp(lower, upper)

This method clamps the matrix components to a given range [min,max], modifying the matrix in-place. The caller matrix instance is returned.

> a = matrix(3,3,{1,2,3,4,5,6,7,8,9})
> = a
1 2 3
4 5 6
7 8 9
# Matrix of size [3,3] in row_major [0xe56a30 data= 0xe56f40]
> a:clamp(3,6)
> = a
3 3 3
4 5 6
6 6 6
# Matrix of size [3,3] in row_major [0xe56a30 data= 0xe56f40]

matrix = m:adjust_range(min, max)

This method modifies in-place the matrix components, interpolating the values to be in the given range [min,max]. The caller matrix is returned.

> a = matrix(3,3,{1,2,3,4,5,6,7,8,9})
> a:adjust_range(3,6)
> = a
3 3.375 3.75
4.125 4.5 4.875
5.25 5.625 6
# Matrix of size [3,3] in row_major [0x25cca30 data= 0x25ccf40]
> = a:adjust_range(0,1)
0 0.125 0.25
0.375 0.5 0.625
0.75 0.875 1
# Matrix of size [3,3] in row_major [0x25cca30 data= 0x25ccf40]
> = a:adjust_range(1,9)
1 2 3
4 5 6
7 8 9
# Matrix of size [3,3] in row_major [0x25cca30 data= 0x25ccf40]

string = m:get_major_order()

Returns a string with the major order: row_major or col_major

boolean = m:is_transposed()

Returns a boolean which indicates if the matrix is a transposition of another or not.

boolean = m:is_contiguous()

Indicates if the matrix internal data is contiguous at memory.

matrix = m:contiguous()

Returns a contiguous version of the caller matrix. If the matrix is contiguous, returns itself. Otherwise, returns a copy of the caller. Note that, if the matrix is a slice or a transposition of another matrix, therefore it could be non-contiguous.

Data initialization methods

matrix = m:fill(number)

This is an in-place method which sets all components to a given value.

> a = matrix(2,3):fill(4) -- a 2x3 matrix filled with 4
> = a
4 4 4
4 4 4
# Matrix of size [2,3] in row_major [0x26ff9b0 data= 0x26ffa20]

matrix = m:zeros()

This is equivalent to m:fill(0)

matrix = m:ones()

This is equivalent to m:fill(1)

matrix = m:linear(start=0, step=1)

Initializes the matrix starting at the given index and using the given step. The index and the step is optional.

> m = matrix(3,2,2):linear(1,2)
> = m

# pos [1,1,1]
1 3
5 7

# pos [2,1,1]
9 11
13 15

# pos [3,1,1]
17 19
21 23
# Matrix of size [3,2,2] in row_major [0x149de00 data= 0x149daa0]
> m = matrix(2,2):linear()
> = m
0 1
2 3
# Matrix of size [2,2] in row_major [0x149f110 data= 0x149f1e0]

matrix = m:linspace(a=1, b=m:size())

Initializes the matrix with a linear space distribution. It receives two optional arguments, why default a=1 and b=m:size(). It returns the caller matrix.

> m = matrix(5):linspace(1,20)
> = m
 1           5.75        10.5        15.25       20        
# Matrix of size [5] in row_major [0x291f200 data= 0x291ecd0]

matrix = m:logspace(a=1, b=m:size(), base=10)

Initializes the matrix with a logarithmic distribution between a and b with the given logarithm base. It receives three optional arguments, why default a=1, b=m:size() and base=10. It returns the caller matrix.

> m = matrix(5):logspace(0.001,0.1)
> = m
 0.001       0.00316228  0.01        0.0316228   0.1       
# Matrix of size [5] in row_major [0x291fed0 data= 0x291fd50]

matrix = m:uniform(lower, upper [, random] )

This method initializes the matrix with random integers taken uniformly from the given range of values:

> m = matrix(10):uniform(0,10,random(1234))
> = m
3 6 5 4 8 9 1 7 9 10
# Matrix of size [10] in row_major [0x2716b10 data= 0x2716490]

The random object is optional, but to ensure reproducibility it is recommended.

matrix = m:uniformf(lower=0, upper=1 [, random] )

This method initializes the matrix with random floats taken uniformly from the given range of values:

> m = matrix(2,2):uniformf(-10, 10, random(1234))
> = m
-6.16961 -0.0467267
2.44218 6.35677
# Matrix of size [2,2] in row_major [0x1000e90 data= 0xe47410]

The random object is optional, but to ensure reproducibility it is recommended.

matrix = m:diag(number)

This method sets the matrix diagonal components to a given value, modifying in-place the caller matrix. For any number of dimensions, the diagonal are whose components which positions are equals at all dimensions.

> a = matrix(3,3,3):ones():diag(5)
> = a

# pos [1,1,1]
5 1 1
1 1 1
1 1 1

# pos [2,1,1]
1 1 1
1 5 1
1 1 1

# pos [3,1,1]
1 1 1
1 1 1
1 1 5
# Matrix of size [3,3,3] in row_major [0x1718f10 data= 0x1718d50]

Matrix serialization

string = m:toString( mode='ascii' )

This method returns a Lua string which represents the caller matrix. It receives an optional argument indicating if the matrix data will be stored in ascii or binary format (by default ascii).

> a = matrix(3,5):ones()
> = a
1 1 1 1 1
1 1 1 1 1
1 1 1 1 1
# Matrix of size [3,5] in row_major [0xd80a10 data= 0xd815d0]
> = a:toString()
3 5
ascii row_major
1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 

> = a:toString("ascii")
3 5
ascii row_major
1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 

> = a:toString("binary")
3 5
binary row_major
8Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe98Ffe9

matrix = matrix.fromString(filename)

This method loads a matrix from a Lua string generated by method matrix.toString.

> a = matrix.fromString[[3 5
>> ascii row_major
>> 1 1 1 1 1 1 1 1 1
>> 1 1 1 1 1 1 
>> ]]
> = a
1 1 1 1 1
1 1 1 1 1
1 1 1 1 1
# Matrix of size [3,5] in row_major [0xd983b0 data= 0xdfe5c0]

string = m:to_lua_string(mode='binary')

This method is similar to toString(), but instead of returning a string with matrix data, it also includes Lua code which makes the string loadable. It is useful to serialize the matrix into a file.

> m1  = matrix(2,2):uniformf()
> str = "return " .. m1:to_lua_string()
> m2  = loadstring(str)()
> print(m1 == m2)
true

m:toFilename(filename, mode='ascii')

This method stores a matrix in a given filename. It also receives an optional argument with ascii or binary (by default ascii). It allows to compress the output file using GZIP, if the filename has '.gz' extension.

> a = matrix(3,3)
> a:toFilename("a.mat", "binary")
> a:toFilename("a.mat.gz", "binary")

matrix = matrix.fromFilename(filename [,order] )

This method loads a matrix from a given filename, expecting the format used by matrix.toFilename method. It allows to load compressed files using GZIP, if the filename has '.gz' extension. The second argument is optional, and if given it forces to load the matrix following the given order string, which could be row_major, col_major, or nil. By default it is nil, loading the matrix with the order given by the file content.

> a = matrix.fromFilename("a.mat")
> a = matrix.fromFilename("a.mat.gz")

m:toTabFilename(filename)

This method stores a matrix in a given filename, but without header, and formatting the data to be formatted by lines and spaces, one matrix row per line. It is limited to bi-dimensional matrices. It allows to compress the output file using GZIP, if the filename has '.gz' extension.

> a = matrix(3,3)
> a:toTabFilename("a.mat")
> a:toTabFilename("a.mat.gz")

matrix = matrix.fromTabFilename(filename [,"row_major"] )

This method loads a matrix from a given filename, formatted as done by matrix.toTabFilename. The size of the matrix is computed in a first loop over all the data, so this method needs two passes to load the matrix. It allows to load compressed files using GZIP, if the filename has '.gz' extension. The second argument is optional, and if present it could be row_major o or col_major, indicating in which format you want to load the matrix.

> a = matrix.fromTabFilename("a.mat")
> a = matrix.fromTabFilename("a.mat.gz")

m:toTabStream(file)

This method appends a matrix in a given file object, following the same format mas m:toTabFilename(...). This method allows to append several matrices in the same file object (stream), which could be opened using io.open, and could be a standard text file or a gzipped file. The destination file must be opened for writing, otherwise the method will fail with error:

stdin:1: calling 'toTabStream' on bad self (FILE* expected, got nil)

The following is an example of how this method works:

> a = matrix(3,3)
> b = matrix(4,3)
> dest = io.open("dest.mat.gz", "w")
> a:toTabStream(dest)
> b:toTabStream(dest)
> dest:close()
> aux = matrix.fromTabFilename("dest.mat.gz")
> = aux
 1           2           3         
 4           5           6         
 7           8           9         
 1           2           3         
 4           5           6         
 7           8           9         
 10          11          12        
# Matrix of size [7,3] in row_major [0x29e0240 data= 0x299e2a0]

m:toMMap(filename)

Stores the matrix in a file in a binary machine-dependent format, so it could be loaded using mmap function (matrix.fromMMap). The endianism must be the same between machines where matrix is stored/loaded.

matrix = matrix.fromMMap(filename [,write [,shared]])

Loads the matrix from a file in a binary machine-dependent format, by using the mmap function (matrix.toMMap). The endianism must be the same between machines where matrix is stored/loaded. Two additional boolean arguments are allowed. The second boolean argument indicates if writing is available, by default it is true. Be careful, if writing is set to false, any attempt of writing will throw a segmentation fault. The third boolean argument indicates if the data is shared between different processes, by default it is true. If both arguments are true, any writing will be available to any process which shares this map. Besides, writings will be synchronized in the hard disk (but not instantly). If writing is true, but shared is false, then the memory is mapped as copy-on-write. For more info, see the manual page of mmap function (PROTECT_WRITE, MAP_SHARED and MAP_PRIVATE).

table = m:toTable()

This method returns a plain Lua table (one-dimensional table) which contains the matrix data in row_order, as expected by matrix constructors.

> a = matrix(3,2,{1,2,3,4,5,6})
> = a
1 2
3 4
5 6
# Matrix of size [3,2] in row_major [0x9ddce0 data= 0x9ddd30]
> t = a:toTable()
> = table.concat(t, " ")
1 2 3 4 5 6

Low-level matrix access

These methods allows raw accessing of matrix components.

number = m:size()

This method returns the number of elements in the matrix.

> a = matrix(3,4,{1,2,3,4, 5,6,7,8, 10,11,12,13})
> = a:size()
12

table = m:stride()

This method is similar to m:dim, but returning the stride of the dimension (the offset between elements at each dimension)

> a = matrix(4,3,2)
> = a:stride()
table: 0x23a5fe0
> = table.concat(a:stride(), " ")
6 2 1
> = a:stride(1), a:stride(2), a:stride(3)
6	2	1
> a = matrix.col_major(4,3,2)
> = a:stride(1), a:stride(2), a:stride(3)
1	4	12

number = m:offset()

It returns the offset from data first position. Only sub-matrices has an offset!=0.

> a = matrix(2,3)
> = a:offset()
0
> b = a:slice({2,1},{1,1})
> = b:offset()
3

number = m:raw_get(pos)

It receives a raw position at the underlying data pointer, and returns its value. It is useful to combine stride and offset methods in order to compute the raw position.

> a = matrix(3,2, {1,2,3,4,5,6})
> = a
1 2
3 4
5 6
# Matrix of size [3,2] in row_major [0x144fce0 data= 0x144fd90]
> = a:raw_get(a:offset() + a:stride(1)*1 + a:stride(2)*0), a:get(2,1)
3	3

NOTE! that the strides are multiplied by matrix position minus 1.

m:raw_set(pos, value)

It receives a raw position at the underlying data pointer and a number. The given position is set to given number value. It is useful to combine stride and offset methods in order to compute the raw position.

> a = matrix(3,2, {1,2,3,4,5,6})
> = a
1 2
3 4
5 6
# Matrix of size [3,2] in row_major [0x144fce0 data= 0x144fd90]
> -- equivalent to a:set(2,1, 10)
> a:raw_set(a:offset() + a:stride(1)*1 + a:stride(2)*0, 10)
> = a
1 2
10 4
5 6
# Matrix of size [3,2] in row_major [0x144fce0 data= 0x144fd90]

NOTE! that the strides are multiplied by matrix position minus 1.

Sliding window iterator

For fast and easy matrix traversal, a C++ sliding window object is binded to Lua. It works similarly to dataset.matrix, but circularity and out-of-matrix default values are not supported. The object is constructed using the method sliding_window of matrix, and could be iterated using its method iterate():

> m = matrix(4,2,3):uniformf(-10,10,random(1234)) -- randomly initialized matrix
> for submat in m:sliding_window():iterate() do print(submat) end

# pos [1,1,1]
-6.16961 -0.0467267 2.44218
6.35677 -1.24545 2.24224
# Matrix of size [1,2,3] in row_major [0x253f160 data= 0x253dec0]

# pos [1,1,1]
5.70717 5.4272 5.59952
7.2134 -4.54815 -6.98726
# Matrix of size [1,2,3] in row_major [0x253fa40 data= 0x253dec0]

# pos [1,1,1]
-4.47071 -6.02962 6.03744
6.30326 9.16279 -6.82369
# Matrix of size [1,2,3] in row_major [0x2540230 data= 0x253dec0]

# pos [1,1,1]
7.51865 -7.67724 -2.84365
-9.74185 0.0199025 -0.263331
# Matrix of size [1,2,3] in row_major [0x25409c0 data= 0x253dec0]

It is possible to modify the default behavior giving this parameters to sliding_window method:

• offset: a Lua table with offset applied to the window in each coordinate (starting at 0).
• size: a Lua table with the window size for each coordinate.
• step: a Lua table with the step size at each coordinate (each value must be >= 1).
• numSteps: a Lua table with the number of steps in each coordinate (each value must be >= 1).
• orderStep: a Lua table with the traversal order of coordinates (starting at 1).

> m = matrix(4,2,3):uniformf(-10,10,random(1234))
> for w in m:sliding_window{ step={2,1,1}, size={1,1,2} }:iterate() do print(w) end
# pos [1,1,1]
-6.16961 -0.0467267
# Matrix of size [1,1,2] in row_major [0x9fdb90 data= 0x9cf2d0]

# pos [1,1,1]
-4.47071 -6.02962
# Matrix of size [1,1,2] in row_major [0x9fe0f0 data= 0x9cf2d0]

Manual iteration of the sliding_window is also possible using the following methods:

• matrix = sw:get_matrix([matrix]): returns the matrix generated by the window at its current position. It is possible to pass an optional argument, a destination matrix, so the computation effort is dismissed to constant. NOTE that this matrix must be created by a previous call to get_matrix over the same sliding_window.

• sw:next(): moves the window to the next position.

• sw:is_end(): returns true if the window has finished the matrix traversal.

> m = matrix(4,2,3):uniformf(-10,10,random(1234))
> wo = m:sliding_window{ step={2,1,1}, size={1,1,2} }
> while not wo:is_end() do m=wo:get_matrix(m) print(m) wo:next() end

Fast mathematical operations

This operations uses standard Lua math operators for friendly user interaction, but they work with BLAS API for best performance. However, all this operations return a new instantiated matrix, for best performance it is recommended to use directly the BLAS interface.

The operators binary +, -, *, /, and unary operators -, ^, are implemented as algebraic operations. The + and - operators only work when the matrices has the same sizes:

> a= matrix(3,3,3,{1,2,3,4,5,6,7,8,9,
                   10,11,12,13,14,15,16,17,18,
                   19,20,21,22,23,24,25,26,27})
> = a+a

# pos [1,1,1]
2 4 6
8 10 12
14 16 18

# pos [2,1,1]
20 22 24
26 28 30
32 34 36

# pos [3,1,1]
38 40 42
44 46 48
50 52 54
# Matrix of size [3,3,3] in row_major [0x1196d90 data= 0x1196e40]
> = a-a

# pos [1,1,1]
0 0 0
0 0 0
0 0 0

# pos [2,1,1]
0 0 0
0 0 0
0 0 0

# pos [3,1,1]
0 0 0
0 0 0
0 0 0
# Matrix of size [3,3,3] in row_major [0x1198d80 data= 0x1198a50]

The operator * only works with vectors or bi-dimensional matrices. If needed, you can rewrap the matrix data before the operation. Depending on the dimension of the two matrices, the multiplication could be:

• A dot product between two vectors: when the two matrices are unidimensional vectors, the first one a vector and the second one a column vector:

> a, b = matrix(4,{1,2,3,4}), matrix(4,1,{5,6,7,8})
> = a*b
70
# Matrix of size [1] in row_major [0xfa9230 data= 0xfc2300]

• An outter product between two vectors: when the first matrix is a column vector, and the second matrix is a unidimensional matrix or a bi-dimensional matrix (row or column vector).

> a = matrix(4,{1,2,3,4})
> b = matrix(4,1,{5,6,7,8})
> = b*a
5 10 15 20
6 12 18 24
7 14 21 28
8 16 24 32
# Matrix of size [4,4] in row_major [0x1001940 data= 0x1176a80]

• A matrix-vector product when the first matrix is a bi-dimensional matrix and the second is a vector. The output has the same number of dimensions as the given vector.

> a = matrix(2,2,{1,2,3,4})
> b = matrix(2,{5,6})
> = a*b
17 39
# Matrix of size [2] in row_major [0x105baa0 data= 0xfe80f0]
> b = matrix(1,2,{5,6})
> = a*b
17
39
# Matrix of size [2,1] in row_major [0x107e3c0 data= 0x107fb30]
> b = matrix(2,1,{5,6})
> = a*b
17
39
# Matrix of size [2,1] in row_major [0x10c4700 data= 0x10c6890]

• A matrix-matrix product when the two matrices are bi-dimensional and not vectors.

> a=matrix(3,2,{1,2,3,4,5,6})
> b=matrix(2,4,{1,2,3,4,5,6,7,8})
> = a*b
11 14 17 20
23 30 37 44
35 46 57 68
# Matrix of size [3,4] in row_major [0x1114270 data= 0x11165d0]

A multiplication by a scalar is also possible, if you multiply one matrix by one number.

> a=matrix(3,2,{1,2,3,4,5,6})
> = a*5
5 10
15 20
25 30
# Matrix of size [3,2] in row_major [0x10f2160 data= 0x10d14e0]

The component-wise operator / is allowed for division between matrix and a scalar, or between a scalar and a matrix.

The operator ^ is also allowed only with scalars.

The unary operator - is equivalent to multiply by -1.

matrix = m:scalar_add(number)

Adds to all the components, in-place, the given scalar number. Returns the caller matrix object.

matrix = m:div(scalar)

Produces the computation between the component-wise inversion of the matrix and the given scalar. This operation is done in-place.

> m = matrix(2,2,{1,2,3,4})
> m:div(1)
> = m
1      0.5
0.3333 0.25
# Matrix of size [3,2] in row_major [0x1cf2160 data= 0x10d15e0]

BLAS interface

The most efficient way to do operations if using the BLAS interface directly. All the methods are prepared to adjust the BLAS operations to the given matrices, so you don't need to be worried about strides and sizes.

All of this methods are in-place, so they modify the caller object, and returns it to allow operation sequences.

matrix = m:axpy(alpha, X)

The AXPY operation computes addition of vectors:

Y = alpha * X + Y

The method receives two positional parameters: the alpha scalar and the matrix X. The X and Y matrix sizes must be equals, and the number of dimensions is not a problem. This method interprets all the data as a sequence, calling several times to AXPY BLAS function if necessary:

> a = matrix(4,{1,2,3,4})
> b = matrix(4,{5,6,7,8})
> a:axpy(2.0, b)
> = a
11 14 17 20
# Matrix of size [4] in row_major [0x107e3c0 data= 0x1110970]
> a = matrix(2,2,2,{1,2,3,4,5,6,7,8})
> b = matrix(2,2,2,{9,10,11,12,13,14,15,16})
> a:axpy(1.0, b)
> = a

# pos [1,1,1]
10 12
14 16

# pos [2,1,1]
18 20
22 24
# Matrix of size [2,2,2] in row_major [0xfb1f40 data= 0x1056f00]

matrix = m:gemv{ beta, alpha, Y, X, trans_A}

The GEMV operation computes matrix-vector multiplication:

Y = beta * Y + alpha * op(A) * X

being Y the caller matrix (a vector), A another matrix, and X a vector (unidimensional matrix, or bi-dimensional with one row (or one column)), and beta and alpha are scalars. The op(A) is transposition operation.

The method receives a table with:

• A field, the other matrix.

• X field, the vector.

• alpha field, the scalar

• beta field, the other scalar.

• trans_A field, a boolean which indicates if the A matrix will be transposed or not. It is optional, by default is false.

> a = matrix(3,2,{1,2, 3,4, 5,6})
> b = matrix(2,{7,8})
> c = matrix(3)
> c:gemv{ A=a, X=b, alpha=2, beta=0 }
> = c
46 106 166
# Matrix of size [3] in row_major [0xfbeff0 data= 0xfaf640]

matrix = m:gemm{ beta, alpha, A, B, ... }

The GEMM operation computes matrix-matrix multiplication:

Y = beta * Y + alpha * op(A) * op(B)

being Y the caller matrix (a vector), A another matrix, and B a matrix, and beta and alpha are scalars. The op(A) and op(B) are transposition operations.

The method receives a table with:

• A field, the other matrix.

• B field, the vector.

• alpha field, the scalar

• beta field, the other scalar.

• trans_A field, a boolean which indicates if the A matrix will be transposed or not. It is optional, by default is false.

• trans_B field, a boolean which indicates if the B matrix will be transposed or not. It is optional, by default is false.

> a = matrix(3,2,{1,2, 3,4, 5,6})
> b = matrix(4,2,{7,8, 9,10, 11,12, 13,14})
> c = matrix(3,4):ones()
> c:gemm{ A=a, B=b, alpha=1, beta=1, trans_B=true}
> = c
24 30 36 42
54 68 82 96
84 106 128 150
# Matrix of size [3,4] in row_major [0x1452a20 data= 0x144cbf0]

matrix = m:ger{ X, Y, alpha }

The GER operation computes outter product of vectors:

Z = Z + alpha * X * Y'

being Z the caller matrix (a squared matrix), X and Y two vectors, and beta and alpha are scalars. The Y vector is transposed.

> a = matrix(3,{1,2,3})
> b = matrix(3,{4,5,6})
> c = matrix(3,3):zeros()
> c:ger{ X=a, Y=b, alpha=2 }
> = c
8 10 12
16 20 24
24 30 36
# Matrix of size [3,3] in row_major [0x1f06b20 data= 0x1f18080]

number = m:dot(matrix)

The DOT operation computes the dot-product of two vectors, the caller matrix and a given matrix. It returns a number.

> a = matrix(3,{1,2,3})
> b = matrix(3,{4,5,6})
> = a:dot(b)
32
# Matrix of size [1] in row_major [0x1f4ffe0 data= 0x2076e20]

matrix = m:scal(number)

The SCAL operation computes the multiplication of a matrix by a scalar.

> a = matrix(3,{1,2,3})
> a:scal(4)
> = a
4 8 12
# Matrix of size [3] in row_major [0x1f3b230 data= 0x201e9a0]

matrix = m:copy(matrix)

The COPY operation copies the content of a given matrix in the caller matrix object.

> a = matrix(3,3,{1,2,3,4,5,6,7,8,9})
> b = matrix(3,3):fill(5)
> a:copy(b)
> = a
5 5 5
5 5 5
5 5 5
# Matrix of size [3,3] in row_major [0x1f7e870 data= 0x1f49ef0]
> a = matrix(3,3,{1,2,3,4,5,6,7,8,9})
> b = matrix(2,2,{1,2,3,4})
> c = a:slice({2,1},{2,2})
> c:copy(b)
> = a
1 2 3
1 2 6
3 4 9
# Matrix of size [3,3] in row_major [0x1fb64e0 data= 0x1fbd600]

LAPACK interface

U,S,V = m:svd()

This method computes the Singular Values Decomposition of the caller matrix. The decomposition is computed in col_major order, so the returned matrices will be in col_major.

matrix = m:inv()

Computes the inverse of the caller matrix. Check that your matrix is not singular, otherwise the returned matrix won't be correct. It can work with row_major matrices, but internally they are transformed to col_major. The returned matrix will be in col_major.

matrix = m:pinv()

Computes the pseudo-inverse of the caller matrix, using the SVD method. Internally it uses col_major order, so the returned matrix will be in col_major.

Component-wise operations

This operations are applied in-place and over all the components of the caller matrix. If it is possible, the caller matrix is returned.

matrix = m:tan()

Computes the TAN function of all the components.

matrix = m:tanh()

Computes in-place the TANH function of all the components.

matrix = m:atan()

Computes in-place the ATAN function of all the components.

matrix = m:atanh()

Computes in-place the ATANH function of all the components.

matrix = m:sin()

Computes in-place the SIN function of all the components.

matrix = m:sinh()

Computes in-place the SINH function of all the components.

matrix = m:asin()

Computes in-place the ASIN function of all the components.

matrix = m:asinh()

Computes in-place the ASINH function of all the components.

matrix = m:cos()

Computes in-place the COS function of all the components.

matrix = m:cosh()

Computes in-place the COSH function of all the components.

matrix = m:acos()

Computes in-place the ACOS function of all the components.

matrix = m:acosh()

Computes in-place the ACOSH function of all the components.

matrix = m:abs()

Computes in-place the ABS function of all the components.

matrix = m:complement()

Computes in-place the complement function of all the components: X = 1 - X

matrix = m:log()

Computes in-place the LOG function of all the components.

matrix = m:log1p()

Computes in-place the LOG1P function of all the components.

matrix = m:plogp()

Computes in-place the p*log(p) operation over all components. It is useful to compute entropy related measures.

matrix = m:exp()

Computes in-place the EXP function of all the components.

matrix = m:pow(number)

Computes in-place the POWER of all the components by a given scalar.

matrix = m:sqrt()

Computes the SQRT function of all the components.

matrix = m:cmul(matrix)

Computes in-place a component-wise multiplication between the caller and a given matrix.

Matrix level operations

This operations are applied taking into account all the data at the matrix.

min,argmin = m:min()

Returns the minimum and its position in the matrix.

> a = matrix(3,4,{1,2,3,4,5,6,7,12,9,10,11,8})
> = a:min()
1	1

matrix,matrixInt32 = m:min(dim [, matrix[, matrixInt32]] )

Applies the min operator over the elements of the given dimension, and returns a matrix with the same number of dimensions, but with the size of dimension dim equals 1. The second matrix argument is optional, and if given, the returned matrix will be this second argument.

> a = matrix(3,4,{1,2,3,4,
>>                5,6,7,12,
>>                9,10,11,8})
> = a:min(1)
 1           2           3           4         
# Matrix of size [1,4] in row_major [0x1f06bb0 data= 0x1f06cb0]
> = a:min(2)
 1         
 5         
 8         
# Matrix of size [3,1] in row_major [0x1f07560 data= 0x1f06d90]

max,argmax = m:max()

Returns the maximum and its position in the matrix.

> a = matrix(3,4,{1,2,3,4,5,6,7,12,9,10,11,8})
> = a:max()
12	8

matrix,matrixInt32 = m:max(dim [, matrix[, matrixInt32] ] )

Applies the max operator over the elements of the given dimension, and returns a matrix with the same number of dimensions, but with the size of dimension dim equals 1. The second matrix argument is optional, and if given, the returned matrix will be this second argument.

> a = matrix(3,4,{1,2,3,4,
>>                5,6,7,12,
>>                9,10,11,8})
> = a:max(1)
 9           10          11          12        
# Matrix of size [1,4] in row_major [0x1f05500 data= 0x1f05600]
> = a:max(2)
 4         
 12        
 11        

matrix = m:lt(value)

This method computes in-place a comparison between all the components with the given value, and converts the component in 1 or 0 if it is less than the given value or not. If the value is another matrix, both matrices will be compared component-by-component.

matrix = m:gt(value)

This method computes in-place a comparison between all the components with the given value, and converts the component in 1 or 0 if it is greater than the given value or not. If the value is another matrix, both matrices will be compared component-by-component.

number = m:sum( )

Computes the sum of all the components of the caller matrix, and returns its value.

matrix = m:sum( number [, matrix] )

Receives a number indicating the dimension where the sum must be run, and returns a matrix with each possible sum of the given dimension. The second matrix argument is optional, and if given, the returned matrix will be this argument.

> m = matrix(2,2,{1,2,3,4})
> = m:sum(1)
 4           6         
# Matrix of size [1,2] in row_major [0x19e0620 data= 0x19d2480]
> = m:sum(2)
 3         
 7         
# Matrix of size [2,1] in row_major [0x19e0a40 data= 0x19d3b90]

number = m:norm2()

The NORM2 operation computes the euclidean norm of the caller matrix. It returns a number.

> a = matrix(2,2,2,{1,2,3,4,5,6,7,8})
> = a:norm2()
14.282856941223

Matrix dictionary (hash set): matrix.dict

The class matrix.dict is a C++ binded class, which allows to store several matrices in a hash set. They are indexed by strings, and basic component-wise math operations are implemented, allowing to operate between two matrix.dict instances.

obj = matrix.dict( [table] )

The constructor allows to build an empty matrix dictionary, or could receive a Lua table with the dictionary fields.

> obj = matrix.dict() -- an empty dictionary
> obj = matrix.dict{ a = matrix(2,3):linear() }
> = type( obj('a') )
matrix

matrix = obj(name)

Returns the matrix stored with the given name string, by using the __call metamethod. If the given name is not associated with any matrix, nil will be returned.

matrix = obj:find(name)

Equivalent to the previous one.

obj[name] = matrix

Allows to modify the matrix associated with the given name.

obj = obj:insert(name, matrix)

Equivalent to the previous one. Returns the caller object`

table = obj:keys()

Returns a table will all the string keys of the dictionary.

Lua iterator = obj:iterate()

Retruns a Lua iterator, which allows to perform loops over all the matrices:

> for key,mat in obj:iterate() do print(key) print(mat) end

Lua iterator = pairs(obj)

The same as the previous one.

> for key,mat in pairs(obj) do print(key) print(mat) end

string = obj:to_lua_string(format)

Returns a callable Lua string which builds a serialized copy of the caller object. It receives a format string, which could be ascii or binary.

another = obj:clone()

Returns a deep copy of the caller object.

another = obj:clone_only_dims()

Returns a deep copy of the caller object, but without copying the matrix data content, only cloning the matrix dimension sizes.

Implemented operations

The following list of operations are implemented to be executed over all the contained matrices:

• number = obj:size()

• obj = obj:fill(number)

• obj = obj:ones()

• obj = obj:zeros()

• obj = obj:axpy(number, matrix.dict)

• number = obj:dot(matrix.dict)

• obj = obj:copy(matrix.dict)

• obj = obj:scalar_add(number)

• obj = obj:complement()

• obj = obj:pow(number)

• obj = obj:scal(number)

• obj = obj:inv()

• obj = obj:sqrt()

• obj = obj:exp()

• obj = obj:plogp()

• obj = obj:log1p()

• obj = obj:cos()

• obj = obj:cosh()

• obj = obj:acos()

• obj = obj:acosh()

• obj = obj:tan()

• obj = obj:tanh()

• obj = obj:atan()

• obj = obj:atanh()

• obj = obj:sin()

• obj = obj:sinh()

• obj = obj:asin()

• obj = obj:asinh()

Other kind of matrices

Currently is possible to use complex, double, int32 and char matrices, supporting load and save, matrix structural methods, and some of them also support mathematical operations:

• matrixComplex: fully working matrix type, with almost all the methods described above.
• matrixDouble: partial working matrix type, only allow structural methods (explained at MatFormat section).
• matrixInt32: partial working matrix type, only allow structural methods (explained at MatFormat section).
• matrixChar: partial working matrix type, only allow structural methods (explained at MatFormat section).

In all cases, you could use april_help to ask which methods are available. Complex, Double and Int32 matrices implements a method to_float() which converts the given object to a standard matrix with float numeric precission. The matrixChar type implements a method to_string_table.

matrixComplex

The constructor of a matrixComplex receives a table with complex numbers (see utils section). A complex number uses float single precission resolution for real and imaginary part:

> -- using strings which are converted to complex numbers (slow performance)
> m = matrixComplex(2,2, { "1+1i", "2+2i", "3+2i", "4+1i" })
> = m
        1+1i         2+2i
        3+2i         4+1i
# MatrixComplex of size [2,2] in row_major [0x24d52c0 data= 0x24d4a00]
>
> -- using directly complex numbers
> m = matrixComplex(2,2, { complex(1,1), complex(2,2), complex(3,2), complex(4,1) })
> = m
        1+1i         2+2i
        3+2i         4+1i
# MatrixComplex of size [2,2] in row_major [0x24d6550 data= 0x24d6650]

Besides the standard matrix methods, the matrixComplex implements the following:

• caller = m:conj() computes the conjugate in-place, modifying the caller matrix, and returning the caller matrix instance.
• matrix = m:real() returns the real part of the caller matrixComplex.
• matrix = m:img() returns the imaginary part of the caller matrixComplex.
• matrix = m:abs() returns the modulus of the polar form of matrixComplex.
• matrix = m:angle() returns the angle of the polar form of matrixComplex.
• matrix = m:to_float() converts the caller matrix in a matrix object which has one additional dimension. This additional dimension has always size 2, and keeps the real and imaginary parts of the caller matrixComplex. If the caller matrix is ordered in row_major, then the additional dimension will be the last, otherwise (col_major order), the additional dimension will be the first. Note that the returned matrix and the matrixComplex caller references the same memory pointer.

matrixDouble

matrixDouble is the type of matrices for double data. This kind of matrices don't accept mathematical operations, but yes structural operations as select, slice, etc.

It is also possible to convert this matrix in an standard float matrix using the method to_float([false]), which returns the same matrix but casting the data to float. It receives an optional boolean argument, which indicates if the resulting matrix will be in col_major. If not given, it is taken as false.

> m = matrixDouble(2,3,{1,2,3,4,5,6})
> = m
          1           2           3
          4           5           6
# MatrixDouble of size [2,3] [0x2512c70 data= 0x251ad70]
> = m:to_float(true)
 1           2           3         
 4           5           6         
# Matrix of size [2,3] in col_major [0x25142b0 data= 0x251adf0]

matrixInt32

matrixInt32 is the type of matrices for integer data. This kind of matrices don't accept mathematical operations, but yes structural operations as select, slice, etc.

It is also possible to convert this matrix in an standard float matrix using the method to_float([false]), which returns the same matrix but casting the data to float. It receives an optional boolean argument, which indicates if the resulting matrix will be in col_major. If not given, it is taken as false.

> m = matrixInt32(2,3,{1,2,3,4,5,6})
> = m
          1           2           3
          4           5           6
# MatrixInt32 of size [2,3] [0x2512c70 data= 0x251ad70]
> = m:to_float(true)
 1           2           3         
 4           5           6         
# Matrix of size [2,3] in col_major [0x25142b0 data= 0x251adf0]

matrixChar

matrixChar is the type of matrices for char data. This kind of matrices don't accept mathematical operations, but yes structural operations as select, slice, etc.

Exists an special method, to_string_table(), which converts the matrix in a table of strings, concatenating the chars in row_major order.

> m = matrixChar(2,2, { "h","ola" })
> = m
[1,1] = h
[1,2] = o
[2,1] = l
[2,2] = a
# MatrixChar of size [2,2] [0x12c3310 data= 0x12c3410]
> = unpack(m:to_string_table())
ho	la