Julia: functions and structures-as-functions

Despite the fact that the concept of Julia lacks “classical” object-oriented programming with classes and methods, the language provides abstraction tools, a key role in which is played by the type system and elements of functional programming. Let us consider the second point in more detail.

The concept of functions in Julia is probably the most similar to languages ​​from the Lisp family (to be more precise, the Lisp-1 branches), and functions can be considered at three levels: as subprograms, as abstractions to a certain sequence of actions, and as data representing this abstraction .

Level 1. Functions as routines

The allocation of subprograms and the assignment of their own names has been going on since prehistoric times, when Fortran was considered a high-level language, and C was not there yet.

In this sense, Julia products are standard. “Feature” can be called the fact that syntactically there is no division into procedures and functions. Regardless of whether the subroutine is called to get some value or just to perform some action on the data, it is called a function.

The definition of a function begins with a keyword function, followed by a list of arguments, a sequence of commands in brackets, and the word ends with a definition end:

"""
    sum_all(collection)
Sum all elements of a collection and return the result
"""
function sum_all(collection)
    sum = 0
    for item in collection
        sum += collection
    end
    sum
end

Distinctive syntax is the behavior inherited from Lisp: for a “normal” return of a value from a function, a word is returnnot necessary: ​​the value of the last one calculated before the endexpression is returned . In the example above, the value of the variable will be returned sum. Thus, it returncan be used as a marker of the special behavior of a function:

function safe_division(number, divisor)
    if divisor == 0
        return 0
    end
    number / divisor
end
# то же самое
function safe_division1(number, divisor)
    if divisor == 0
        0 # при попадании в эту ветку до выхода из функции нет других инструкций
    else
        number / divisor
    end
end

For functions with a short definition, there is a shortened syntax similar to a mathematical notation. So, the calculation of the length of the hypotenuse along the length of the legs can be defined as follows:

hypotenuse(a, b) = sqrt(a^2 + b^2)

The “safe” division using the ternary operator can be written as follows:

safe_division(number, divisor) = divisor == 0 ? 0 : number / divisor

As you can see, it is not necessary to specify types for function arguments. Given how the Julia JIT compiler works, duck typing will not always result in poor performance.

As I tried to demonstrate in a previous article , the Julia compiler can infer the type of the return result by the types of input arguments. Therefore, for example, a function safe_divisionfor quick work requires minimal modification:

function safe_division(number, divisor)
    if divisor == 0
        return zero(number / divisor)
    end
    number / divisor
end

Now, if the types of both arguments are known at the compilation stage, the type of the returned result is also unambiguously displayed, since the function zero(x)returns a null value of the same type as its argument (and division by zero, according to IEEE 754 , has a quite representable value in the format of floating-point numbers).

Functions can have a fixed number of positional arguments, positional arguments with default values, named arguments, and a variable number of arguments. Syntax:

# фиксированный список аргументов
function hello(name)
    println("Hello, ", name)
end
# аргумент со значением по умолчанию
# аргументы со значением по умолчанию должны идти после обязательных аргументов
function greeting_d(name, greeting = "Hello")
    println(greeting, ", ", name)
end
# именованный аргумент со значением по умолчанию
# именованные аргументы отделяются точкой с запятой и идут после
# обязательных аргументов и опциональных позиционных аргументов
function greeting_kw(name; greeting = "Hello")
    println(greeting, ", ", name)
end
# аргумент greeting хоть и именованный, но должен обязательно быть предоставлен при вызове
function greeting_oblkw(name; greeting)
    println(greeting, ", ", name)
end
# функция с неопределённым числом аргументов
# все аргументы, начиная со второго, в теле интерпретируются как кортеж names
function greeting_nary(greeting, names...)
    print(greeting)
    for name in names
        print(", ", name)
    end
    print('\n')
end
julia> hello("world")
Hello, world
julia> greeting_d("world")
Hello, world
julia> greeting_d("Mr. Smith", "How do you do")
How do you do, Mr. Smith
julia> greeting_kw("Mr. Smith")
Hello, Mr. Smith
julia> greeting_kw("mom", greeting = "Hi")
Hi, mom
julia> greeting_oblkw("world")
ERROR: UndefKeywordError: keyword argument greeting not assigned
Stacktrace:
 [1] greeting_oblkw(::String) at ./REPL[23]:3
 [2] top-level scope at none:0
julia> greeting_oblkw("mom", greeting = "Hi")
Hi, mom
julia> greeting_nary("Hi", "mom", "dad", "everyone")
Hi, mom, dad, everyone

Level 2. Functions as data

The function name can be used not only in direct calls, but also as an identifier with which the procedure for obtaining the value is associated. For example:

function f_x_x(fn, x)
    fn(x, x)
