Functions

Syntax

To improve the structure and ensure reusability of pieces of our program we can use the function command. We have already used functions when, for example, calling typeof(input). Other examples of functions that can be found on any common calculators are sin(x) or exp(x). The syntax to define our own functions is the following:

function foo(input)
    # function body that is executed when foo is called

    return output
end

So if we want to define a function which prints out and returns $\sin(\cos(x))$ we can write

function sincos(x)
    result = sin(cos(x))
    println("sin(cos($x)) = $result")
    
    return result
end

If we wish to specify multiple inputs and outputs we can do so as well:

function sincos(x, y)
    result1 = sin(cos(x))
    result2 = sin(cos(y))
    
    println("sin(cos($x)) = $result1")
    println("sin(cos($y)) = $result2")
    
    return result1, result2
end

We can call a function with multiple outputs via out = foo(input) and access the output at index $i$ through out[i]. We can also write (assuming two outputs) out1, out2 = foo(input) and out1, _ = foo(input) if we only need one of the outputs.

julia> x = 1; y = 1.5;
julia> res1, res2 = sincos(x, y)
sin(cos(1)) = 0.5143952585235492
sin(cos(1.5)) = 0.07067822452613834
(0.5143952585235492, 0.07067822452613834)

julia> res1
0.5143952585235492

julia> res = sincos(x, y);

julia> res[1]
0.5143952585235492

Let us practice this syntax by revisiting loops:

function mysum1(V)
    s = zero(eltype(V))

    for i in eachindex(V)
        s += V[i]
    end

    return s
end

function mysum2(V)
    s = zero(eltype(V))

    for v in V
        s += v
    end

    return s
end

V = rand(100_000)
mysum1(V) ≈ mysum2(V)
isapprox(mysum1(V),  mysum2(V); atol=1e-10, rtol=1e-10)

Anonymous functions

Functions can also be created anonymously, that is, without giving a name. We call these functions anonymous functions and they are especially feasible when we want to use a function as an argument. The following snippet creates an anonymous function that takes one argument x and returns $x^2 + 1$:

julia> x -> x^2 + 1
#1 (generic function with 1 method)

Unfortunately, this function can not be accessed again, since we do not have a name/variable which could be accessed. But we can apply a value right away and e.g. evaluate this anonymous function for $x = 3$:

julia> (x -> x^2 + 1)(3)
10

This is not a good application for an anonymous function but there are applications where they are indeed quite helpful. Throughout this workshop, we will occasionally use anonymous functions as function arguments. Let us take a look at the function map(f, c) that allows us to transform a collection (e.g. a vector) c by applying function f to every element. The following example applies $x \mapsto x^2 + 1$ to every element of the vector [1, 2, 3, 4, 5]:

julia> map(x -> x^2 + 1, [1, 2, 3, 4, 5])
5-element Vector{Int64}:
  2
  5
 10
 17
 26

Element-wise operations and input specifications

As always, we can use the dot operation . to evaluate an array of inputs element-wise. Define the scalar function

function sincos(x::Float64)
    return sin(cos(x))    
end

and run

julia> x = ones(2, 3);

julia> sincos.(x)
2×3 Matrix{Float64}:
 0.514395  0.514395  0.514395
 0.514395  0.514395  0.514395

Moreover, we can assign values to inputs in the function definition. If the caller does not specify the input, these values will be used instead.

function sincos(x::Float64=(0.5 * pi))
    return sin.(cos.(x))
end

We can now call this function via

julia> sincos()
6.123233995736766e-17

julia> sincos(0.0)
0.8414709848078965

Note that the dot operation can also be used for built-in functions like sin or cos.

Call by reference

Julia functions do not copy the input but directly operate on the input data. This means that changing values of the input in the function body will also change this data for the function caller. Whenever we define a function which will modify the input, we should indicate this with a ! behind the function name:

function sincos!(x)
    x .= sin.(cos.(x))

    return x
end

Calling this function leads to

julia> x = ones(2);
julia> println("Function value is ", sincos!(x))
Function value is [0.5143952585235492, 0.5143952585235492]

julia> x
2-element Vector{Float64}:
 0.5143952585235492
 0.5143952585235492

function sincos1!(x)
    x .= sin.(cos.(x))

    return x
end

function sincos2!(x)
    x = sin.(cos.(x))

    return x
end

julia> x = ones(2);

julia> sincos1!(x);

julia> x
2-element Vector{Float64}:
 0.5143952585235492
 0.5143952585235492

julia> x = ones(2);

julia> sincos2!(x);

julia> x
2-element Vector{Float64}:
 1.0
 1.0

function sincos1!(x)
    println("Address input: ", pointer_from_objref(x))
    x .= sin.(cos.(x))
    println("Address output: ", pointer_from_objref(x))

    return x
end

function sincos2(x)
    println("Address input: ", pointer_from_objref(x))
    x = sin.(cos.(x))
    println("Address output: ", pointer_from_objref(x))

    return x
end

julia> x = ones(2);

julia> y = sincos1!(x);
Address input: Ptr{Nothing} @0x00007f3325795b90
Address output: Ptr{Nothing} @0x00007f3325795b90

