在 Julia 中使用元编程优化递归函数 [英] Optimizing a recursive function with metaprogramming in Julia

问题描述

按照 this answer 的方法，我试图了解究竟发生了什么以及表达式和生成的函数在 Julia 中是如何工作的元编程的概念内.

Following the approach of this answer I am trying to understand what happens exactly and how expressions and generated functions work in Julia within the concept of metaprogramming.

目标是使用表达式和生成函数优化递归函数(具体示例，您可以查看上面提供的链接中回答的问题).

The goal is to optimize a recursive function using expressions and generated functions (for a concrete example you can have a look at the question answered in the link provided above).

考虑以下修改后的斐波那契函数，其中我想计算斐波那契数列直到 n 并将其乘以一个数字 p.

Consider the following modified fibonacci function, in which I want to compute the fibonacci series up to n and multiply it by a number p.

简单的递归实现是

function fib(n::Integer, p::Real)
    if n <= 1
        return 1 * p
    else
        return n * fib(n-1, p)
    end
end

作为第一步，我可以定义一个返回表达式而不是计算值的函数

As a first step, I could define a function which returns an expression instead of the computed value

function fib_expr(n::Integer, p::Symbol)
    if n <= 1
        return :(1 * $p)
    else
        return :($n * $(fib_expr(n-1, p)))
    end
end

其中，例如返回类似

julia> ex = fib_expr(3, :myp)
:(3 * (2 * (1myp)))

通过这种方式，我得到了一个完全扩展的表达式，它取决于分配给符号myp的值.这样我就再也看不到递归了，基本上我是元编程:我创建了一个函数来创建另一个函数"(在这种情况下，我们称之为表达式).我现在可以设置 myp = 0.5 并调用 eval(ex) 来计算结果.但是，这比第一种方法慢.

In this way I get an expression which is fully expanded and depends on the value assigned to the symbol myp. In this way I do not see the recursion anymore, basically I am metaprogramming: I created a function that creates another "function" (in this case we call it expression though). I can now set myp = 0.5 and call eval(ex) to compute the result. However, this is slower than the first approach.

不过，我能做的是通过以下方式生成参数函数

What I can do though, is to generate a parametric function in the following way

@generated function fib_gen{n}(::Type{Val{n}}, p::Real)
    return fib_expr(n, :p)
end

神奇的是，调用 fib_gen(Val{3}, 0.5) 可以完成任务，而且速度非常快.
那么，发生了什么事?

And magically, calling fib_gen(Val{3}, 0.5) gets things done, and is incredibly fast.
So, what is going on?

据我了解，在第一次调用 fib_gen(Val{3}, 0.5) 时，参数函数 fib_gen{Val{3}}(...) 被编译，其内容是通过fib_expr(3, :p)得到的完全展开的表达式，即3*2*1*p和p 替换为输入值.之所以这么快，是因为 fib_gen 基本上只是一系列乘法，而原始 fib 必须在堆栈上分配每个递归调用慢一点，我说的对吗?

To my understanding, in the first call to fib_gen(Val{3}, 0.5), the parametric function fib_gen{Val{3}}(...) gets compiled and its content is the fully expanded expression obtained through fib_expr(3, :p), i.e. 3*2*1*p with p substituted with the input value. The reason why it is so fast then, is because fib_gen is basically just a series of multiplications, whereas the original fib has to allocate on the stack every single recursive call making it slower, am I correct?

为了给出一些数字，这是我的简短基准使用 BenchmarkTools.

To give some numbers, here is my short benchmark using BenchmarkTools.

julia> @benchmark fib(10, 0.5)
...
mean time: 26.373 ns
...

julia> p = 0.5
0.5

julia> @benchmark eval(fib_expr(10, :p))
...
mean time: 177.906 μs
...

julia> @benchmark fib_gen(Val{10}, 0.5)
...
mean time: 2.046 ns
...

我有很多问题:

• 为什么第二种情况这么慢?
• ::Type{Val{n}} 究竟是什么意思?(我从上面链接的答案中复制了它)
• 由于 JIT 编译器，有时我会迷失在 compile-time 和 run-time 发生的事情中，就像这里的情况一样......

• Why the second case is so slow?
• What exactly is and means ::Type{Val{n}}? (I copied that from the answer linked above)
• Because of the JIT compiler, sometimes I am lost in what happens at compile-time and at run-time, as it is the case here...

此外，我尝试根据

@generated function fib_tot{n}(::Type{Val{n}}, p::Real)
    if n <= 1
        return :(1 * p)
    else
        return :(n * fib_tot(Val{n-1}, p))
    end