end
julia> f_x_x(+, 3)
6 # +(3, 3) = 3+3 = 6
julia> f_x_x(*, 3)
9 # *(3, 3) = 9
julia> f_x_x(^, 3)
27 # ^(3, 3) = 3^3 = 27
julia> f_x_x(log, 3)
1.0 # log(3, 3) = 1

The "classic" functions that take a functional argument are map, reduceand filter.

map(f, x...)applies the function fto the values ​​of all elements from x(or tuples of i-elements) and returns the results as a new collection:

julia> map(cos, [0, π/3, π/2, 2*π/3, π])
5-element Array{Float64,1}:
  1.0                  
  0.5000000000000001   
  6.123233995736766e-17
 -0.4999999999999998   
 -1.0
julia> map(+, (2, 3), (1, 1))
(3, 4)

reduce(f, x; init_val)"Reduces" the collection to a single value, "expanding" the chain f(f(...f(f(init_val, x[1]), x[2])...), x[end]):

function myreduce(fn, values, init_val)
    accum = init_val
    for x in values
        accum = fn(accum, x)
    end
    accum
end

Since it is not really determined in what order the array will pass during the reduction, or whether it will be called fn(accum, x)or fn(x, accum), the reduction will give a predictable result only with commutative or associative operators, such as addition or multiplication.

filter(predicate, x)returns an array of elements xthat satisfy the predicate predicate:

julia> filter(isodd, 1:10)
5-element Array{Int64,1}:
 1
 3
 5
 7
 9
julia> filter(iszero, [[0], 1, 0.0, 1:-1, 0im])
4-element Array{Any,1}:
    [0] 
   0.0  
    1:0 
 0 + 0im

Using higher-order functions for operations on arrays instead of writing a loop has several advantages:

1. the code is getting shorter
2. map()or reduce()show the semantics of the operation being performed, then you still need to understand the semantics of what is happening in the cycle
3. map() allows the compiler to understand that operations on array elements are independent of data, which allows you to apply additional optimizations

Level 3. Functions as Abstractions

Often in map()or filter()you need to use a function that has not been assigned its own name. Julia in this case allows you to express the abstraction of the operations on the argument, without entering your own name for this sequence. Such an abstraction is called an anonymous function , or a lambda function (since in the mathematical tradition such functions are denoted by the letter lambda). The syntax for this view is:

# именованная функция
square(x) = x^2
# анонимный аналог
x -> x^2
# именованная функция
hypot(a, b) = sqrt(x^2 + y^2)
# анонимный аналог - если аргументов больше одного, вокруг них ставятся скобки,
# чтобы не спутать с кортежем из переменной и анонимной функции от одной переменной
(x, y) -> sqrt(x^2 + y^2)
# именованная функция
fortytwo() = 42
# анонимная функция
() -> 42
julia> map(i -> map(x -> x^i, 1:5), 1:5)
5-element Array{Array{Int64,1},1}:
 [1, 2, 3, 4, 5]         
 [1, 4, 9, 16, 25]       
 [1, 8, 27, 64, 125]     
 [1, 16, 81, 256, 625]   
 [1, 32, 243, 1024, 3125]

Both named and anonymous functions can be assigned to variables and returned as values:

julia> double_squared = x -> (2 * x)^2
#17 (generic function with 1 method)
julia> double_squared(5)
100

Variable scope and lexical closures

Normally, they try to write functions in such a way that all the data necessary for the calculation are obtained through formal arguments, i.e. any variable names that occur in the body are either the names of formal arguments or the names of the variables entered inside the body of the function.

function normal(x, y)
    z = x + y
    x + y * z
end
function strange(x, y)
    x + y * z
end

You normal()can say about a function that in its body all variable names are related , i.e. if we everywhere (including the argument list) replace “x” with “m” (or any other identifier), “y” with “n”, and “z” with “sum_of_m_and_n”, the meaning of the expression will not change. In the function, the strange()name z is unbound , i.e. a) the meaning may change if this name is replaced by another; and b) the correctness of the function depends on whether a variable with the name “z” was defined at the time the function was called.

Generally speaking, normal()everything is not so clean with a function either:

1. What happens if a variable named z is defined outside the function?
2. The + and * characters, in fact, are also unrelated identifiers.

With point 2, nothing can be done but to agree - it is logical that the definitions of all the functions used in the system must exist, and we hope that their real meaning corresponds to our expectations.

Point 1 is less obvious than it seems. The fact is that the answer depends on where the function is defined. If it is defined globally, then zinside it normal()will be a local variable, i.e. even if there is a global variable, zits value will not be overwritten. If the definition of the function is inside the code block, then if there is an earlier definition in this block, it zwill be the value of the external variable that will be changed.

