iconEuler Home

The Matrix Language of Euler

The matrix language is the core of the numerical part of Euler. It
allows generating tables of functions without the need to write loops.

Of course, it can also be used to compute expressions with vectors and
matrices. But this is explained in detail in the introduction to
linear algebra.

See: 04 - Introduction to Linear Algebra

First of all, a matrix can be defined element by element. Separate the
elements with ",", and the rows with ";".

Let us enter a 3x3 matrix and assign it to the variable A.
>shortformat; A:=[1,2,3;4,5,6;7,8,9]
          1           2           3 
          4           5           6 
          7           8           9 
Incomplete rows will be filled with 0.
>[1,2,3;3,4]
          1           2           3 
          3           4           0 
A matrix with one row is a row vector. The matrix language makes a
difference between row and column vectors.
>v:=[1,2,3]
[ 1  2  3 ]
The : operator creates equal spaced vectors. By default the
increment is 1.
>1:3
[ 1  2  3 ]
But the increment can also be any other value.
>1:0.2:2
[ 1  1.2  1.4  1.6  1.8  2 ]
Even negative.
>10:-1:1
[ 10  9  8  7  6  5  4  3  2  1 ]
For a narrower output use shortformat, shortestformat, or fracformat.
In the following command, we use four places for the output with no
decimal fraction. The function intrandom generates integer random
numbers.
>format(4,0); intrandom(10,10,6), shortformat;
  3   1   2   6   3   1   1   2   5   2 
  6   1   2   4   2   1   6   6   2   4 
  3   4   5   6   2   6   6   2   4   1 
  2   3   5   3   4   5   1   2   3   3 
  3   1   6   3   1   2   5   4   2   4 
  3   3   5   5   5   3   5   4   1   1 
  5   5   3   2   3   2   1   6   1   4 
  4   6   5   5   4   2   2   1   2   1 
  3   3   4   5   3   5   1   3   3   1 
  4   2   2   2   5   2   6   4   4   5 
For a column vector, use a 1xn matrix.
>[1;2;3]
          1 
          2 
          3 
Or transpose a row vector.
>[1,2,3]'
          1 
          2 
          3 

Basic Rule

The most important feature of the matrix language is that all
operators and functions are applied element for element to the matrix
or vector.

Since v was 1,2,3, and v*v multiplies each element of v with itself,
we get all elements squared. So the result w has the property

01 - Introduction to the Matrix Language

Note that this is different from the matrix product or the inner
product!
>v*v
[ 1  4  9 ]
This is the same as v^2.
>v^2
[ 1  4  9 ]
Not only operators, but almost all mathematical functions are
applied element for element to a matrix or vector.
>sqrt(v)
[ 1  1.4142  1.7321 ]
Always remember that * is not the matrix product.
>A*A
          1           4           9 
         16          25          36 
         49          64          81 
The matrix product uses the "." operator.
>A.A
         30          36          42 
         66          81          96 
        102         126         150 
Let us define a smaller matrix.
>B:=[1,2;3,4]
          1           2 
          3           4 
B^(-1) is not the inverse of B! By the mapping rule, it will
invert all elements of B.
>B^(-1)
          1         0.5 
    0.33333        0.25 
To compute the inverse, we have to use the function inv.
>inv(B)
         -2           1 
        1.5        -0.5 
Alternatively, there is the function matrixpower.
>matrixpower(B,-1)
         -2           1 
        1.5        -0.5 
We can print this matrix using fractional output.
>fracprint(inv(B));
       -2         1 
      3/2      -1/2 
Or we can go use the fracformat for all outputs.
>fracformat; inv(B), shortformat;
                 -2                   1 
                3/2                -1/2 
Let us check the inverse.
>B.inv(B)
          1           0 
          0           1 
Euler is a numerical program. So we do not get the exact inverse.
But to see the error, we have to choose a longer format.
>longestformat; B.inv(B), shortformat;
     0.9999999999999998  1.110223024625157e-016 
