abstractcontrol

> It's still silly to compare Julia's tensors to PyTorch's because they are generic in dimension, layout, and type unlike PyTorch's and these seem to be exactly the pieces that you're ragging on. This means that your beef against PyTorch is fine, but it doesn't necessarily apply to Julia. I am curious if you can go into more detail about Julia's implementation instead of talking about possible Julia problems through an explanation of PyTorch.

Well, no, I can't go into more detail. It is the responsibility of those who take up the challenge I made, so why don't you go into detail instead, in the context of a GPU map kernel? I am genuinely curious about this, but not so much that I actually want to dive deep into the language.

Apart from that, I do like the rest of your reply. It is true that best performance lies somewhere between full inlining and full heap allocation, and that full inlining will kill compile times. One of the points I want to make with Spiral is that starting from the position of full inlining and then adding boxing and join points is a lot easier and much less error prone than starting from the other end and inlining by hand. I feel that this is important other programmers realize along with the fact that sophisticated compiler heuristics for when to inline are the only alternative.

> This sounds like a good idea. We should take this over to our community and see whether it would work well in the language.

I know that Odesky thinks that pattern matching on types directly only really works on languages with simple type systems such as Spiral's. If a language uses abstract types then that becomes more difficult. I definitely see a language having both intensional and extensional (parametric) polymorphism as a significant challenge. I'd love to see it met because extensional polymorphism is a lot easier on the IDEs and the way Spiral does it will make it difficult to create good tooling for it.

If Julia could support intensional polymorphism then inlining guarantees would just fall out from that. It is quite powerful - the inlining guarantees Spiral gives does not just apply to top level functions, but any function passed as an argument and captured by other functions.

> Not really. Julia is dynamic and interactive. C++ with its full templating is very powerful but can allow its compile-times to explode because its main use case is not with a REPL (and templated C++ code compile times are...).

If C++'s templating system was a car, then it would be missing wheels and a steering wheel. And possibly the brakes. How can one possibly do staging without join points, and all the other things Spiral has? As a language it misses several key concepts that would actually make such a system work well.

With the last paragraph I half agree with you.

It is true that right now Spiral uses a compile-run loop and acts as a sophisticated code generator, but it could be made into caching interpreter with some work. It probably will have to be if it gets bigger for the sake of compile times and library support.

Getting instant feedback in the REPL is hardly the domain of dynamic languages anymore - F# can do it well for example and so can various other static functional languages. Even C# can do it these days.

> To me, Spiral looks like a cool alternative to C++/Rust/Go in the scientific computing domain, and package ecosystems that mix Spiral+Julia could be quite powerful.

Julia does a lot of multi-stage compilation. It probably has had a lot of work done on its JIT by now. And language mixing, especially of languages with different strengths and weaknesses is a great idea. Staging gives the user access to it in a more direct manner so why not bring it out?

Spiral is quite a small language so if JITing abilities are flexible enough, it might be a good idea to implement it in Julia as a DSL and merge it with main language. That would combine all the advantages of Julia and Spiral as languages.

I cannot change my attitude otherwise I would never get anything done. The only principle that I hold dear, and the only way I can truly be fair is that if I am actually wrong I will change my mind. I would not be so rough on Julia if it didn't present itself as a direct competitor to Spiral. Since it presents itself as a numerical computation language capable of replacing C++ and dominating machine learning and GPU programming I will hold it to the same standard I would take Spiral.

Any numerical computation language should have tensors generic in dimension, layout and type. That is the point I am trying to hammer home. If you think that Julia's tensors are better than PyTorch's and feel that my equal treatment of them is unfair - I honestly feel that is the wrong way to think about this. Rather than fighting to be the king of a very small hill, it would be better to try climbing a much larger mountain instead. Look inward and be unhappy at what the language has now.

Envy Spiral's wonderful tensors. Wouldn't Julia be a much better language if it had its inlining guarantees? And maybe add to that the intensional polymorphism that would allow it to implement them.

Why if that happened, I'd lose my reason for working on the language and have to do something else like study reinforcement learning for example. How horrible would that be? :)

How would one enable this in Visual Studio?

I might be persuaded to see things from your point of view if you can show that Julia's tensors meet my criteria. I might be wrong about this - I do not know much about Julia and it might very well be capable of this. But these are the kinds of things that if it were capable of, it would advertise. If it can do it I'll apologize and change the line to show Julia in a more positive light.

