Basics¶
In this section, the user will find informations about basic manipulations on the matrix object : creation (and destruction), etc. matrix is defined within the namespace castor meaning that in order to call the functions, their name should be preceeded by castor::. In order to use directly their name, the user should add
using namespace castor;
in the preamble of the .cpp file. A minimum main file would then look like this:
#include "castor/matrix.hpp"
using namespace castor;
int main()
{
matrix<double> A;
// Write your code here
return 0;
}
For some more advanced features, see the corresponding section Advanced.
The user may also refer to the examples in the demo/* subdirectories of the main directory of the castor project.
Matrix creation (and destruction)¶
There are two main paths to create a matrix object. The first way is to call one of the constructors explicitly. If no value is specified, the object will be filled with zeros. The second way it to call a builder (see From builder and Display).
From matrix constructor¶
Remark: In the following, we will use a lot the disp function which is meant to produce a formatted output of the content of a matrix object.
The code below creates an empty matrix of type int.
matrix<int> A;
disp(A,2);
Matrix 0x0 of type 'i' (0 B):
-empty-
Remark: If the using namespace castor clause was not added in the preamble as explained at the beginning of this section, the code above becomes
castor::matrix<int> A;
castor::disp(A,2);
// etc.
By passing the value “2” as a second argument to disp, we can see that it is a matrix of size 0x0 of integer type whose data size is 0 Byte.
By passing as argument a single value of type T, a singleton matrix is created.
matrix<int> A(M_PI);
disp(A,2);
Matrix 1x1 of type 'i' (4 B):
3
Here, A has been declared as an matrix of integers but a double containing the value of pi was passed as argument. As a consequence, it was cast to an int, thus the value 3. Note that 4 Bytes is the size of an integer in C++ when standard compilation options are used.
Next, we initialize a matrix using an initialization-list. By passing a single list of values as argument to the constructor, a line-matrix is created. By passing a list of a list, a matrix is created whose number of lines is the number of elements in the outer list and the number of columns is the number of elements in the inner lists. Please be careful that the number of elements in the inner list should be the same for all. These two options are illustrated below.
matrix<float> A({1,2,3,4,5}); // matrix of floats
matrix<> B({{1,2,3},{4,5,6}}); // matrix of doubles
disp(A,2);
disp(B,2);
Matrix 1x5 of type 'f' (20 B):
1.00000 2.00000 3.00000 4.00000 5.00000
Matrix 2x3 of type 'd' (48 B):
1.00000 2.00000 3.00000
4.00000 5.00000 6.00000
Finally, it is possible to create a matrix by giving its dimensions and a fill-value. By default, the matrix is filled with 0s. In the example below, we create a 2x3 matrix filled with the value 4, then we modify one of the entries.
matrix<> A(2,3,4.);
disp(A,2);
A(1,2) = -0.5;
disp(A,2);
Matrix 2x3 of type 'd' (48 B):
4.00000 4.00000 4.00000
4.00000 4.00000 4.00000
Matrix 2x3 of type 'd' (48 B):
4.00000 4.00000 4.00000
4.00000 4.00000 -0.50000
Please refer to the constructors list in the class matrix description.
From builder¶
We describe now some of the useful builders for the matrix class.
The code below creates a 2x3 matrix of double filled with zeros.
matrix<> A = zeros(2,3);
disp(A);
1.0000 1.0000 1.0000
1.0000 1.0000 1.0000
Remark: this is equivalent to calling explicitly the matrix constructor.
The code below creates a 1x10 matrix of double initialized with linear spaced values :
matrix<> A = linspace(0,1,10);
disp(A,2);
Matrix 1x10 of type 'd' (80 B):
0 0.11111 0.22222 ... 0.77778 0.88889 1.00000
For a 2x2 random matrix of float, use
matrix<float> A = rand<float>(2);
disp(A,2);
Matrix 2x2 of type 'f' (16 B):
0.84019 0.39438
0.78310 0.79844
This last result may differ depending on your random number generator.
Notes :
Matrices and vectors are objects of the matrix template class. A vector is considered as a (1xn) size by default or (nx1).
The template argument of class matrix is double by default. It is possible to specify type both for matrix constructors and builders.
Clear a matrix¶
If for some reason the content of a matrix needs to be cleared (for example, free some RAM), there are to possibilities. The first solution (the clean one) is to call the clear function.
auto A=rand(5);
disp(A,2);
clear(A);
disp(A,2);
Matrix 3x3 of type 'd' (72 B):
0.84019 0.39438 0.78310
0.79844 0.91165 0.19755
0.33522 0.76823 0.27777
Matrix 0x0 of type 'd' (0 B):
-empty-
The second one is to assign an empty matrix in place of an existing one.
auto A=rand(3);
disp(A,2);
A = {};
disp(A,2);
Matrix 3x3 of type 'd' (72 B):
0.84019 0.39438 0.78310
0.79844 0.91165 0.19755
0.33522 0.76823 0.27777
Matrix 0x0 of type 'd' (0 B):
-empty-
Display¶
A very useful function is the disp function. It produces a formatted output of the content of a matrix object with additional informations. Let us create a 2x2 random matrix.
auto A = rand<>(2);
The simplest call to disp displays the raw content without additional informations.
disp(A);
0.84019 0.39438
0.78310 0.79844
In fact, this is equivalent to calling disp(A,1). The second (optional) argument determines the level of informations to be displayed. disp(A,0) will produce the same output as disp(A) but no end-of-line character is added to the output. disp(A,2) will add informations on the size, the type and the RAM storage of the matrix, as illustrated before.
The third argument to disp is the output stream (std::ostream) which by default is the standard output std::cout. Finally, the user may specify two additional arguments which are the number of lines and columns which need to be displayed. By default, their value is 3 meaning that the first 3 and last 3 element of each direction are displayed.
auto A = rand(10);
disp(A,2);
Matrix 10x10 of type 'd' (800 B):
0.84019 0.39438 0.78310 ... 0.76823 0.27777 0.55397
0.47740 0.62887 0.36478 ... 0.71730 0.14160 0.60697
0.01630 0.24289 0.13723 ... 0.10881 0.99892 0.21826
...
0.53161 0.03928 0.43764 ... 0.73853 0.63998 0.35405
0.68786 0.16597 0.44010 ... 0.89337 0.35036 0.68667
0.95647 0.58864 0.65730 ... 0.81477 0.68422 0.91097
Now we modifiy a little bit the format.
disp(A,2,std::cout,4,2);
Matrix 10x10 of type 'd' (800 B):
0.84019 0.39438 ... 0.27777 0.55397
0.47740 0.62887 ... 0.14160 0.60697
0.01630 0.24289 ... 0.99892 0.21826
0.51293 0.83911 ... 0.29252 0.77136
...
0.23828 0.97063 ... 0.51254 0.66772
0.53161 0.03928 ... 0.63998 0.35405
0.68786 0.16597 ... 0.35036 0.68667
0.95647 0.58864 ... 0.68422 0.91097
Note : disp can also display the content of other variables :
disp("pi value is :");
disp(M_PI);
pi value is
3.14159
Size and indexing¶
We describe now a few useful functions to begin manipulating matrices.
Size, length, numel¶
The size function returns the two-element vector containing the number of rows and columns in the matrix. The result is also a matrix.
matrix<> A = eye(3,4);
disp(size(A),2)
Matrix 1x2 of type 'm' (16 B):
3 4
If a dimension is specified, size returns only the length of the specified dimensions :
matrix<> A = eye(3,4);
disp(size(A, 1));
disp(size(A, 2));
3
4
The length and numel functions returns respectively the maximum length and the number of elements in the matrix :
matrix<> A = eye(3,4);
disp(length(A));
disp(numel(A));
4
12
Accessing elements¶
The elements of a matrix can be accessed using either linear or bilinear indexing.
Linear indexing consists in accessing the n-th element of the matrix in the order the data is stored. Since matrix uses a row-major layout, the rows of the matrix are concatenated one after the other.
matrix<> A = {{1,2,3},
{4,5,6},
{7,8,9}};
disp(A(5));
6
The 5-th element of A thus holds the value 6 (the 6-th holds 7, etc.). Linear indexing is particularly useful when accessing the elements of a one-dimensional matrix (a vector).
Remark: The index numbering follows the C/C++ convention meaning that the indexes start a 0 and ends at n-1 where n would be a dimension of the matrix (see size).
Bilinear indexing is the natural way to access the elements of a matrix.
disp(A(1,2));
6
Help¶
The function help allows you to display the documentation of a function at runtime. You have to give the complete path to the header file matrix.hpp in the documentationFiles variable.
documentationFiles =
{
"/complete/path/to/matrix.hpp"
};
help("size");
============================ DOCUMENTATION ============================
Help on "size":
Size of array.
S = size(A) for m-by-n matrix A returns the two-element vector [m,n]
containing the number of rows and columns in the matrix.
S = size(A,dim) returns the lengths of the specified dimensions dim.
Example(s):
matrix<> A = {{1,2,3},{4,5,6}};
disp(size(A));
disp(size(A,1));
disp(size(A,2));
See also:
length, numel.
=======================================================================
Basic operations¶
The matrix class is designed to be as easy of use as Matlab or Numpy arrays. As a consequence, many operators have been overloaded. We will describe here some basic manipulations with few of all of the available operators. They can be discovered at Operators.
First, we create two matrices A and B, we multiply the first one by M_PI and we add them.
// create two 'double' matrices
auto A = eye(2);
auto B = eye(2);
disp(A,2);
disp(B,2);
A *= M_PI;
disp(A,2);
auto C = A + B;
disp(C,2);
Matrix 2x2 of type 'd' (32 B):
1.00000 0
0 1.00000
Matrix 2x2 of type 'd' (32 B):
1.00000 0
0 1.00000
Matrix 2x2 of type 'd' (32 B):
3.14159 0
0 3.14159
Matrix 2x2 of type 'd' (32 B):
4.14159 0
0 4.14159
Then, we create an orthogonal matrix L and we compute D = C - L*L'*C. Orthogonal matrices are such that the matrix product with their transpose yields the identity matrix. Therefore, the result should be a null-matrix.
double theta = 0.2;
matrix<> L = {
{std::cos(theta),-std::sin(theta)},
{std::sin(theta),std::cos(theta)}
};
auto D = C - mtimes(L,mtimes(transpose(L),C));
disp(D,2);
Matrix 2x2 of type 'd' (32 B):
0 0
0 0
We obtain the expected result. Note that the matrix-matrix product is computed using mtimes. Indeed, the * operator does not compute the matrix-matrix product but a term-by-term product.
auto A = ones(2);
auto B = matrix<>(2,2,2.);
disp(A*B,2);
Matrix 2x2 of type 'd' (32 B):
2.00000 2.00000
2.00000 2.00000