-2.220446049250313e-016                       1 
Of course, Maxima can be used for these linear algebra problems. We
can define the matrix for Euler and Maxima with &:=, and then use it
in symbolic expressions.
>B &:= [1,2;3,4]; &invert(B)
                             [ - 2   1  ]
                             [          ]
                             [  3     1 ]
                             [  -   - - ]
                             [  2     2 ]

The result can be evaluated numerically in Euler. For more
information about Maxima, see the introduction to Maxima.
>&invert(B)()
         -2           1 
        1.5        -0.5 
Euler has also a powerful function xinv, which makes a bigger effort,
and gets more exact results. Note, that with &:= the matrix B will be
symbolic in symbolic expressions, and numerical in numerical
expressions.
>longestformat; B.xinv(B), shortformat;
                      1                       0 
                      0                       1 
To transpose a matrix we use the apostroph.
>v'
          1 
          2 
          3 
So we can compute matrix A times vector b.
>A.v'
         14 
         32 
         50 
Note that v is still a row vector. So v'.v is different from
v.v'.
>v'.v
          1           2           3 
          2           4           6 
          3           6           9 
v.v' computes the norm of v squared for row vectors v. The result is a
1x1 vector, which works just like a real number.
>v.v'
[ 14 ]
There is also the function norm (along with many other function of
Linear Algebra).
>norm(v)^2
14

Sub-Matrices and Matrix-Elements

To access a matrix element, use the bracket notation.
>A[2,2]
5
Euler can access a complete line of a matrix.
>A[2]
[ 4  5  6 ]
In case of row vectors, this returns an element of the vector.
>v[2]
2
To make sure, you get the first row for a 1xn and a mxn matrix,
specify all columns using an empty second index.
>A[2,], v[1,]
[ 4  5  6 ]
[ 1  2  3 ]
If the index is a vector of indices, Euler will return the
corresponding rows of the matrix.

Here we want the first and second row of A.
>A[[1,2]]
          1           2           3 
          4           5           6 
We can even reorder A this way. To be precise, we do not change
A here, but compute a reordered version of A.
>A[[3,2,1]]
          7           8           9 
          4           5           6 
          1           2           3 
The index trick works with columns too.

This example selects all rows of A and the second and third
column.
>A[1:3,2:3]
          2           3 
          5           6 
          8           9 
For abbreviation ":" denotes all row or column indices.
>A[:,3]
          3 
          6 
          9 
Alternatively, leave the first index open.
>A[,2:3]
          2           3 
          5           6 
          8           9 
We can also get the last line of A.
>A[-1]
[ 7  8  9 ]
But note that A[0] is not defined! The result would be an empty
matrix, or an error, depending on a global flag.

If the size of the two matrices involved in an operation is different,
or one of the operands is a number, Euler does the natural thing,
and expands both operands to the same matrix size.
>2*v
[ 2  4  6 ]
Likewise, A*v will expand the row vector v, by duplicating it three
times. So, v gets the same dimensions as the 3x3-matrix A.
Effectively, the columns of A are multiplied with the corresponding
elements of v.

Again, remember that * is not the matrix product!
>A*v
          1           4           9 
          4          10          18 
          7          16          27 
The same can be done with the rows of A.
>A*[1;2;1]
          1           2           3 
          8          10          12 
          7           8           9 
The order of multiplication does not matter here.
>[1;2;1]*A
          1           2           3 
          8          10          12 
          7           8           9 
The transpose of v is a column vector. So v'*v will return the
matrix v[i]*v[j].
>v'*v
          1           2           3 
          2           4           6 
          3           6           9 
v*v' is the same, but both are not equal to v.v'!
>v.v'
[ 14 ]
Now let us change elements of A by assigning a submatrix of A to some
value. This does in fact change the stored matrix A.
>A[1,1]:=4
          4           2           3 
          4           5           6 
          7           8           9 
We can also assign a value to a row of A.
>A[1]:=[-1,-1,-1]
         -1          -1          -1 
          4           5           6 
          7           8           9 
We can even assign to a submatrix.
>A[1:2,1:2]:=[5,6;7,8]
          5           6          -1 
          7           8           6 
          7           8           9 