The reason why I started talking about inlining of monads here is because ultimately they are just functions, and Spiral's tensors are not part of the language, but just functions as well. So if it cannot do monads then it probably cannot do tensors to satisfaction either.

It was pointed out to me that tensor slicing could serve as substitute for partial application. What I am wondering next is the degree of flexibility Julia's tensors have in the context of GPU kernels. Let me illustrate this with a short example. The following is how the inplace map Cuda kernel is made in Spiral. It seems simple, but I'll go at it from the top to show the characteristics and flexibility that Spiral's tensors have.

    met map' w f in out = 
        inl in, out = zip in, zip out
        assert (in.dim = out.dim) "The input and output dimensions must be equal."
        inl in = flatten in |> to_dev_tensor
        inl out = flatten out |> to_dev_tensor
        inl in_a :: () = in.dim
        
        inl blockDim = 128
        inl gridDim = min 64 (divup (s in_a) blockDim)

        w.run {
            blockDim gridDim
            kernel = cuda // Lexical scoping rocks.
                grid_for {blockDim gridDim} .x in_a {body=inl {i} ->
                    inl out = out i
                    inl in = in i
                    out .set (f in.get out.get)
                    }
            }

The first line is `inl in, out = zip in, zip out`. Since the inputs to to the map function can be arbitrary data structures, the `zip` iterates over it and zips the tensors in it into a single one.

In Spiral the tensors have a tuple of arrays layout by default, and it is absolutely trivial to switch between AOT and TOA representations. In fact you can mix and match tensors that have both layouts internally. Zipping tensors in Spiral just merges them into one with the TOA layout though the subtensors could be AOT.

Dealing with these layout transformations is a significant issue in numerical computation that Spiral solves completely. The reason why I have an impression that Julia cannot do this is not because I know for sure. Rather I do not think this can be done within the confines of a parametric type system. Rather I think it would require an intensional type system like the Spiral has.

    inl in = flatten in |> to_dev_tensor
    inl out = flatten out |> to_dev_tensor

`flatten` checks that the tensor is contiguous and flattens it into a 1d tensor. The `to_dev_tensor` is just some formality that I could not get rid of. It takes the pointers out from the references holding them which allows the tensors to be captured lexically before being passed into the kernel. Without this a type error would occur.

    kernel = cuda // Lexical scoping rocks.
        grid_for {blockDim gridDim} .x in_a {body=inl {i} ->
            inl out = out i
            inl in = in i
            out .set (f in.get out.get)
            }

The actual kernel itself is quite small and trivial. `cuda` is just a parser abbreviation for `inl {blockDim gridDim} ->`. In other words it is not a special language feature, but a regular function.

When I was doing this in F#, I had to not pass in every single argument of the function by hand - meaning the tensor pointers, their offsets, their sizes, their dimension ranges for each and every tensor, I also had to make wrapper classes on top of that. If a map operation needed more arguments, then I would have needed to copy paste a kernel, adjust it and the wrapper and repeat the testing process all over again, all this with very little type safety.

Spiral is not quite on the level of Futhark in terms of ergonomics, but it is a vast improvement over what came before. Note that there is no need for a type annotations anywhere - the function is as generic as it could possibly be in a statically typed language. Because of inlining guarantees all of this has zero overhead.

Compared to Spiral, where do Julia and its tensors stand on the axis of generality? Would it be possible to write a simple map as simply as this?

I recall an example where Julia had to use macros in order to achieve loop unrolling. I'll dig it up if you want. This is something can be done quite naturally in Spiral due to its type system.

I find that people often say macros are about syntax, but what they really are is some combination of a parser and a partial evaluator - in other words, their main use seems to doing compile-time function evaluation.

Spiral does not have arbitrary CTFE, but can perform any functionally-pure, non side-effecting computation in its type system plus memoization of function calls.

One thing of practical interest that makes me doubt Julia's inlining capabilities would be monads.

https://github.com/pao/Monads.jl/blob/master/doc/index.rst

"Monads.jl provides a powerful, if relatively slow"

The key word here is slow. I haven't seen a single language as of yet apart from Spiral that can inline all their overheads away.

They are pervasive in Haskell, and Haskell cannot do it. Scala wants to do it, but it can't. F# and Ocaml cannot optimize their overheads away. If Julia could then it would advertise it. This is something that is subject of academic research in the field of partial evaluation.

