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
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
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
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] ^
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 ]
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
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):
Using an expression does also work. If n is a column vector, Euler will plot all functions.
>plot2d("x^n",a=0,b=1);
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
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):
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
>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):
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
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