Before jumping into a discussion about how all of this works, lets look at an example of the DSL. In this example, we will define the sobel operator, then process an images with it. For a great overview of the sobel operator, check out this article.

In the DSL, there are things to worry about: Functions and Chains. A function is a single unit of work that takes an \((x,y)\) coordinate and an arbitrary number of inputs. For example, suppose we have a function \(Grad(x,y)\) that returns the magnitude of the gradient of two images \(I_1\) and \(I_2\) at the point \((x,y)\). We might denote this function with mathematical notation as:

\[ Grad(x,y) = \sqrt{I_1(x,y)^2 + I_2(x,y)^2} \]

In the DSL I have defined, we would denote this operation in a similar manner, sans syntactic differences:

// create a new function named grad
// Function::new takes a number and a lambda as an argument.
// The number indicates how many inputs the function has
// This lambda is always called with (x,y) coordinate values
// and an array of inputs of the length specified.
let grad = Function::new(2, |x, y, inputs| {
    // first we pull out references to the InputExpressions representing our inputs
    let input0 = &inputs[0];
    let input1 = &inputs[1];

    // compute the squares using the input expressions
    // Notice that x and y are both treated like functions.
    // This is essentially a hack to get around the way I've stored the AST
    let t1 = input0(x(), y()) * input0(x(), y());
    let t2 = input1(x(), y()) * input1(x(), y());

    // Compute the sum and the square root of the sum
    // The last expression generated by this lambda is result of the function we are defining
    // The Box::new trick is needed, again, because of the way I've store the AST
    Box::new(SqrtExpr::new(t1 + t2))
});

This isn't the most beautiful way to build a representation of our function, but it works and I learned a lot implementing the magic that makes it work. Each function stores is a syntax tree representing the expression that the function computes. The syntax tree defined by the grad function looks something like this:

There is also helper function that can be used to generate functions which perform a convolution on a single image with a kernel matrix For example, to generate a function that takes a single image as input and returns the convolution with the horizontal sobel matrix, use the following code:

let sobel_x = [[-1, 0, 1],
               [-2, 0, 2],
               [-1, 0, 1]];

let sobel_x = Function::gen_3x3_kernel(sobel_x);

Notice that these functions are defined in a purely functional, mathematical sense. They do not mutate their inputs, nor do they store any state, nor are they coupled to any particular inputs.

There's lots of little details which may be interesting to discuss, but I'm only going to discuss the compilation method. First, we need some slightly more rigorous definitions of things:

• All of the user defined functions take two 64 bit integer values and a list of inputs. They return a 64 bit integer value.
• Because functions themselves have no concept of chaining (the gradient function doesn't call sobel_x and sobel_y), the "inputs" to a user defined function can be thought of as function-pointers which will eventually be resolved to real functions, although this is not how they are implemented.
• All of the system-defined expressions (eg. SqrtExpr) take a 64 bit integer and return a 64 bit integer.

The compilation strategy for the DSL is very simple: Every DSL function is compiled into a function with a signature that would look something like this in C:

inline int64_t function(int64_t x, int64_t y, image inputs[], size_t num_inputs);

The array of image inputs provided here is not equivalent to the list of the inputs given to the DSL function. The inputs given to the DSL function are resolved to other compiled functions (using the chain) during code generation, so our generated grad function will directly call the sobel_x and sobel_y functions.

Since every value is a 64 bit integer, the code generation for an expression essentially just involves spitting out adds and multiplies for integers. The generated grad code roughly corresponds to:

inline int64_t grad(int64_t x, int64_t y, image inputs[], size_t num_inputs) {
    int64_t partial1 = sobel_x(x, y, inputs, num_inputs) * sobel_x(x, y, inputs, num_inputs);
    int64_t partial2 = sobel_y(x, y, inputs, num_inputs) * sobel_x(x, y, inputs, num_inputs);
    int64_t partial3 = partial1 + partial2;

    return core_isqrt(partial3);
}

A driver function is injected into the module. This function performs some bookkeeping tasks, then just loops over the pixels in the output image, calling the appropriate function (whichever was last in the chain) for every pixel:

for (int x = 0; x < output.width; x++) {
  for (int y = 0; y < output.height; y++) {
    int64_t res = function(x, y, inputs, num_inputs);

    /* output image at x, y */ = (uint8_t) res;
  }

Image inputs (the actual images we are processing), are passed to each function. When the compiler reaches an ImageSource in the function chain, it emits a call to a function which returns the pixel in the image at a given \((x,y)\) coordinate.

For anyone interested, I've dumped the entire LLVM IR module for an unoptimized sobel chain here. Some of the code is generated from the file core.c in the github repo for the project, if you need some hints to figure out what's going on here. The entry point is the function jitfunction. There's lots of other interesting little idiosyncrasies in this code but I don't have space and you don't have time to read about all of them.

Anyone who knows a little bit about computers and performance is probably hurting a little bit thinking about how this might perform. You've noticed all of the function calls, don't these have lots of overhead? Uou've noticed that I'm computing the sobel_x and sobel_y values twice in the gradient function. Don't worry, it isn't quite so bad.

Anyone who knows a fair amount about computers and performance noticed that inline keyword and is wondering if I'm somehow relying on function inlining to extract performance from this technique. The answer is yes. Every generated function is marked with an LLVM attribute AlwaysInline which, when combined with the appropriate LLVM optimization passes, guarantees that these functions will always be inlined into their caller.

Overall, this project was extremely enjoyable. I had yet another opportunity to fiddle with LLVM, which is always lots of fun (but sometimes very painful). I learned a little bit about image processing and some of the challenges that arise when shuffling pixels around. Finally, I learned a little bit of Rust. I have only one thing to say about Rust: Rust is an amazing language. Go learn Rust.