Speaking as a functional programmer, as a language Julia feels barely functional to me. It has first class functions, but it does not have pervasive partial application and syntax that makes such a style comfortable like functional languages tend to have. It feels more like an imperative language.

And speaking as a systems programmer, I know that without inlining guarantees you cannot have optimized functional constructs. It is not just monads. By tracking all the factions and their bodies to exactness and forgoing all heap allocation unless explicitly stated, it becomes much easier to do language interop. It is possible to make functions that capture variables from lexical scope and then pass as arguments to a function that passes them through to the GPU. Julia talks about speed, but barely talks about inlining so I am dubious about its claims.

I could go on, but I think this is a large enough sample of my views on Julia. If Julia wants me to use it, then it needs to speak more to my values as a programmer.

Never stumbled upon it in my research. Thanks for the reference.

Yes, it is, but there exist an alternative approach called automatic differentiation which is more imperative and flexible in nature than symbolic differentiation.

The idea is to keep a tape (you can think of it as a list) and then record all the operations on it as you step forward through the program. Then at the end you execute the operations backwards from the tape.

I take this approach in Spiral's ML library.

I can't actually program all the time. Sometimes I need to study machine learning or think about design. During those times I use commits much like journal entries.

I actually keep a separate journal and sometimes paste from it. I've been using LibreOffice Writer for it and decided to drop it recently because it would take so long to save a file - like 5s or more. Around every 3 month I'd fill in about 1000 pages of it and it was forcing me to move to a new file every time that happened. Now I just use VS Code for this sort of thing. Raw text is the best after all, but it is too bad I can't paste images into the journal anymore.

The reason why I am using random projections in the latest test is because I am testing an algorithm that iteratively calculates the inverse Cholesky factor of the covariance matrix and am testing it on Mnist images. The cov matrices made from raw Mnist images are non-invertible, but projecting them to a much smaller dimension allows me to actually test the algorithm on non-synthetic data.

I do not actually need more than I have, but I'll keep your link in mind if I ever need random projections though.

Before I started work on Spiral I spent way too long, like a year working on a ML library before throwing in the towel. I think I skimmed through the source of DeepML over a year ago and my impression was that the author was struggling against the limitations of F#'s type system. There were a bunch of places where he was doing the equivalent of telling the compiler to fuck itself by downcasting to `System.Object`.

The way these projects go is that you start of with an AD library and then when the type system and the language starts getting in the way you either do a lot, and I mean a LOT of engineering with inferior tools or you make something better.

The 'make something better' is always picked and can take various forms. Generally, people would not make a decision to work on a language for a ML library for almost a year before they write down the first line for it. They instead take the middle road of starting with an AD library and gradually extending it.

First they realize that they cannot really express tensors properly with the confines of the language so they make everything symbolic, then they build engines to run such symbolic AST, then they realize that they need JITs and other optimizers to make everything run fast. Tensorflow and PyTorch are currently at this stage. There are more stages after this.

Of course the task of working with ASTs is what a compiler does. I think it is a great pity that Tensorflow and PyTorch are written in C++ under the hood. C++ is probably the worst language I can imagine for working on compilers, so I applaud the authors of DeepML for picking F# instead. Had the TF team picked statically typed functional language for this they could have cut the size of TF by over 10x and saved themselves over a million lines of code.

