Examples

Examples

Let's show some examples of how use Dolang to compile simple Julia functions.

Basic example

First, let's construct some equations, arguments, and parameters

using Dolang
eqs = [
    :(foo(0) = log(a(0))+b(0)/x(-1)),
    :(bar(0) = c(1)+u*d(1))
]
args = [(:a, -1), (:a, 0), (:b, 0), (:c, 0), (:c, 1), (:d, 1)]
params = [:u]
defs = Dict(:x=>:(a(0)/(1-c(1))))
targets = [(:foo, 0), (:bar, 0)]

ff = FunctionFactory(
    eqs, args, params, targets=targets, defs=defs, funname=:myfun
)
FunctionFactory
            name: myfun
  # of equations: 2
       # of args: 6
     # of params: 1
    Has targets?: true

Now that we have our FunctionFactory we can generate some code (warning, the generated code is not pretty):

code = make_function(ff)
begin
    @generated function myfun{D}(::Dolang.TDer{D}, V::AbstractVector, p)
            ff = Dolang.FunctionFactory{Array{Tuple{Symbol,Int64},1},Array{Symbol,1},Dict{Symbol,Expr},DataType}(Expr[:(_foo__0_ = log(_a__0_) + _b__0_ / (_a_m1_ / (1 - _c__0_))), :(_bar__0_ = _c__1_ + _u_ * _d__1_)], Tuple{Symbol,Int64}[(:a, -1), (:a, 0), (:b, 0), (:c, 0), (:c, 1), (:d, 1)], Symbol[:u], Symbol[:_foo__0_, :_bar__0_], Dict(:x=>:(a(0) / (1 - c(1)))), :myfun, Dolang.SkipArg, Dolang.IncidenceTable(Dict(2=>Dict(:d=>Set([1]),:bar=>Set([0]),:c=>Set([1])),1=>Dict(:a=>Set([0, -1]),:b=>Set([0]),:c=>Set([0]),:foo=>Set([0]),:x=>Set([-1]))), Dict(:a=>Set([0, -1]),:b=>Set([0]),:d=>Set([1]),:bar=>Set([0]),:c=>Set([0, 1]),:foo=>Set([0]),:x=>Set([-1])), Dict(0=>Set(Symbol[:a, :b, :bar, :c, :foo]),-1=>Set(Symbol[:a, :x]),1=>Set(Symbol[:d, :c]))))
            Dolang.func_body(ff, Der{D})
        end
    @generated function myfun!{D}(::Dolang.TDer{D}, out, V::AbstractVector, p)
            ff = Dolang.FunctionFactory{Array{Tuple{Symbol,Int64},1},Array{Symbol,1},Dict{Symbol,Expr},DataType}(Expr[:(_foo__0_ = log(_a__0_) + _b__0_ / (_a_m1_ / (1 - _c__0_))), :(_bar__0_ = _c__1_ + _u_ * _d__1_)], Tuple{Symbol,Int64}[(:a, -1), (:a, 0), (:b, 0), (:c, 0), (:c, 1), (:d, 1)], Symbol[:u], Symbol[:_foo__0_, :_bar__0_], Dict(:x=>:(a(0) / (1 - c(1)))), :myfun, Dolang.SkipArg, Dolang.IncidenceTable(Dict(2=>Dict(:d=>Set([1]),:bar=>Set([0]),:c=>Set([1])),1=>Dict(:a=>Set([0, -1]),:b=>Set([0]),:c=>Set([0]),:foo=>Set([0]),:x=>Set([-1]))), Dict(:a=>Set([0, -1]),:b=>Set([0]),:d=>Set([1]),:bar=>Set([0]),:c=>Set([0, 1]),:foo=>Set([0]),:x=>Set([-1])), Dict(0=>Set(Symbol[:a, :b, :bar, :c, :foo]),-1=>Set(Symbol[:a, :x]),1=>Set(Symbol[:d, :c]))))
            Dolang.func_body!(ff, Der{D})
        end
    @generated function myfun(::Dolang.TDer{0}, V::AbstractVector, p)
            ff = Dolang.FunctionFactory{Array{Tuple{Symbol,Int64},1},Array{Symbol,1},Dict{Symbol,Expr},DataType}(Expr[:(_foo__0_ = log(_a__0_) + _b__0_ / (_a_m1_ / (1 - _c__0_))), :(_bar__0_ = _c__1_ + _u_ * _d__1_)], Tuple{Symbol,Int64}[(:a, -1), (:a, 0), (:b, 0), (:c, 0), (:c, 1), (:d, 1)], Symbol[:u], Symbol[:_foo__0_, :_bar__0_], Dict(:x=>:(a(0) / (1 - c(1)))), :myfun, Dolang.SkipArg, Dolang.IncidenceTable(Dict(2=>Dict(:d=>Set([1]),:bar=>Set([0]),:c=>Set([1])),1=>Dict(:a=>Set([0, -1]),:b=>Set([0]),:c=>Set([0]),:foo=>Set([0]),:x=>Set([-1]))), Dict(:a=>Set([0, -1]),:b=>Set([0]),:d=>Set([1]),:bar=>Set([0]),:c=>Set([0, 1]),:foo=>Set([0]),:x=>Set([-1])), Dict(0=>Set(Symbol[:a, :b, :bar, :c, :foo]),-1=>Set(Symbol[:a, :x]),1=>Set(Symbol[:d, :c]))))
            Dolang.func_body(ff, Der{0})
        end
    @generated function myfun!(::Dolang.TDer{0}, out, V::AbstractVector, p)
            ff = Dolang.FunctionFactory{Array{Tuple{Symbol,Int64},1},Array{Symbol,1},Dict{Symbol,Expr},DataType}(Expr[:(_foo__0_ = log(_a__0_) + _b__0_ / (_a_m1_ / (1 - _c__0_))), :(_bar__0_ = _c__1_ + _u_ * _d__1_)], Tuple{Symbol,Int64}[(:a, -1), (:a, 0), (:b, 0), (:c, 0), (:c, 1), (:d, 1)], Symbol[:u], Symbol[:_foo__0_, :_bar__0_], Dict(:x=>:(a(0) / (1 - c(1)))), :myfun, Dolang.SkipArg, Dolang.IncidenceTable(Dict(2=>Dict(:d=>Set([1]),:bar=>Set([0]),:c=>Set([1])),1=>Dict(:a=>Set([0, -1]),:b=>Set([0]),:c=>Set([0]),:foo=>Set([0]),:x=>Set([-1]))), Dict(:a=>Set([0, -1]),:b=>Set([0]),:d=>Set([1]),:bar=>Set([0]),:c=>Set([0, 1]),:foo=>Set([0]),:x=>Set([-1])), Dict(0=>Set(Symbol[:a, :b, :bar, :c, :foo]),-1=>Set(Symbol[:a, :x]),1=>Set(Symbol[:d, :c]))))
            Dolang.func_body!(ff, Der{0})
        end
    function myfun(V::AbstractVector, p)
        myfun(Dolang.Der{0}, V::AbstractVector, p)
    end
    function myfun!(out, V::AbstractVector, p)
        myfun!(Dolang.Der{0}, out, V::AbstractVector, p)
    end
    begin
        function myfun(::Dolang.TDer{0}, V::AbstractArray, p)
            out = Dolang._allocate_out(eltype(V), 2, V)
            begin
                nrow = size(out, 1)
                for _row = 1:nrow
                    __out__row = view(out, _row, :)
                    __V__row = Dolang._unpack_obs(V, _row)
                    myfun!(Dolang.Der{0}, __out__row, __V__row, p)
                end
                return out
            end
        end
        function myfun(V::AbstractArray, p)
            out = Dolang._allocate_out(eltype(V), 2, V)
            begin
                nrow = size(out, 1)
                for _row = 1:nrow
                    __out__row = view(out, _row, :)
                    __V__row = Dolang._unpack_obs(V, _row)
                    myfun!(Dolang.Der{0}, __out__row, __V__row, p)
                end
                return out
            end
        end
    end
    begin
        function myfun!(::Dolang.TDer{0}, out, V::AbstractArray, p)
            begin
                expected_size = Dolang._output_size(2, V)
                if size(out) != expected_size
                    msg = "Expected out to be size $(expected_size), found $(size(out))"
                    throw(DimensionMismatch(msg))
                end
            end
            begin
                nrow = size(out, 1)
                for _row = 1:nrow
                    __out__row = view(out, _row, :)
                    __V__row = Dolang._unpack_obs(V, _row)
                    myfun!(Dolang.Der{0}, __out__row, __V__row, p)
                end
                return out
            end
        end
        function myfun!(out, V::AbstractArray, p)
            begin
                expected_size = Dolang._output_size(2, V)
                if size(out) != expected_size
                    msg = "Expected out to be size $(expected_size), found $(size(out))"
                    throw(DimensionMismatch(msg))
                end
            end
            begin
                nrow = size(out, 1)
                for _row = 1:nrow
                    __out__row = view(out, _row, :)
                    __V__row = Dolang._unpack_obs(V, _row)
                    myfun!(Dolang.Der{0}, __out__row, __V__row, p)
                end
                return out
            end
        end
    end
    (myfun, myfun!)
