Engee documentation
Notebook

Unit matrix

In this script we will describe several ways to create unit matrices and show how to compare these matrices in terms of their memory footprint.

Task description

Create a matrix whose main diagonal is filled with units (the other elements are either zero or missing). Assume that such a matrix can be not only square, but ideally also multidimensional.

Using loops

A unit matrix can be created using loops as follows.

In [ ]: 
Pkg.add(["LinearAlgebra"])
In [ ]: 
n = 14
m = 11

eye = zeros( m, n )

for i in 1:n 
    if m < n
        if i <= m
            eye[i,i] = 1
        end
    else
        eye[i,i] = 1
    end
end
In [ ]: 
gr()
heatmap( eye, yflip=:true, c=:Reds )
Out[0]:

Through a short cycle

Putting units along the main diagonal can be done in a compressed cycle.

In [ ]: 
eye_array = zeros( 5, 5 );

[eye_array[ i,i ] = 1 for i in 1:minimum( size(eye_array) )];

heatmap( eye_array, yflip=:true, c=:Purples )
Out[0]:

This code can be adapted for the multidimensional case by using the function ndims, which returns the number of dimensions of the original matrix, whose main diagonal we fill with units.

In [ ]: 
eye3D = zeros(5,4,3)

[eye3D[ repeat( [i], ndims(eye3D) )... ] = 1 for i in 1:minimum( size(eye3D) )];

display( eye3D )

5×4×3 Array{Float64, 3}:
[:, :, 1] =
 1.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0

[:, :, 2] =
 0.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0

[:, :, 3] =
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0

But this custom code will quickly become difficult to maintain if the system requirements grow. It is better to use additional libraries, especially since they are also built into the Engee environment.

LinearAlgebra libraries: object I

The most economical way to create such a matrix is to use the command I.

In [ ]: 
using LinearAlgebra
I
Out[0]:
UniformScaling{Bool}
true*I

This object occupies very little memory, behaves like a unit matrix, can be dominated by any object, and can be converted to a matrix using the Matrix() command.

In [ ]: 
Matrix( 5I, 3, 3)
Out[0]:
3×3 Matrix{Int64}:
 5  0  0
 0  5  0
 0  0  5

Create a rectangular matrix of Boolean type numbers with units on the main diagonal:

In [ ]: 
Matrix{Bool}(I, 4, 10)
Out[0]:
4×10 Matrix{Bool}:
 1  0  0  0  0  0  0  0  0  0
 0  1  0  0  0  0  0  0  0  0
 0  0  1  0  0  0  0  0  0  0
 0  0  0  1  0  0  0  0  0  0

Create, for comparison, a unit matrix of size 5 by 5.

In [ ]: 
eye_I = Matrix{Bool}(I, 5, 5)
Out[0]:
5×5 Matrix{Bool}:
 1  0  0  0  0
 0  1  0  0  0
 0  0  1  0  0
 0  0  0  1  0
 0  0  0  0  1

LinearAlgebra library: command Diagonal

Another interesting way is provided by the function Diagonal. It can be used to position a vector along the main diagonal of a square matrix. Other elements are not stored in memory, the matrix is displayed as a sparse array (SparseArray).

In [ ]: 
using LinearAlgebra
eye_Diagonal = Diagonal( repeat([1], 5) )
Out[0]:
5×5 Diagonal{Int64, Vector{Int64}}:
 1  ⋅  ⋅  ⋅  ⋅
 ⋅  1  ⋅  ⋅  ⋅
 ⋅  ⋅  1  ⋅  ⋅
 ⋅  ⋅  ⋅  1  ⋅
 ⋅  ⋅  ⋅  ⋅  1

LinearAlgebra library: function diagm

Another function allows you to position vectors along any diagonal of a matrix, not just the main diagonal.

In [ ]: 
eye_diagm = diagm( 0 => repeat([1], 5) )
Out[0]:
5×5 Matrix{Int64}:
 1  0  0  0  0
 0  1  0  0  0
 0  0  1  0  0
 0  0  0  1  0
 0  0  0  0  1

Let's visualise this matrix:

In [ ]: 
heatmap( diagm( 0 => repeat([1], 5), 3=>repeat([1], 5) ), yflip=:true, c=:Blues )
Out[0]:

SparseArrays Library

You can construct a unit sparse matrix directly using the library SparseArrays as follows:

In [ ]: 
using SparseArrays: spdiagm
eye_spdiagm = spdiagm(0 => repeat([1], 5))
Out[0]:
5×5 SparseArrays.SparseMatrixCSC{Int64, Int64} with 5 stored entries:
 1  ⋅  ⋅  ⋅  ⋅
 ⋅  1  ⋅  ⋅  ⋅
 ⋅  ⋅  1  ⋅  ⋅
 ⋅  ⋅  ⋅  1  ⋅
 ⋅  ⋅  ⋅  ⋅  1

Such a matrix can also be visualised without any problems.

In [ ]: 
heatmap( eye_spdiagm, yflip=:true, c=:Greens )
Out[0]:

But, as we will see, it occupies more memory than all previous variants.

In [ ]: 
varinfo( r"eye_array|eye_I|eye_Diagonal|eye_diagm|eye_spdiagm" )
Out[0]:
name size summary
eye_Diagonal 88 bytes 5×5 Diagonal{Int64, Vector{Int64}}
eye_I 65 bytes 5×5 Matrix{Bool}
eye_array 240 bytes 5×5 Matrix{Float64}
eye_diagm 240 bytes 5×5 Matrix{Int64}
eye_spdiagm 288 bytes 5×5 SparseArrays.SparseMatrixCSC{Int64, Int64} with 5 stored entries

Conclusion

We have compared several ways of creating unit matrices. The matrix created using the standard object I.