end

但是速度很慢

julia> @benchmark fib_tot(Val{10}, 0.5)
...
mean time: 4.601 μs
...

我在这里做错了什么?是否可以将 fib_expr 和 fib_gen 组合在一个函数中?

What am I doing wrong here? Is it even possible to combine fib_expr and fib_gen in a single function?

我意识到这更像是一本专着而不是一个问题，尽管我阅读了 元编程部分几次，我很难掌握所有内容，尤其是像这个这样的应用示例.

I realize this is more a monograph rather than a question, however, even though I read the metaprogramming section few times, I am having a hard time to grasp everything, in particular with an applied example such as this one.

推荐答案

专着回应:

首先从普通"宏开始会更容易.我会放宽你使用的定义:

It will be easier to start with "normal" macros first. I'll relax the definition you used a bit:

function fib_expr(n::Integer, p)
    if n <= 1
        return :(1 * $p)
    else
        return :($n * $(fib_expr(n-1, p)))
    end
end

这不仅允许为 p 传递符号，例如整数文字或整个表达式.鉴于此，我们可以为相同的功能定义一个宏:

That allows to pass in more than just symbols for p, like integer literals or whole expressions. Given this, we can define a macro for the same functionality:

macro fib_macro(n::Integer, p)
    fib_expr(n, p)
end

现在，如果在代码中的任何位置使用 @fib_macro 45 1，在编译时它将首先被一个长嵌套表达式替换

Now, if @fib_macro 45 1 is used anywhere in the code, at compile time it will first be replaced by a long nested expression

:(45 * (44 * ... * (1 * 1)) ... )

然后正常编译——变成一个常量.

and then compiled normally -- to a constant.

这就是宏的全部内容，真的.在编译期间替换语法；并且通过递归，这可以是在编译​​和计算表达式函数之间任意长的变化.而对于本质上是不变的，但写起来很乏味的东西，它非常有用:一个很好的例子是 Base.Math.@evalpoly.

That's all there is to macros, really. Replacing syntax during compile time; and by recursion, this can be an arbitrarily long alteration between compiling, and evaluating functions on expressions. And for things that are essentially constant, but tedious to write otherwise, it is very useful: a bood example example is Base.Math.@evalpoly.

但它的问题是你不能检查只在运行时才知道的值:你不能实现 fib(n) = @fib_macro n 1，因为在编译时，n 是代表参数的符号，而不是您可以发送的数字.

But it has the problem that you cannot inspect values which are only known at runtime: you can't implement fib(n) = @fib_macro n 1, since at compile time, n is a symbol representing the parameter, and not a number you can dispatch on.

下一个最佳解决方案是使用

The next best solution to this would be to use

fib_eval(n::Integer) = eval(fib_expr(n, 1))

这可行，但每次调用时都会重复编译过程——这比原始函数的开销要大得多，因为现在在运行时，我们对表达式执行整个递归树，然后在结果上调用编译器.不好.

which works, but will repeat the compilation process every time it is called -- and that is much more overhead than the original function, since now at runtime, we perform the whole recursion on the expression tree and then call the compiler on the result. Not good.

所以我们需要一种方法来混合运行时和编译时间.输入 @generated 函数.这些将在运行时在 type 上分派，然后像定义函数体的宏一样工作.

So we need a way to intermingle runtime and compile time. Enter @generated functions. These will at runtime dispatch on a type, and then work like a macro defining the function body.

首先关于类型调度.如果我们有

First about type dispatch. If we have

f(x) = x + 1

并有一个函数调用f(1)，大概会发生以下情况:

and have a function call f(1), about the following will happen:

1. 确定参数的类型(Int)
2. 查阅函数的方法表以找到最佳匹配方法
3. 为特定的 Int 参数类型编译方法体，如果之前没有这样做的话
4. 根据具体参数评估编译的方法

1. The type of the argument is determined (Int)
2. The method table of the function is consulted to find the best matching method
3. The method body is compiled for the specific Int argument type, if that hasn't been done before
4. The compiled method is evaluated on the concrete argument

如果我们随后输入 f(1.0)，同样的情况将再次发生，为 Float64 编译一个新的、不同的专用方法，基于相同的函数体.

If we then enter f(1.0), the same will happen again, with a new, different specialized method being compiled for Float64, based on the same function body.

现在，Julia 具有可以将数字用作类型的特殊功能.这意味着上面概述的调度过程也将适用于以下函数:

Now, Julia has the peculiar feature that you can use numbers as types. That means that the dispatch process outlined above will also work on the following function:

g(::Type{Val{N}}) where N = N + 1

