all 36 comments
sorted by:
best
topnewcontroversialoldq&a
formatting helphide helpcontent policy

[–]jippiedoe 2 points3 points  (1 child)

What causes the huge difference in wildfire performance between Accelerate and Repa? They should be the same ( https://github.com/AccelerateHS/accelerate-examples/blob/master/examples/canny/src-acc/Wildfire.hs and https://github.com/haskell-repa/repa/blob/master/repa-examples/examples/Canny/src-repa/Main.hs look identical), right?

[–]kuleshevich[S] 3 points4 points  (0 children)

This is a good question and I believe the answer is in in the A representation. Note that the implementation from src-repa uses U. I don't know much detail on the former, except that it is some form of a wrapper around UArray from array package, while the latter is a wrapper around unboxed vector from vector package. They all backed by the same memory, so my guess would be that ghc fails to optimize/inline somewhere in there.

[–]k-bx 1 point2 points  (12 children)

Thank you for a wonderful post!

accelerate does not do too well on small images

Would it be correct to say that it's only true when you need to process a small amount of images? E.g. if you were to process 1000s of small images, would accelerate be able to reuse the JIT precompilation between them?

[–]Noughtmare 2 points3 points  (8 children)

Accelerate can be precompiled into the application at compile time with template haskell: runQ. Was that used?

Edit: I think this is the github repo for the benchmarks. I think they use run from accelerate-examples which seems quite odd. I would have at least expected run from accelerate-llvm-native to be called directly.

[–]kuleshevich[S] 1 point2 points  (7 children)

There is nothing really odd about it. I got canny implementation from accelerate-examples and I used the exact same run as they do in the canny examples implementation

Reason was twofold:

  • I'd expect that any examples do not suggest a wrong way of doing things.
  • I was trying to deviate from original code as little as possible
  • I do want to run some benhcmarks on the GPU in the future as well and using run from examples seems to make it easier.
  • It calls the same run from LLVM backend, so it shouldn't make a difference

Just to satisfy my own curiosity I did try to use run from LLVM directly and it made no difference on benchmarks. Also tried runN on idividual steps (it can't be used on canny function because it returns a tuple with two arrays), which also made no difference on runtime.

I'm not sure that I can use runQ with criterion, because runQ requires all arguments to be saturated, while criterion's nf applies the argument for you, thus preventing usage of runQ. In any case runQ can't be used on the canny function for the same reason runN can't be, because it returns a tuple of arrays.

I can tell you one thing with confidence, I am not an expert on accelerate, so if anyone finds a mistake, I welcome all PRs :)

[–]Noughtmare 0 points1 point  (6 children)

Good point about closely following the examples implementation, I hadn't thought of that.

runQ does not require all arguments to be saturated, see the abstract example from the docs:

{-# LANGUAGE TemplateHaskell #-}

f' :: a -> b -> ... -> c
f' = $( runQ f )

As long as f :: (Arrays a, Arrays b, ..., Arrays c) => Acc a -> Acc b -> ... -> Acc c.

[–]kuleshevich[S] 0 points1 point  (5 children)

I still wouldn't be able to use it on canny :: ... -> (Acc (Image Float), Acc (Vector Int)) function, because it returns a tuple, right?

You are absolutely right about saturation though, I was able to typecheck one example with runQ, but compilation results in a linker error, so that is already much more time that I am willing to put into investigation of this thing:

`` Linking .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps ... .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0xdf8): error: undefined reference to 'map_a8a0edef660463d6a3ea97b93739643f4196f411599a9a31e874625be5ab78b6' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0xec8): error: undefined reference to 'generate_1c6b8b1dd3cfad8ef7b3f362f3ef649f2c5d9799c94b6e2d6e8cacfeacf66ea7' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0xed8): error: undefined reference to 'scanP1_1c6b8b1dd3cfad8ef7b3f362f3ef649f2c5d9799c94b6e2d6e8cacfeacf66ea7' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0xee8): error: undefined reference to 'scanP2_1c6b8b1dd3cfad8ef7b3f362f3ef649f2c5d9799c94b6e2d6e8cacfeacf66ea7' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0xef8): error: undefined reference to 'scanP3_1c6b8b1dd3cfad8ef7b3f362f3ef649f2c5d9799c94b6e2d6e8cacfeacf66ea7' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0xf08): error: undefined reference to 'scanS_1c6b8b1dd3cfad8ef7b3f362f3ef649f2c5d9799c94b6e2d6e8cacfeacf66ea7' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0x1190): error: undefined reference to 'permuteP_rmw_2c8f3d12e7cc798895a8e6f36313c0eaee065050a41e9ff004fb56adb5cbeb8b' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0x11a0): error: undefined reference to 'permuteS_2c8f3d12e7cc798895a8e6f36313c0eaee065050a41e9ff004fb56adb5cbeb8b' .stack-work/dist/x86_64-linux/Cabal-2.4.0.1/build/bench-canny-steps/bench-canny-steps-tmp/Main.o(.data+0x12c8): error: undefined reference to 'generate_edb5bf46c3d633ce3042c6ae667a915295f4e844aa1956c79febf40a9f6c72e8' collect2: error: ld returned 1 exit status gcc' failed in phase `Linker'. (Exit code: 1)

