Working with Numerical Arrays

Get started handling multidimensional arrays in Julia with this simple guide.

by Martin D. Maas, Ph.D

@MartinDMaas

Last updated: 2021-10-15

Manipulating Multidimensional Arrays, such as Vectors and Matrices.

Julia, like the best mathematical languages, has native support for vectors, matrices, and multidimensional arrays, which are fundamental objects of computational science and engineering.

In the examples in this post, we will be relying almost exclusively on what we can do with the Julia Base syntax.

Initializing Arrays

Direct Input

The most convenient way to input a matrix is by using whitespace-separated columns, and semicolons for rows, as follows:

A = [1 2 3; 1 2 4; 2 2 2]     # 3×3 Matrix{Float64}

For easier readability, you can also resort to newlines:

A = [1 2 3; 
     1 2 4; 
     2 2 2]     

As for inputting vectors, each element can be separated by either commas or semicolons. However, note that separating entries with whitespace will result in a 1x3 matrix, which is a different type of entity for Julia than a vector.

b1 = [4.0, 5, 6]                # 3-element Vector{Float64}
b2 = [4.0; 5; 6]                # 3-element Vector{Float64}
m1 = [4.0 5 6]                  # 1×3 Matrix{Float64}

Julia also supports arrays of non-numerical types such as Strings, or even arrays of ‘Any’, which could include strings and numbers, and can be initialized as:

A = ["Hello", 1, 2, 3]
  4-element Vector{Any}:
  "Hello"
  1
  2
  3 

Uninitialized (undef) Arrays

Array in Julia can be constructed in various ways, other than direct initialization of arrays from given values.

For examples, for performance reasons, it is wise to initialize arrays of a given type (and size), without specifying any values. That is, not even initializing a new array with zeros.

To do this, we can employ keywords such as Vector{T}, Matrix{T} or Array{T}, where T is a type (see, for example, Numeric Types in the documentation). Among numeric types, arguably Float64 (double-precision real), ComplexF64 (double-precision complex), and Int64 are among the most frequently used in numerical computing.

The reason to declare uninitialized arrays is to be able to fill them later – for example, with the use of a for loop.

What happens under the hood when we declare an undef array, is that a certain portion of memory gets reserved (or allocated) for this specific use. As the computer is not even filling that chunk of memory with zeros, we are saving some time.

The following are some examples of undef arrays of type Float64:

n = 5
A1 = Array{Float64}(undef,n,n)          # 5×5 Matrix{Float64}
A2 = Matrix{Float64}(undef,n,n)         # 5×5 Matrix{Float64}

V1 = Array{Float64}(undef,n)            # 5-element Vector{Float64}
V2 = Vector{Float64}(undef,n)           # 5-element Vector{Float64}

For non-numerical types like String or Any, we can resort to

A = Array{String}(undef,n)
A = Array{Any}(undef,n)

Empty Arrays

Empty arrays can be a useful starting point when we are in a situation where it is hard or impossible to know which array sizes we need in advance. An empty array can be later grown dynamically, and filled with values. (see the “Dynamic Arrays” section below in this page).

To initialize an empty array, it is perfectly valid to use \( n = 0 \) in the above expressions.

v = Array{Float64}(undef,0)

There is a shorthand for this expression:

v = Float64[]

A possible source of errors would be to confuse this array with an empty array of “Any” type, which is initialized as follows:

v = []    # Same as Any[], and you can't change this type easily later (gotcha!)

If we later fill this array dynamically with Float values, its type would remain fixed at “Any”, which could lead to bad performance – or even errors, if we plan on passing this array to a function which requires arrays of type Float, for example.

Generators

Julia provides many functions for constructing and initializing special kinds of arrays. For example, for null matrices we can use the ‘zeros’ function, and for matrices filled with ones we have the ‘ones’ function.

A = ones(8,9)
B = zeros(8,9)

The documentation includes many other examples, such as matrices with different types of random entries, boolean values, and so forth, such as:

C = rand(6,6)

In the case of the identity matrix, placing an integer right before a letter ‘I’ will do the job, but this requires the LinearAlgebra package:

using LinearAlgebra
M = 5I + rand(2,2)

 2×2 Matrix{Float64}:
 5.50162   0.462804
 0.498287  5.30439

Array Operations

Now that we know several ways of inputting arrays, we should take a look at how we can operate with them.

Getting the Size of an Array

The size of an array can be obtained with the ‘size’ function. Importantly, the elements of s=size(A) are treated as an array of integers, one for each dimension, even in the one-dimensional case. That is, we should always employ a syntax of the form s[dim].

A = rand(6)
s = size(A)
print(s)        # (6,)
print(s[1])     # 6

Matrix multiplication