end

In order to actually call methods of myfun we need to call eval on our code:

eval(code)
(myfun, myfun!)

And now we can call it

V = [1.5, 2.5, 1.0, 2.0, 3.0, 1.1]
p = [0.5]
myfun(V, p)

We can call the vectorized version:

myfun([V V]', p)

Or the mutating one:

out = zeros(2)
myfun!(out, V, p)
out

We can evaluate the Jacobian (first derivative matrix)

myfun(Der{1}, V, p)

... or the Hessian

myfun(Der{2}, V, p)

... or higher order derivatives

myfun(Der{3}, V, p)
myfun(Der{10}, V, p)

Note that the Hessian is returned in as a SparseMatrixCSC where the columns are ordered with the first variable increasing faster. In this example, the (1, 4) element of the Hessian will be the second partial derivative of the first equation (the 1) with fourth and first variables in args.

For all higher order derivatives the return value is a Vector{Dict{NTuple{N,Int}, Float64}}, where N is the order of derivative. The keys of each dict will be non-increasing so mixed partials are only computed and stored once. Notice that for the order 3 derivatives an entry for (1, 1, 3) appears, but not for (1, 3, 1) or (3, 1, 1). The user of these routines is required to handle the equality of partials.

With dispatch

Now lets consider an example that leverages the dispatch argument. In this case we will use the convenience method for make function with signature: make_function(::Vector{Expr}, ::AbstractVector, ::AbstractVector). We also will not show the compiled code

eqs = [
    :(sin(x(0)) + exp(2*x(1))),
    :(y(0) / (2 * (1 - β))),
]
# NOTE: Arguments can be pre-normalized and contain unicode
variables = [(:x, 0), (:y, 0), :_x__1_, :β]
to_diff = 1:3
code = make_function(eqs, variables, to_diff, dispatch=Int, name=:my_int_fun)
eval(code)
(my_int_fun, my_int_fun!)

Let's check the methods of my_int_fun and my_int_fun!:

methods(my_int_fun)
methods(my_int_fun!)

Notice that because we set the dispatch argument to Int, all methods of either the mutating or allocating function require that an Int is passed to direct dispatch.

Let's try calling our function

V = [π, 5.0, 0.8]
p = [0.98]
my_int_fun(V, p)  # doesn't work

The above fails with

ERROR: MethodError: no method matching my_int_fun(::Array{Float64,1}, ::Array{Float64,1})

Now if we call that method where the first argument is an instance of Int, it will work

V = [π, 5.0, 0.8]  # hide
p = [0.98]  # hide
my_int_fun(1, V, p)

The actual integer we pass doesn't matter

my_int_fun(10_000_001, V, p)

We can still evaluate derivatives

my_int_fun(Der{1}, 1, V, p)

Notice that the order of derivative comes first, then the dispatch argument, then all other function arguments.

Higher order derivatives also work

my_int_fun(Der{2}, 10, V, p)
my_int_fun(Der{5}, 10, V, p)
my_int_fun(Der{10}, 10, V, p)

With grouped args

Finally we show one more example of how to use an associative mapping for the FunctionFactory args field to produce a function with grouped arguments.

eqs = [
    :(sin(x(0)) + exp(2*x(1))),
    :(y(0) / (2 * (1 - β))),
]
args = Dict(
    :a => [(:x, 0), (:y, 0)],
    :b => [(:x, 1)]
)
params = [:β]
ff = FunctionFactory(eqs, args, params, funname=:grouped_args)
code = make_function(ff)
eval(code)
(grouped_args, grouped_args!)

Let's take a look at the methods for this function

methods(grouped_args)

Notice that each of theses routines have arguments for a, b, and p instead of just the V and p we saw in previous examples.

When we evaluate these methods we need to be sure that a has two elements, b has one, and p has one:

a = [π, 5.0]
b = [0.8]
p = [0.98]
grouped_args(a, b, p)

Notice we can evaluate a method where a and b are vectorized. Here p will be repeated as necessary

grouped_args([a a+1]', [b b+1]', p)

We can also evaluate a partially vectorized version of the function where only a is a matrix and b is a vector. In this case both b and p will be repeated

grouped_args([a a+1]', b, p)

Non-allocating versions of the function are also available

out = zeros(2)
grouped_args!(out, a, b, p)
@show @allocated grouped_args!(out, a, b, p)
@show out

As of today (2017-06-13) only first order derivative code has been implemented for functions with grouped args. Trying to evaluate a derivative will result in an error that looks like this:

grouped_args(Der{1}, a, b, p)

Trying to evaluate any higher order derivative will result in an error that looks like this:

julia> grouped_args(Der{2}, a, b, p)
ERROR: MethodError: no method matching equation_block(::Dolang.FunctionFactory{Dict{Symbol,Array{Tuple{Symbol,Int64},1}},Array{Symbol,1},Dict{Symbol,Any},DataType}, ::Type{Dolang.Der{2}})