Moreover, some shortcuts are allowed.
>A[1:2,1:2]:=-1
         -1          -1          -1 
         -1          -1           6 
          7           8           9 
A warning: Indices out of bounds return empty matrices, or an error
message, depending on a system setting. The default is an error
message. Remember, however, that negative indices may be used to
access the elements of a matrix counting from the end.
>A[4]
Index 4 out of bounds!

Error in :
A[4]
    ^

Other Matrix Functions

To build a matrix, we can stack one matrix on top of another. If both
do not have the same number of columns, the shorter one will be
filled with 0.
>v_v
          1           2           3 
          1           2           3 
Likewise, we can attach a matrix to another side by side, if both
have the same number of rows.
>A|v'
         -1          -1          -1           1 
         -1          -1           6           2 
          7           8           9           3 
If they do not have the same number of rows the shorter matrix
is filled with 0.

There is an exception to this rule. A real number attached to
a matrix will be used as a column filled with that real number.
>A|1
         -1          -1          -1           1 
         -1          -1           6           1 
          7           8           9           1 
To get the size of A, we can use the following functions. The
function length returns the maximum of the rows and columns. It is
useful for row or column vectors.
>C=zeros(2,4); rows(C), cols(C), length(C), size(C),
2
4
4
[ 2  4 ]
There are many other functions, which generate matrices.
>ones(2,2)
          1           1 
          1           1 
Also a matrix of random numbers can be generated with random (uniform
distribution) or normal (Gauß distribution).
>random(2,2)
   0.021243     0.46542 
    0.18176     0.14257 
Here is another useful function, which restructures the elements of a
amtrix into another matrix.
>redim(1:9,3,3)
          1           2           3 
          4           5           6 
          7           8           9 
With the following function, we can use this and the dup function to
write a rep() function, which repeats a vector n times.
>function rep(v,n) := redim(dup(v,n),1,n*cols(v))
Let us test.
>rep(1:3,5)
[ 1  2  3  1  2  3  1  2  3  1  2  3  1  2  3 ]
The function multdup doublicates elements of a vector.
>multdup(1:3,5), multdup(1:3,[2,3,2])
[ 1  1  1  1  1  2  2  2  2  2  3  3  3  3  3 ]
[ 1  1  2  2  2  3  3 ]

Examples

We compute i*j for i,j from 1 to 10. The trick is to multiply 1:10
with its transpose. The matrix language of Euler automatically
generates a table of values.
>fixedformat(4,0); (1:10)*(1:10)', shortformat;
  1   2   3   4   5   6   7   8   9  10 
  2   4   6   8  10  12  14  16  18  20 
  3   6   9  12  15  18  21  24  27  30 
  4   8  12  16  20  24  28  32  36  40 
  5  10  15  20  25  30  35  40  45  50 
  6  12  18  24  30  36  42  48  54  60 
  7  14  21  28  35  42  49  56  63  70 
  8  16  24  32  40  48  56  64  72  80 
  9  18  27  36  45  54  63  72  81  90 
 10  20  30  40  50  60  70  80  90 100 
Using the matrix language, we can generate several tables of several
functions at once.

In this example, we compute t[j]^i for i from 1 to n. We get a matrix,
where each row is a table of t^i for one i. I.e., the matrix has the
elements

01 - Introduction to the Matrix Language

This can be plotted in one plot in one command. The plot command plots
one function for each row of the matrix.
>n=(1:10)'; t=0:0.01:1; plot2d(t,t^n):

01 - Introduction to the Matrix Language

Using an expression does also work. If n is a column vector, Euler
will plot all functions.
>plot2d("x^n",a=0,b=1);

Operators and the Matrix Language

Euler has comparison operators, like "==", which checks for equality.

We get a vector of 0 and 1, where 1 stands for true.
>t=(1:10)^2; t==25
[ 0  0  0  0  1  0  0  0  0  0 ]
From such a vector, "nonzeros" selects the non-zero elements.