Two dimensional arrays (or matrices) are a fundamental part of Julia. For example, you can resort to a Matlab-style sytax for matrix-matrix multiplication:

A * B

A Matrix and a Vector can be also multiplied with the * operator.

A * v

Element-wise multiplication

In case you need to multiply the elements of an n-dimensional array on an element-wise fashion, you can resort to the dot operator, which will broadcast the scalar multiplication operator on an element-wise fashion, just as we discussed in the previous post of this tutorial series:

A .* B

Dot product

You can do dot products by calling the ‘dot’ function

v = rand(1000)
w = rand(1000)
z = dot(v,w)

Alternatively, you can resort to a typical linear algebra notation:

z = v'w

Array slicing

Array slicing is a very useful operation in Linear Algebra, which consists in extracting a subset of elements from an array and packaging them as another array. In Julia, we can resort to a very straightforward syntax for this operation.

For example, the code below will extract the odd indices in both dimensions from a 6x6 matrix, and store the results in a new 3x3 matrix.

A = rand(6,6)           # 6×6 Matrix{Float64}
B = A[1:2:5,1:2:5]      # 3×3 Matrix{Float64}

Logical Indexing

A very useful way of selecting array elements is to resort to a boolean array, in a process which is often referret to as logical indexing. This can be easily done in Julia by resorting to broadcasting the logical operations of interest.

For example, lets generate a matrix of random numbers and set to zero all numbers in the array which are below 0.5:

A = rand(6,6)
A[ A .< 0.5 ] .= 0

Note that in this example we also needed to broadcast the assignment operator as well.

Backslash operator

Just like in Matlab, Julia has a built-in operator to solve matrices. For square matrices, it will try to solve the linear system, while for rectangular matrices, it will seek for a least squares solution.

One ‘gotcha’ that you will probably encounter sooner or later is that, as 1xN matrices are not the same as N-element vectors. We have to be careful and always employ vectors in the right-hand side of an equation.

b1 = [4.0, 5, 6]                # 3-element Vector{Float64}
b2 = [4.0; 5; 6]                # 3-element Vector{Float64}
m1 = [4.0 5 6]                  # 1×3 Matrix{Float64}

x=A\b1                          # Solves A*x=b
x=A\b2                          # Solves A*x=b  
x=A\m1                          # Error!!

Dynamic Arrays

Arrays in Julia are Dynamic (i.e. growable, resizable) unless defined otherwise. We can change this default behavior by resorting to the StaticArrays.jl package, which can enhance performance for certain type of operations.

But lets go about how to use dynamic arrays in Julia.

Resizing 1D Arrays with Push and Pop

One-dimensional Arrays can be resized with functions such as Push and Pop. In our example below we are using numerical arrays, but these functions can be also applied to non-numerical arrays.

A = Float64[]       # Equivalent to A=Array{Float64}(undef,0)
push!(A, 4)         # Adds the number 4 at the end of the array
push!(A, 3)         # Adds the number 3 at the end of the array
v = pop!(A)         # Returns 3 and removes it from A

There are also functions such as pushfirst!, and popfirst! which work on the beginning of the array as opposed to the end. Additionally, the function splice!(A,i) and deleteat!(A,i) will incorporates (resp. delete) an element at a given position i.

Multi-Dimensional Array Concatenation

A common operation is to concatenate arrays, for example, to add rows to an existing matrix. In Julia, there are a few ways to do this.

For example, we can resort to hcat, vcat, and cat functions, to concatenate in the horizontal, vertical, or along any given dimension.

For example,

A = [4 5 6] 
B = [6 7 8] 

M1 = vcat(A, B)
2×3 Matrix{Int64}:
 4  5  6
 6  7  8

M2 = hcat(A, B)
1×6 Matrix{Int64}:
 4  5  6  6  7  8

The cat function, in turn, can be used to acomplish the same results as above. Importantly, cat will also work for n-dimensional arrays, as we can concatenate along any dimension. For example.

M1 = cat(A, B, dims=1)
M2 = cat(A, B, dims=2)
M3 = cat(A, B, dims=3)

For 2D matrices, we can also resort to the same syntax that we used to input a matrix, using brackets and spaces and semi-colons as alternatives to hcat and vcat, respectively. For example,

M1 = [A; B]
2×3 Matrix{Int64}:
 4  5  6
 6  7  8

M2 = [A B]
1×6 Matrix{Int64}:
 4  5  6  6  7  8

Conclusion

We have covered the basics of Array manipulation, using the Julia Base syntax.

There is, of course, more array functionality in several packages, such as LinearAlgebra.jl, StaticArrays.jl, SparseArrays.jl, which will be the subject of future posts. Stay tuned!