这有点棘手.请记住，类型本身就是 Julia 中的值:Int isa Type.

That's a bit tricky. Remember that types are themselves values in Julia: Int isa Type.

这里，Val{N} 是每个 N 的所谓 singleton 类型，只有一个实例，即 Val{N}() -- 就像 Int 是一个有很多实例的类型 0, -1, 1, -2, ....

Here, Val{N} is for every N a so-called singleton type having exactly one instance, namely Val{N}() -- just like Int is a type having many instances 0, -1, 1, -2, ....

Type{T} 也是一个单例类型，它的单个实例类型为 T.Int 是 Type{Int}，Val{3} 是 Type{Val{3}}-- 事实上，两者都是它们类型的唯一值.

Type{T} is also a singleton type, having as its single instance the type T. Int is a Type{Int}, and Val{3} is a Type{Val{3}} -- in fact, both are the only values of their type.

因此，对于每个 N，都有一个类型 Val{N}，它是 Type{Val{N}}.因此，g 将为每个单独的 N 分派和编译.这就是我们可以将数字作为类型发送的方式.这已经允许优化:

So, for each N, there is a type Val{N}, being the single instance of Type{Val{N}}. Thus, g will be dispatched and compiled for each single N. This is how we can dispatch on numbers as types. This already allows for optimization:

julia> @code_llvm g(Val{1})

define i64 @julia_g_61158(i8**) #0 !dbg !5 {
top:
  ret i64 2
}

julia> @code_llvm f(1)

define i64 @julia_f_61076(i64) #0 !dbg !5 {
top:
  %1 = shl i64 %0, 2
  %2 = or i64 %1, 3
  %3 = mul i64 %2, %0
  %4 = add i64 %3, 2
  ret i64 %4
}

但请记住，它需要在第一次调用时对每个新的 N 进行编译.

But remember that it requires compilation for each new N at the first call.

(而 fkt(::T) 只是 fkt(x::T) 的缩写，如果你不使用 x在体内.)

(And fkt(::T) is just short for fkt(x::T) if you don't use x in the body.)

最后来生成函数.它们是对上述调度模式的轻微修改:

Finally to generated functions. They work as a slight modification of the above dispatch pattern:

1. 确定参数的类型(Int)
2. 查阅函数的方法表以找到最佳匹配方法
3. 方法体被视为一个宏，并且使用 Int 参数类型作为参数调用，如果以前没有这样做的话.生成的表达式被编译成一个方法.
4. 根据具体参数评估编译的方法

1. The type of the argument is determined (Int)
2. The method table of the function is consulted to find the best matching method
3. The method body is treated as a macro and called with the Int argument type as a parameter, if that hasn't been done before. The resulting expression is compiled into a method.
4. The compiled method is evaluated on the concrete argument

此模式允许更改函数被分派的每种类型的实现.

This pattern allows to change the implementation for each type which the function is dispatched on.

对于我们的具体设置，我们希望调度代表斐波那契数列参数的 Val 类型:

For our concrete setting, we want to dispatch on the Val types representing the arguments of the Fibonacci sequence:

@generated function fib_gen{n}(::Type{Val{n}}, p::Real)
    return fib_expr(n, :p)
end

你现在看到你的解释是完全正确的:

You now see that your explanation was exactly right:

在第一次调用 fib_gen(Val{3}, 0.5) 时，参数函数fib_gen{Val{3}}(...) 被编译，其内容是完整的通过fib_expr(3, :p)得到的扩展表达式，即3*2*1*pp 替换为输入值.

in the first call to fib_gen(Val{3}, 0.5), the parametric function fib_gen{Val{3}}(...) gets compiled and its content is the fully expanded expression obtained through fib_expr(3, :p), i.e. 3*2*1*p with p substituted with the input value.

我希望整个故事也回答了您列出的所有三个问题:

I hope that the whole story has also answered all three of your listed questions:

1. 使用 eval 的实现每次都会复制递归，加上编译的开销
2. Val 是将数字提升为类型的技巧，而 Type{T} 是仅包含 T 的单例类型——但我希望这些例子很有帮助
3. 编译时间不是在执行之前，因为 JIT —— 它是每次第一次编译方法时，因为它被调用了.

1. The implementation using eval replicates the recursion every time, plus the overhead of compilation
2. Val is a trick to lift numbers to types, and Type{T} the singleton type containing only T -- but I hope the examples were helpful enough
3. Compile time is not before execution, because of JIT -- it is every time a method gets compiled first time, because it get's called.

这篇关于在 Julia 中使用元编程优化递归函数的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！