[ML](https://en.wikipedia.org/wiki/Standard_ML) derivatives like it were made for that sort of thing and I cannot imagine doing Spiral in a non-ML styled language.

Indeed. But F# cannot actually do this, and neither can even dependently typed languages to the degree Spiral can. In Spiral types act much like values meaning you can bind them to variables and you can reflect on them without cost like with values, though if you tried something like `(type 1) + 2` you'd just get a error during the last phase of compilation. It is rather useful and without it, type inference would be impossible in Spiral.

I do not think there is a language with capabilities similar to Spiral. If it existed it would have saved me a lot of time. The closest you could get would would be [Scala LMS](https://scala-lms.github.io/) and the work done on staging such as in [context of Ocaml](http://ocamllabs.io/iocamljs/staging2018.html).

Personally I think macro-style staging in statically typed languages is an atrocity and want nothing to do with it. I kind of like the type driven version as is done in LMS, but compared to Spiral the LMS is less expressive.

A lot of work is being done on macros/partial evaluation/staging in Scala and I assume it will continue to be done as the authors are aware of the need for this in a language - without it a functional language would be dog slow, but my impression is that they are having significant difficulties making it fit with Scala's much more complex type system. One consequence of that is that they seemingly rewrite their macro system every few years. This time they will surely get it right despite being at it for decade(s).

The reason why what Spiral has works so well is because its type system is so simple, so staging can be completely intervowen with how it does type inference.

My view based on seeing the last 50 years of language design is that staging is really difficult to bolt onto an already established design. It tends to degenerate into having two languages - the core language and a scripting language on top of it. C++'s template system is the worst example of this. Another example which actually embraces this dichotomy would be [Terra](http://terralang.org/) which could be described as bolting Lua on top of C. I rather dislike this sort of design.

A language needs to be designed around staging from the ground up and it needs to be accepted that in the absence of magic, the programmer will have to mold his style to fit with the new paradigm.

If a language is to have staging it should be everywhere and at all times much like regular type inference is in F# and other ML styled languages.

When I looked at Julia last time I do not think it had the capability to partially apply tensors. In Spiral indexing into a 3d tensor would give you back a 2d tensor. It might be possible to do this with function in Julia, but that would depend on its inlining capabilities and unlike Spiral it does not provide guarantees with regards to that.

Its type system is definitely less powerful than Spiral's - some of the things you could express naturally in Spiral would require macros in Julia.

All in all, the language strikes me as a better Matlab which I think it succeeds at, but Spiral is made to be a better F#. Some of the tradeoffs I had to make in language design means that I could not succeed at this across the board.

The fact that it has a more powerful type system makes it take an immediate hit in terms of ergonomics, and the language is less pleasant to code in than F# because the IDE is not providing immediate and constant feedback.

On the other hand for systems programming, I can't imagine what C (and even C++) could possibly do better than it except maybe compile times. Despite that, Spiral is also nothing ilke Rust - it does not try to check that pointers are safe. It is rather a extremely expressive static language that has no type annotations anywhere and looks like a dynamically typed functional language.

Author here. The language was made for Cuda programming and until something comes along to displace it, I have no plans of moving from it.

That having said, the language is quite small and it would be possible to create a new backend in about 500 lines of code so if you really want it, you could probably do it yourself. I'd help you integrate it into the language. I know nothing about SPIRV at the moment though.

> It does not seem that there is any significant correlation between memory size and allocation cost. Nor should there be - it's just bookkeeping, after all.

I guess this explains your lack of concern about memory pooling. I am going to have to do some research to see if OpenCL does pooling natively, but I can assure you that this is not the behavior on Cuda devices. I can't really get your example to compile right now, but I'll look into it later.

Also you are right that is strange that allocations are so slow on Cuda.

But if it turns out that OpenCL is doing pooling behind the scenes, this is going to be an issue for you when you decide to do a Cuda backend because the performance profile will start to get dominated by allocations.

> This would have the same effect as embedding assembly into C code, i.e. it basically makes everything nearby off-limits for optimisation.

This is true, but that would also be the case if you were using an FFI. Being able to embed C code would shorten that path a bit by not requiring the user to compile a code fragment to a .dll before invoking it.

> In fact, Futhark is a bit of a counter-reaction to the run-time-code-generating embedded languages that were popular in the FP community some years ago.

Hmmm...I never heard about this. But Then again, relatively speaking I haven't been programming that long.

> Have you considered looking at something like Obsidian[2]? It is a much lower level functional language that directly exposes GPU concepts. While Obsidian is Haskell-specific, there is a successor called FCL[3] that is standalone. Unlike Futhark, these languages do not depend on aggressive and expensive compiler transformations, and therefore may be a better fit for embedding.

I can't find any documentation for Obsidian, but there is a paper for FCL in the Github, so I'll look at that.

> Do you think some of the pain would be alleviated if Futhark had a simple FFI that would allow calling hand-written primitives when they exist?

Yeah, definitely. You should go a bit further and allow embedding C code for performance oriented users. I do not foresee using it, but somebody is going to need it eventually I guarantee it. I can envision some better alternatives than C, such as the language I am working on, but right now it is still incomplete and won't be for some time.

> first, memory allocations on GPUs are not unusually slow

The time it takes to allocate a chunk is linear in its size which is quite slow. I am not sure why that is, but maybe it is faster using OpenCL instead of Cuda? What were your timings for allocating and disposing raw memory plotted against size?

> If used as a library, a Futhark program does not dispose everything once it stops running, but only once the library is unloaded.

That is interesting. Can Futhark be used as a library apart from Haskell (in which it is written)?

> For example, if a Futhark function returns an array, that array still lives on the GPU, and if you use it as an argument to another Futhark function, there will have been no traffic (except bookkeeping stuff) between CPU and GPU.

But still, Futhark will probably not be able to optimize away all the intermediates. And more to the point, some programs like neural nets do in fact accumulate intermediates by necessity. If a particular Futhark program is run multiple times, the memory in those intermediates should be held in a pool.

> Second, I don't see why an interpreter would be any better at managing memory than the current Futhark runtime system.

It could potentially allow memory pool to be shared amongst multiple Futhark programs.

The idea is not to turn Futhark into an interpreted language per se, it would still be a compiled language, but to instead add an extra layer that would allow easier communication with other languages. I do not have a concrete vision of how this should be done.

It goes back to what you mentioned about using Futhark as a library. I am expecting a negative answer that it can be used as a library from anything other than Haskell, but if you were to go more in the direction I am suggesting, instead of making backends for C#, Java and such what you could do is make something that will allow Futhark to be used as a library.

Thinking about to some of the C examples that I have seen, I do not think users will appreciate having massively bloated code files needed to compile the stuff dumped into their projects folders by Futhark.

Not to mention, C# will need to be compiled to C which will result in more temporary files. Futhark as an embedded language could take responsibility for managing all of that.

I do not feel like writing an entire essay on language integration right now - most likely I will open an issue 6-12 months from now on Futhark's Github repo just to talk about it, but let me start with agreeing everything what DannyBee said and adding a few thoughts of my own.

1) Let me just say that language integration is a very serious issue. It is not just between completely different languages, but also between modules written in the same language. Once you start using algebraic datatypes to emulate language features the main language lacks, you essentially step into a dynamic sublanguage and have to deal with the friction caused by crossing module boundaries. This friction is both on the programmer's side – he has to deal with writing boilerplate for crossing the boundaries, and on the computer's side which has to do marshaling which is terrible for performance and just nasty.

2) This friction in the case of Futhark will be magnified manifold as it is a completely different language. Since Futhark is not a general purpose language, a realistic use case for it to be called indirectly directly by other languages who will generate code for it. This is generally how high-speed anything is used today. There is a collection of highly optimized, assembly written routines (such as BLAS) bundled into a library and they are called from very slow high-level languages such as Python.