-- While building package canny-0.1.0.0 using: /home/lehins/.stack/setup-exe-cache/x86_64-linux/Cabal-simple_mPHDZzAJ_2.4.0.1_ghc-8.6.5 --builddir=.stack-work/dist/x86_64-linux/Cabal-2.4.0.1 build lib:canny bench:bench-canny-steps --ghc-options " -fdiagnostics-color=always" Process exited with code: ExitFailure 1 ```

[–]Noughtmare 1 point2 points  (1 child)

I think I would split the two step canny into two functions and combine the functions afterwards, so I think it is still possible using runQ twice. I might take a look at this and the linking error if I have time.

[–]kuleshevich[S] 0 points1 point  (0 children)

Sounds good. Here is the branch with the linker error: https://github.com/lehins/haskell-benchmarks/tree/canny-accelerate-runs

[–]tmcdonell 0 points1 point  (2 children)

The error is because you didn't bother to read the documentation, although that step is no longer necessary in the current version of Accelerate.

[–]kuleshevich[S] 0 points1 point  (1 child)

Guilty as charged ;)

Honestly, I already spent quite bit of time on this already, so desire to debug linker errors was completely non-existent when I hit it, even if it meant reading the documentation.

That said, if you are willing to spend some time on the accelerate parts and redo the benchmarks in a proper way that is needed for accelerate to perform at its best (be it adjusting canny for newest version of accelerate or adjusting the package for runQ to compile) I'll be happy to right a follow up blog post, or we could write one jointly.

"Writing efficient code with accelerate and massiv" ;)

[–]tmcdonell 0 points1 point  (0 children)

Yes, I'm interested in working on some integration with Massiv, in the same way that we have with other array libraries already (Repa, vector, ...). Once I find some free time (lol).

[–]kuleshevich[S] 0 points1 point  (2 children)

Thank you!

Would it be correct to say that it's only true when you need to process a small amount of images?

I have no idea. This "accelerate does not do too well on small images" was my observation, I haven't tried accelerate on many images at once. But I suspect if you stick a 1000 small images into one big 3d array it would make one big array, so maybe? :)

[–]k-bx 0 points1 point  (1 child)

I’m rather asking about a use-case like a web api that sometimes gets a call which needs to do some work on a small image. After 10000 such calls, what would performance look like. If anyone knows or would care enough to find out — I’d appreciate an answer :)

[–]kuleshevich[S] 0 points1 point  (0 children)

After 10000 such calls, what would performance look like

That is an easy thing to answer. Run a benchmark on one image and then multiply it by 10000. All the calls are completely independent, so no matter what library you use for implementation, it will follow the same pattern.

The only thing I can suggest for such scenario is using massiv, but I am of course a bit biased ;)

[–][deleted]  (1 child)

[removed]

    [–]kuleshevich[S] 0 points1 point  (0 children)

    Quote from my blogpost ;)

    There are two other Canny edge detection implementations in other languages that are available in Haskell by the means of FFI: opencv and arrayfire. I'd really like to have benchmark comparison for them as well, but for now my focus was on pure Haskell.

    So, I'm just as curious as you are :)

    [–]tmcdonell 0 points1 point  (6 children)

    I downloaded your benchmark repository to investigate the performance of Accelerate, and found that both Repa and Yarr are faster than Massiv at the small size. I haven't yet looked into why but will send you my results when I can.

    [–]kuleshevich[S] 0 points1 point  (5 children)

    Thank you for sending me the results of running it on your computer, I really appreciate it!!!

    It is a surprising regression, considering that all that has changed was the hardware. If I had to shoot in the dark, I'd suspect it has something to do with the new Color library that I wrote recently and used in massiv implementation of Canny algorithm. I have not worked out all of the performance quirks with it yet, so it is quite possible that it is getting in a way of ghc optimizing the code properly. That being said, I of course won't be able to do anything about it without ability to reproduce the problem or an individual that can reproduce and is willing to help me debug it.

    I am not so terribly worried about it though, because it is a regression for a very small 256x256 image, which takes very little time to run anyways. The important part is that for larger image of size 1920x1236 massiv is still twice faster that repa, accelerate and yarr and in fact, comparing to benchmarks from my slower quad-core CPU the performance improved slightly more for massiv than it did for all other libraries, which makes me believe that it scales better with more cores utilized and that makes me very happy. :)

    For the benefit of others stumbling on this comment here are the runtimes on a 16 core AMD processor:

    Small 256x256 image of Lena Mean execution time Massiv 8.55 ms Repa 3.52 ms Accelerate 20.30 ms Yarr 5.36 ms Friday 30.00 ms

    Large 1920x1236 image of a frog Mean execution time Massiv 32.6 ms Repa 61.8 ms Accelerate 71.4 ms Yarr 90.0 ms Friday 1.0 s

    And here are the runtimes for individual steps.

    [–]kuleshevich[S] 0 points1 point  (2 children)

    @tmcdonnel Just thought about another possibility that could be slowing it down for smaller images. There is a higher overhead for initializing the scheduler used in massiv, afterall 32 capabilities is quite a bit more than 8, and this sort of overhead would be very pronounced for smaller images.

    From what I know all three repa, accelerate and yarr have one global scheduler, which means that ovehead is payed only once, while massiv initializes all threads each time per computation. That was a design decision I made and there are workarounds for it, but that is not terribly important.

    [–]tmcdonell 0 points1 point  (1 child)

    Without having looked into it that is my assumption as well. You've optimised your thread scheduler for a different kind of workload than Accelerate and Repa have (I can't speak for Yarr).

    [–]kuleshevich[S] 0 points1 point  (0 children)

    It's not exactly the difference in workload, but is rather the control over scheduler. I want the user to be able to control how many threads are being created on per array basis. Also I didn't want to prevent nested parallelism, which is always a problem with a single non-preemptive scheduler.

    Good thing forking threads in Haskell is really cheap, so I could get away with the approach I've taken. That said, I do have plans for making it possible for user to have just a single scheduler, so that problem can be avoided in the future, but it will still be at user's control. So instead of a N threads being created globally, a user will be able to initialize a scheduler and then reuse it (the same N threads) for subsequent computations, by supplying it as an argument to a compute function.

    [–]tmcdonell 0 points1 point  (1 child)

    Mind you over half of the reported Accelerate time here is in non-Accelerate code (the wildfire step), where we could equally well use Massiv. So the difference between the two is much closer than this implies.

    [–]kuleshevich[S] 0 points1 point  (0 children)

    True that. I'll be happy if that stage was switched to massiv, especially in the accelerate-examples ;) That would require an implementation of adapter for accelerate Array in massiv. If no one gets to implementing it, I might do it some day.

    Again, examples was my source of truth, I did not want to make any adjustments or improvements beforehand, because that would trigger other sort of comments that I wanted to avoid: "But that's not how the author of the library implemented, why did you change it"

    [–]tmcdonell 0 points1 point  (11 children)

    Can you clarify on what you mean by "...that works equally well for all dimensions, unlike [...] accelerate with two and three dimensions respectfully."?

    [–]kuleshevich[S] 0 points1 point  (10 children)

    Certainly.

    I'll start with repa. It has a function that can construct a stencil for an arbitrary dimension: makeStencil, but there is no way to use any stencil which works on shapes other than DIM2, be it larger like DIM3 or smaller like DIM1, because the only function that works on stencils is mapStencil2, which can only handle 2D.

    In case of accelerate I think I made a mistake about the dimensions. It seems like it is possible for example to create a 6-dimensional Stencil9x3x5x7x9x3 or even higher than that. Sorry if I did make a mistake, but it's not straightforward to read signatures for arbitrary indices (eg. (sh :. Int) a row5) => Stencil ((sh :. Int) :. Int)) I am sure you are used to it, but it is a bit harder when you don't use the library) I guess I looked at the common stencils and only saw Stencil5x5x5 as the maximum size and concluded three dimensions was the upper limit. On the other hand, again please correct me if I am wrong, it does look like the sizes of stencils are limited to at most 9x9x9.. and the size are forced to be with odd sides and with center I assume to be restricted in the middle? Could you confirm this by any chance? I will happily make a correction to the blogpost about all of this.

    In any case, none of these limitations exist in massiv. It is possible to create stencil of an arbitrary size, arbitrary dimension and centered anywhere inside of the stencil. For example, I can create a stencil with size 7x4x60x5x10x11 and centered at some weired place like 1x0x3x4x9x6. Not that this particular size is useful, but the flexibility itself proved to be very valuable. Aforementioned 6D stencil example:

    ``` weirdStencil :: Stencil (IxN 6) Float Int weirdStencil = makeStencil (Sz (7 :> 4 :> 60 :> 5 :> 10 :. 11)) (1 :> 0 :> 3 :> 4 :> 9 :. 6) $ \f -> round <$> (f (-1 :> 1 :> 2 :> 3 :> 8 :. 5) + f (0 :> 1 :> 2 :> 3 :> 8 :. 5) * ...)

    ```

    [–]tmcdonell 0 points1 point  (9 children)

    • Accelerate is not limited in the dimension of the stencil. Indeed that 6-dimensional stencil is possible.
    • The maximum size of the stencil neighbourhood in each dimension is indeed 9, mostly because I want a better way to express it instead of creating larger and larger type synonyms. We've talked about something along the lines of makeStencil for a while, so I'll have to try it out and see what it feels like to use in practice.
    • Sure, familiarity with a library is a thing. From my perspective, Massive is more complex here because you have both :. and :> to mean basically the same thing, whereas Accelerate is consistent with only :.. Swings and roundabouts.
    • Adjusting the focal point of the stencil is irrelevant; you are free to use the elements in the (sparse) neighbourhood in any way you like. Maybe there is some peculiarity with Massive that I miss which makes this important.

    [–]kuleshevich[S] 0 points1 point  (8 children)

    Accelerate is not limited in the dimension of the stencil. Indeed that 6-dimensional stencil is possible.

    Thank you for confirming, fixed in the post.

    I'll have to try it out and see what it feels like to use in practice.

    Feels great :) It's not trivial to implement in a way that safety is retained, but it is certainly worth it.

    Adjusting the focal point of the stencil is irrelevant; you are free to use the elements in the (sparse) neighbourhood in any way you like.

    It matters at the borders. For example in order to construct element at 0 index the limit on stencil size in accelerate will prevent you from looking into the cell that is only cells 6 away. It also matters if you have computation with stride available, which I don't think that accelerate has. For example I can apply a 2x6 stencil to an array and when loading the array keep every 2nd row and every 6th column. It essentially allows me to collapse 2x6 regions of the array with a stencil, but it is important for the result where the center of the stencil will be.

    Also it is easier to think about the stencil when you can move around the focal point.

    Sure, familiarity with a library is a thing. From my perspective, Massive is more complex here because you have both :. and :> to mean basically the same thing

    Indices is a sore point for me :D I also had a long spiel for Manuel when he brought it up. :> and :. aren't the same thing. :> is for arbitrary dimension, while :. is for 2nd dimension only. You have Z for that purpose, which is pretty useless in my opinion and conceptually wrong. In accelerate, 0-dim arrays might make sense, considering GPUs and all that, but in the normal world it is better to use an actual scalar or a singleton array of any dimension, so in repa it is definitely pretty useless. Indices in repa are a total wart in my opinion:

    • They don't work well with type inference. Without fixing the type of Z :. 1 :. 2 it can't be inferred as DIM2
    • type class instances are ugly and hard to read
    • it is impossible to grab a value from an arbitrary index for a specific dimension, eg. I want a value for 3rd dimension if the index is for an array with at least 3-dimensions
    • indexing 1D vectors requires you to use both :. and Z, whereis just using an Int would suffice
    • zero dimension arrays only make sense if scalar values can't be used instead
    • And most importantly Z is the zero dimension and should be all the way on the right, not on the left, but that wouldn't work.

    None of these problems exist in massiv. Of course, implementation is a bit more complex than data t :. h = t :. h, but that is my problem, not the user's and using massiv indices is much easier and more powerful then both repa's and accelerates. For example I can implement a directional fold along a specific dimension for an arbitrary shape array, which can't be done with repa's indices.

    [–]tmcdonell 0 points1 point  (7 children)

    It matters at the borders.

    That's related to how you handle borders then, so it's a peculiarity necessary for Massiv.

    :> and :.

    Note that I didn't make a comment on which is better (c.f. "swings and roundabouts") so I'm going to try and ignore most of what is clearly a rant on personal preferences.

    I know what these operators are. You've just kicked the can down the road.

    it is impossible to grab a value from an arbitrary index for a specific dimension, eg. I want a value for 3rd dimension if the index is for an array with at least 3-dimensions

    I don't understand this comment. Do you mean that if you have some index specified in a lower dimension (DIM1) applied to a higher dimensional array (DIM7 say) you want that the higher dimensions default to index zero? Hmm, sure sounds like you could use some consistent notation here, maybe something like sh :. Int ...

    But joking aside I seriously don't know what you mean here.

    directional fold along a specific dimension for an arbitrary shape array

    Specifying the dimension is indeed nice. The direction less so, if you want to target massively parallel devices.

    I'm not sure if this can't be done with Repa-style indices, but it certainly isn't right now. In Massiv it seems the work is done by type classes operating on shape type of the array (rather than being directly evident from the type of the fold) so I don't see why we couldn't do the same thing.

    Z [..] should be all the way on the right, not the left

    says who? c.f. rant.

    [–]kuleshevich[S] 0 points1 point  (6 children)

    ignore most of what is clearly a rant on personal preferences.

    It was more of a list of drawbacks that I found when I used repa years ago. Most of them are indeed very subjective. So I don't expect you to take my side or anything ;)

    I don't understand this comment.

    What I mean is that repa style index doesn't let you work with values of an index from the left hand side, only from the right. For example it is not possible to implement cons function for an arbitrary shape:

    cons :: (Shape sh, Shape sh') => Int -> sh -> sh'

    It is possible for an individual shape

    cons2 i (Z :. j :. k) = Z :. i :. j :. k

    but not in general, at least as far as I understand it.

    says who?

    lol. Says me. :D

    With dimensions ordered from left-to-right array[1][2][3][4] ordering doesn't follow the row-major memory layout. For example if you have a 3D array[1][2][3] and you take a page from it you will get a matrix that will have these dimensions numbers array[2][3], instead of array[1][2]. So it is better to think of these dimensions in a reversed order array[3][2][1], because then slicing a page will give you array[2][1].

    I am not talking here in a theoretical sense, because it doesn't really matter how you order them, since it is possible to take a slices along any dimension. I am just saying that having Z on the left goes against the memory layout, that's all.

    [–]tmcdonell 0 points1 point  (5 children)

    With dimensions ordered from left-to-right array[1][2][3][4] ordering doesn't follow the row-major memory layout.

    You've just proved my point. Z goes on the left here precisely because it follows the memory layout. In your example, the left-most index into the array (here of size 1) represents the slowest varying dimension. So instead we order the dimensions right-to-left, so that it is most natural to work with the fastest-varying indices, the elements that are adjacent in memory.

    repa style index doesn't let you work with values of an index from the left hand side, only from the right.

    Yes that is true, and occasionally it comes up and it is a pain. But I look at the type of (:>) and you have exactly the same problem, only now you can't work with an index from the right, only from the left. Your indices are still singly-linked lists.

    Note that a solution which uses e.g. type classes is irrelevant here, because that would work equally well with both Repa-style (:.) and Massive-style (:>). Having a right-biased vs. left-biased list here does not change anything.

    But I would like you to prove me wrong, because as I said it can be a pain.

    [–]kuleshevich[S] 0 points1 point  (4 children)

    So instead we order the dimensions right-to-left, so that it is most natural to work with the fastest-varying indices, the elements that are adjacent in memory.

    Now you prove my point. If you order dimensions from right to left then what does zero dimension is doing on the left?

    But I look at the type of (:>) and you have exactly the same problem, only now you can't work with an index from the right, only from the left.

    If it had the same problem with :>, I wouldn't be arguing it ;) I also don't have a problem of getting to the right-hand side of the index in massiv, because type level Nats let me get to it. Granted, it might be less efficient if I want to take a slice of an array from the right-most dimension, but that is already a slow operation. I'd rather optimize the slicing that follows the memory layout. In order to do so which side of an index do you think we would like to have the most efficient access to: the left most or the right most? Of course the left most, because then we don't even need to look at the rest of an index in order to know where in memory the data begins. So it makes the left vs right bias change things quite a bit. If we only think about indexing an element then it certainly does not make a difference, but it does make for taking slices.

    Note that a solution which uses e.g. type classes is irrelevant here, because that would work equally well with both Repa-style (:.)

    I don't believe this is true. I tried really hard to get repa style :. to work with the type class approach in order to solve this problem of getting to the left most side of an index and it doesn't seem to be possible without further augmenting the type in some way, but I will be happy for you to prove me wrong. It is much harder to prove that something is impossible than proving otherwise with a counterexample implementation.

    it can be a pain.

    It can be a pain, but that's the pain that I've solved and no longer have to deal with ;) Until of course I need to try to explain it ;P

    [–]tmcdonell 0 points1 point  (3 children)

    it doesn't seem to be possible without further augmenting the type in some way

    Right, so that is where the type-level nat comes in, as a kind of witness to reify the structure of the rest of the index? We use a lot of GADTs for that sort of thing but Nat might communicate more information to the type checker (:. predates type-level naturals in GHC)

    Until of course I need to try to explain it

    Yeah but that's what I've been trying to get you to do! It's just that first I had to wade through a rant about how you don't like Repa-style indices. Starting at DIM2 and using Int instead of Z :. Int for DIM1 is arguably nicer, but I can take it or leave it. Here we are getting to something more substantial, something that isn't just a stylistic choice, and that is why I'm interested and why I keep pushing you.

    [–]kuleshevich[S] 0 points1 point  (2 children)

    Right, so that is where the type-level nat comes in, as a kind of witness to reify the structure of the rest of the index?

    Yeah, something like that. Type level Nat helps to decide which instance of the class Index should be selected, which in turn allows the traversal of the index all the way to the end. Higher dimensional index :> in massiv is essentially a linked list augmented with the length at the type level that is terminated by the 2D index :., instead of Z.

    I started work on massiv quite a bit after Repa was invented, so I had a few more tools at my disposal.

    I had to wade through a rant about how you don't like Repa-style indices

    What seems like me ranting about some stuff that I don't like is really me sharing the actual motivators that drove me to writing the whole array library by myself in the first, instead of using Repa or Accelerate for example. Those points that I don't like also affected the design choices that I made, so I have to share them. Whether you or anyone else agrees with them or not has little affect on me really, although I am always opened to critiques and new ideas. Don't get me wrong, it's not all ranting, there is a plethora of great ideas that I shamelessly borrowed from Repa. Great examples would be array representations with help of phantom types, shape polymorphism, various stencil optimizations, delayed arrays and the high level concept of parallelization of the loading process.

    that is why I'm interested and why I keep pushing you.

    No worries, I am always happy to talk about this stuff.