If the function body contains the name of an external variable, then this name is associated with the value that existed in the environment where the function was created. If the function itself is exported from this environment (for example, if it is returned from another function as a value), then it “captures” the variable from the internal environment, which is no longer accessible in the new environment. This is called lexical closure.

Closures are mainly useful in two situations: when you need to create a function according to the given parameters and when you need a function that has some internal state.

Consider the situation with a function that encapsulates an internal state:

function f_with_counter(fn)
    call_count = 0
    ncalls() = call_count
    # invoke() будет подсчитывать, сколько раз она вызвана
    # получить это значение можно, вызвав ncalls()
    function invoke(args...)
        call_count += 1
        fn(args...)
    end
    # возвращаем обе функции как кортеж с именованными полями
    # call_count теперь будет извне недоступна напрямую,
    # но экспортированные invoke() и call_count() имеют к нему доступ и могут менять
    (call = invoke, call_count = ncalls)
end
julia> abscount = f_with_counter(abs)
(call = getfield(Main, Symbol("#invoke#22")){typeof(abs)}(abs, Core.Box(0)), call_count = getfield(Main, Symbol("#ncalls#21"))(Core.Box(0)))
julia> abscount.call_count()
0
julia> abscount.call(-20)
20
julia> abscount.call_count()
1
julia> abscount.call(im)
1.0
julia> abscount.call_count()
2

Case study: all the same polynomials

In a previous article, the presentation of polynomials as structures is considered. In particular, one of the storage structures is a list of coefficients, starting with the youngest. To calculate a polynomial pat a point, it was xproposed to call a function evpoly(p, x)that calculates a polynomial according to Horner's scheme.


Representation of a polynomial in the form of a structure does not fully correspond to its intuitive understanding as a mathematical function. But with the help of returning the functional value, polynomials can also be specified directly as functions. So it was:

struct Polynomial
    degree::Int
    coeff::NTuple{N, Float64} where N
    function Polynomial(v::T where T<:Union{Vector{<:Real},
                                            NTuple{<:Any, <:Real}})
        # в случае пустого массива / кортежа в аргументе возвращаем P(x) ≡ 0
        coeff = isempty(v) ? (0.0,) : tuple([Float64(x) for x in v]...)
        # возврат значения - специальным оператором new
        # аргументы - перечисление значений полей
        return new(length(coeff)-1, coeff)
    end
end
"""
    evpoly(p::Polynomial, z::Real)
Evaluate polynomial `p` at `z` using the Horner's rule
"""
function evpoly(p::Polynomial, z::Real)
    ans = p.coeff[end]
    for idx = p.degree:-1:1
        ans = p.coeff[idx] + z * ans
    end
    return ans
end

We transform this definition into a function that takes an array / tuple of coefficients and returns the actual function that computes the polynomial:

function Polynomial_as_closure(v::T where T<:Union{Vector{<:Real},
                                        NTuple{<:Any, <:Real}})
    # в случае пустого массива / кортежа в аргументе возвращаем P(x) ≡ 0
    if isempty(v)
        return x::Real -> 0.0
    end
    coeff = tuple(map(float, v)...)
    degree = length(coeff) - 1
    function evpoly(z::Real)
        ans = coeff[end]
        for idx = degree:-1:1
            ans = coeff[idx] + z * ans
        end
        return ans
    end
    evpoly
end
julia> p = Polynomial_as_closure((0, 1, 1)) # x² + x
(::getfield(Main, Symbol("#evpoly#28")){Tuple{Float64,Float64,Float64},Int64}) (generic function with 1 method)
julia> p(1) # ура, можно обойтись без evpoly()!
2.0
julia> p(11)
132.0

Similarly, you can write a function for the interpolation polynomial.

An important question: was there something that was lost in the new definition in the previous definition? Unfortunately, yes - defining a polynomial as a structure gave hints for the compiler, and for us it was possible to overload arithmetic operators for this structure. Unfortunately, Julia does not provide for functions of such a powerful type system.

Fortunately, in this case we can take the best of both worlds, since Julia allows you to create so-called callable structs. Those. You can specify a polynomial as a structure, but be able to call it as a function! To the definitions of structures from the previous article, you just need to add:

function (p::Polynomial)(z::Real)
    evpoly(p, z)
end
function (p::InterpPolynomial)(z::Real)
    evpoly(p, z)
end

Using functional arguments, you can also add an external constructor of an interpolation polynomial for a certain function constructed from a set of points:

function InterpPolynomial(fn, xsample::T) where {T<:Union{Vector{<:Real},
                                                          NTuple{<:Any, <:Real}}}
    InterpPolynomial(xsample, map(fn, xsample))
end

Conclusion

Opportunities borrowed from functional programming in Julia give a more expressive language compared to a purely imperative style. Representation of structures in the form of functions is a way of more convenient and natural recording of mathematical concepts.