Futhark today is not fit for such a purpose.

* It would be difficult to partition the program written for it into separate pieces. For very simple programs, at a minimum it will generate 2k lines of code (in the C backend).

* This is compounded by the fact that it does not link to the aforementioned optimized libraries, but generates all the code internally. You could then imagine using Futhark intermittently – calling those fast libraries in the main language and using Futhark for the rest since writing code in Futhark is much more convenient compared to C, but then you would need to partition the program and will immediately run into the code bloat issues in the first bullet point.

* Futhark is a very high level language and takes on all the responsibility for managing memory on itself. That even further clashes with the idea of partitioning the program. Memory allocations are extremely slow on the GPU, and in addition to that, they block the whole device meaning they are not asynchronous.

* A minor point of friction compared to the above is that Futhark supports OpenCL which has minor market share instead of Cuda.

3) Based on the above, I question the current integration strategy by Futhark of making backend for different languages. It currently has a C and a Python backend, and an OpenCL CPU backend, and F# is planned, and you can imagine many different backends…

What might be worth trying instead would be to make Futhark an embedded interpreted language. This is not as crazy as it seems – it would define a natural API point for other languages to access it and the interpreter would be responsible for managing GPU memory. It would be a much better model than disposing everything once the program stops running and would allow for efficient intermingling of multiple Futhark programs that could reside in memory. Right now, that sort of thing would be very high friction.

I do not have much advice on how this could be accomplished and would no doubt require much design work, but I am going to try something like that in my own language at some point. I had this crucial insight when I was trying to use it from the language it was written in and realized that it is actually very difficult.

Scala, Clojure and F# in particular had the master stroke of latching themselves to already established ecosystems. Languages targeting the GPU cannot use that strategy directly and will need to be more inventive.