In this case, we get the indices of all elements greater than
50.
>nonzeros(t>50)
[ 8  9  10 ]
Of course, we can use this vector of indices to get the corresponding
values in t.
>t[nonzeros(t>50)]
[ 64  81  100 ]
As an example, let us find all squares of the numbers 1 to 1000,
which are 5 modulo 11 and 3 modulo 13.
>t:=1:1000; nonzeros(mod(t^2,11)==5 && mod(t^2,13)==3)
[ 4  48  95  139  147  191  238  282  290  334  381  425  433  477  524
568  576  620  667  711  719  763  810  854  862  906  953  997 ]
Euler is not completely effective for integer computations. It uses
double precision floating point internally. However, it is often very
useful.

We can check for primality. Let us find out, how many squares plus 1
are primes.
>t:=1:1000; length(nonzeros(isprime(t^2+1)))
112

Mapping Functions to Matrices

Almost all functions in Euler work for matrix and vector input
too, whenever this makes sense. 

E.g., the sqrt function computes the square root of all elements
of the vector or matrix.
>sqrt(1:3)
[ 1  1.4142  1.7321 ]
So you can easily create a table of values. This is one way to
plot a function (the alternative uses an expression).
>x=1:0.01:5; y=log(x)/x^2; plot2d(x,y):

01 - Introduction to the Matrix Language

For the following, we assume, that your already are acquainted
with the programming language contained in Euler. 

Functions will only work for vectors, if they return the result
using a simple expression, which works for matrices.
>function f(x) := x^3-x;
>f(0:0.1:1)
[ 0  -0.099  -0.192  -0.273  -0.336  -0.375  -0.384  -0.357  -0.288
-0.171  0 ]
But if conditional statements or other more complicated statements are
involved, the function must be mapped. The easiest way is to use the
keyword "map" in the function definition.

Let us define

01 - Introduction to the Matrix Language

>function map f(x:real scalar) ...
  if x>0 then return x^3
  else return x^2
  endif;
endfunction
>f(-1:0.5:1)
[ 1  0.25  0  0.125  1 ]
The following method works with all functions.
>fmap(-1:0.5:1)
[ 1  0.25  0  0.125  1 ]
The plot functions, as well as other functions use map by default. It
can be switched off with maps=0.

So the following would work, even if f was not mapped.
>plot2d("f(x)",-1,1):

01 - Introduction to the Matrix Language

Again, remember that sqrt is not a matrix function. it simply maps to
all elements of the matrix.
>sqrt(B)
          1      1.4142 
     1.7321           2 
To apply a function to a matrix, we need to compute its eigenvalue
decomposition. The function matrixfunction does this.
>C:=matrixfunction("sqrt",B)
       0.55369+0.46439i        0.80696-0.21243i 
        1.2104-0.31864i         1.7641+0.14575i 
Check it.
>C.C
                   1+0i                    2+0i 
                   3+0i                    4+0i 

Sorting and Shuffling

The function sort sorts a row vector.
>sort([5,6,4,8,1,9])
[ 1  4  5  6  8  9 ]
It is often necessary to know the indices of the sorted vector in the
original vector. This can be used to reorder another vector in the
same way.

Let us shuffle a vector.
>v=shuffle(1:10)
[ 6  4  7  10  5  2  1  9  8  3 ]
The indices contain the proper order of v.
>{vs,ind}=sort(v); v[ind]
[ 1  2  3  4  5  6  7  8  9  10 ]
Note that sort returns two results, and thus cannot be used in
functions easily, since multiple return values count for multiple
arguments. 

You need to store the result to variable first, or use extra brackets.
>sum((sort(v)))
55
This works for string vectors too.
>s=["a","d","e","a","aa","e"]
a
d
e
a
aa
e
>{ss,ind}=sort(s); ss
a
a
aa
d
e
e
As you see, the position of double entries is somewhat random.
>ind
[ 4  1  5  2  6  3 ]
The function unique returns a sorted list of unique elements of a
vector.
>intrandom(1,10,10), unique(%)
[ 9  6  4  5  10  4  4  3  6  4 ]
[ 3  4  5  6  9  10 ]
This works for string vectors too.
>unique(s)
a
aa
d
e

Euler Home