julia> x = ones(2);

julia> println("Address caller: ", pointer_from_objref(x));
Address caller: Ptr{Nothing} @0x00007f332671cab0

julia> y = sincos1!(x);
Address input: Ptr{Nothing} @0x00007f332671cab0
Address output: Ptr{Nothing} @0x00007f332671cab0

julia> y = sincos2(x);
Address input: Ptr{Nothing} @0x00007f332671cab0
Address output: Ptr{Nothing} @0x00007f33266e2310

Multiple dispatch

One might have observed that since we did not specify any data types, we were able to call functions using vectors and scalars. However, if we call sincos1!(1.0) we see that this might not always be the best idea. Some functions should only be called with a certain data type. We can specify the data type of input and output in the following way:

function sincos1!(x::Array{Float64, 1})
    x .= sin.(cos.(x))

    return x
end

In the same way, we can define functions that have the same name but which perform different operations depending on the data type. For example, we can define the function sincos(x) which evaluates $\sin(\cos(x))$ when $x$ is a matrix or vector and returns a vector. In order to rearrange a matrix $M$ to a vector $m$ we can use m = vec(M). Then, we get:

function sincos(x::Array{Float64, 1})
    println("My input is a vector.")

    return sin.(cos.(x))
end

function sincos(x::Array{Float64, 2})
    println("My input is a matrix.")

    return vec(sin.(cos.(x)))
end

Calling these functions gives

julia> x = ones(2, 2);

julia> sincos(x)
My input is a matrix.
4-element Vector{Float64}:
 0.5143952585235492
 0.5143952585235492
 0.5143952585235492
 0.5143952585235492

julia> x = ones(4);

julia> sincos(x)
My input is a vector.
4-element Vector{Float64}:
 0.5143952585235492
 0.5143952585235492
 0.5143952585235492
 0.5143952585235492

Parametric types for functions

Notice that it is very restrictive to tell our method to only accept inputs of type Float64. It makes perfect sense to evaluate our function at $x = 1$, where $x$ can be an integer. In fact, it is considered to be good practice if we make function inputs as general as possible. Just as for structs, we can use parametric types to make the input more general. If we want the input to be a real number (i.e., a subtype of Real), we can write this as

function sincos(x::T) where T<:Real
    return sin(cos(x))
end

Now, T can be any subtype of Real, that is, we can write

julia> sincos(1)
0.5143952585235492

julia> sincos(1.0)
0.5143952585235492

julia> sincos(1//2)
0.7691963548410085

struct Point{T}
    x::T
    y::T
end

struct Point{T<:Real}
    x::T
    y::T
end

function sincos(x::T) where T<:Real
    return sin(cos(x))
end

function sincos(x::Point{T}) where T<:Real
    norm = sqrt(x.x^2 + x.y^2)
    
    return sin(cos(norm))
end

julia> p = Point{Float64}(1.0, 2.0)
Point{Float64}(1.0, 2.0)

julia> sincos(p)
-0.578813455551511

julia> pim = Point{Complex}(1im, 2.0)
Point{Complex}(0 + 1im, 2.0 + 0.0im)

julia> sincos(pim)
ERROR: MethodError: no method matching sincos(::Point{Complex})
You may have intended to import Base.sincos
Closest candidates are:
  sincos(::T) where T<:Real at REPL[3]:1
  sincos(::Point{T}) where T<:Real at REPL[4]:1
Stacktrace:
 [1] top-level scope
   @ REPL[9]:1

Constructors

Note from the previous exercise, that it might be convenient if every object of type Point computes and stores the norm. Of course this can be done by defining

struct PointFull{T<:Real}
    x::T
    y::T
    norm
end

and creating objects of type PointFull via

julia> p = PointFull{Float64}(1.0, 2.0, sqrt(1.0^2 + 2.0^2))
PointFull{Float64}(1.0, 2.0, 2.23606797749979)

julia> p = PointFull{Float32}(1.0, 2.0, sqrt(1.0^2 + 2.0^2))
PointFull{Float32}(1.0f0, 2.0f0, 2.23606797749979)

julia> p = PointFull{Int64}(1, 2, sqrt(1 + 2))
PointFull{Int64}(1, 2, 1.7320508075688772)

Note that this is quite tedious since we need to copy paste the same formula for the norm every time we construct on object. Conveniently, Julia provides us with constructors, which are functions that are called whenever we create on object of our struct. The syntax is the following:

struct PointFull1{T<:Real}
    x::T
    y::T
    norm

    function PointFull1(x::T, y::T) where {T<:Real}
        norm = sqrt(x^2 + y^2)
        new{T}(x, y, norm)
    end
end

Now, we can call

julia> PointFull1(1.0, 2.0)
PointFull1{Float64}(1.0, 2.0, 2.23606797749979)

julia> PointFull1(1, 2)
PointFull1{Int64}(1, 2, 2.23606797749979)

julia> PointFull1(1im, 2im)
ERROR: MethodError: no method matching PointFull1(::Complex{Int64}, ::Complex{Int64})
Stacktrace:
 [1] top-level scope
   @ REPL[20]:1