4 ms·

The real performance win seems to be reducing the number of allocations (2 vs. 30 thousand).

doublesCs 6y ago
That's not uncommon in HPC. I'm surprised if this was inherent to StaticArrays.jl though. Probably some interop edge case, I'm imagining. Shall we summon Chris Rackauckas?
- ChrisRackauckas 6y ago
  It worked, I answered on Discourse because someone pinged me for an answer, saying "you've been summoned" haha.
marmaduke 6y ago
So the real question is why a library with the word "static" in its title does 30 thousand allocations..
- dnautics 6y ago
  I believe the "static" refers to static, aka, compile-time dimension of the array (yes, julia is a compiling language), versus an array that has a run-time, mutable array dimensionality.
- ChrisFoster 6y ago
  "Static" is a reference to the size of the arrays being fixed and encoded in the type. There's various allocation strategies, for instance StaticArrays.SizedArray uses a normal heap-allocated Array as the storage. That said, the backing storage for the particular static array types used in this example (SMatrix,SVector) is an immutable Tuple which should generally be stack allocated if the library and optimizer are working together as intended. So the allocations here are a bit of a surprise; they may indicate a lack of inlining in something we expected to be inlined.
- tkoolen 6y ago
  1) a mistake in the benchmark procedure explained in https://discourse.julialang.org/t/grassmann-jl-a-b-3x-faster-than-julias-staticarrays-jl/41451/2 https://discourse.julialang.org/t/grassmann-jl-a-b-3x-faster..., and 2) the dynamically-sized result vector containing 10000 statically-sized 5-vectors still needs to be allocated.
spacedome 6y ago
This has been my general experience writing numerical linear algebra methods in Julia. Optimizing a specific algorithm, eliminating unnecessary allocations is the first thing I do, and can give large performance gains, especially for iterative methods. It can be a bit annoying, and can result in much less readable code (eg having to explicitly write things like mul!(C, A, B)). It ends up looking a lot more like the FORTRAN it is meant to replace, and if you aren't careful you lose the ability to use generic types, though multiple dispatch is a great solution. The worst case I have found is iteratively calling linear algebra solvers (Eig, SVD, LU) which do not have any option for preallocating and reusing work arrays. The only way to do this now is to call the BLAS/LAPACK routines directly with ccall, which is a huge pain. There was an attempt to fix this with PowerLAPACK.jl, but it seems abandoned, so for the moment writing optimized methods in FORTRAN/C and calling from Julia is sometimes still preferable. Julia does quite well though, so this only seems necessary for things like core NLA libraries, for example I would not try rewriting ARPACK in Julia for a performance gain. The gain in flexibility for writing these in Julia is absolutely huge though, and I definitely recommend it. The Grassmann.jl library has many examples of what a language like Julia makes possible, or DifferentialEquations.jl.
- Xcelerate 6y ago
  I love Julia and have used it for years, but this is really my only complaint with my language. I used to do HPC work and ended up spending a lot of time tracking down where tiny bits of memory were being allocated in for loops. A lot of my routines ended up having “blocks” of memory to pass in as arguments, which as you mention, moves the language a bit closer to C or Fortran IMO.
  - snicker7 6y ago
    Minimizing allocations (and deallocations) and optimizing memory layout are just things you have to do for HPC in any language. Though I think that it is easier in Julia than in C/Fortran.
- etik 6y ago
  There is a recent effort [1] to provide low-level support for faster operations by transforming user code to take advantage of a compiler's instruction set, memory packing, etc. This is being expanded upon to essentially provide a Julia-native BLAS. Some of the benchmarks are even competitive with or beat Intel MKL (calibrate that statement appropriately to your level of trust in benchmarks). I wouldn't count out a Julia ARPACK implementation just yet. [1] LoopVectorization: https://github.com/chriselrod/LoopVectorization.jl https://github.com/chriselrod/LoopVectorization.jl Announcement post and discussion: https://discourse.julialang.org/t/ann-loopvectorization/32843 https://discourse.julialang.org/t/ann-loopvectorization/3284...

That's not uncommon in HPC. I'm surprised if this was inherent to StaticArrays.jl though. Probably some interop edge case, I'm imagining. Shall we summon Chris Rackauckas?

ChrisRackauckas 6y ago
It worked, I answered on Discourse because someone pinged me for an answer, saying "you've been summoned" haha.

marmaduke 6y ago

So the real question is why a library with the word "static" in its title does 30 thousand allocations..

dnautics 6y ago
I believe the "static" refers to static, aka, compile-time dimension of the array (yes, julia is a compiling language), versus an array that has a run-time, mutable array dimensionality.
ChrisFoster 6y ago
"Static" is a reference to the size of the arrays being fixed and encoded in the type. There's various allocation strategies, for instance StaticArrays.SizedArray uses a normal heap-allocated Array as the storage. That said, the backing storage for the particular static array types used in this example (SMatrix,SVector) is an immutable Tuple which should generally be stack allocated if the library and optimizer are working together as intended. So the allocations here are a bit of a surprise; they may indicate a lack of inlining in something we expected to be inlined.
tkoolen 6y ago
1) a mistake in the benchmark procedure explained in https://discourse.julialang.org/t/grassmann-jl-a-b-3x-faster-than-julias-staticarrays-jl/41451/2 https://discourse.julialang.org/t/grassmann-jl-a-b-3x-faster..., and 2) the dynamically-sized result vector containing 10000 statically-sized 5-vectors still needs to be allocated.

spacedome 6y ago

This has been my general experience writing numerical linear algebra methods in Julia. Optimizing a specific algorithm, eliminating unnecessary allocations is the first thing I do, and can give large performance gains, especially for iterative methods. It can be a bit annoying, and can result in much less readable code (eg having to explicitly write things like mul!(C, A, B)). It ends up looking a lot more like the FORTRAN it is meant to replace, and if you aren't careful you lose the ability to use generic types, though multiple dispatch is a great solution. The worst case I have found is iteratively calling linear algebra solvers (Eig, SVD, LU) which do not have any option for preallocating and reusing work arrays. The only way to do this now is to call the BLAS/LAPACK routines directly with ccall, which is a huge pain. There was an attempt to fix this with PowerLAPACK.jl, but it seems abandoned, so for the moment writing optimized methods in FORTRAN/C and calling from Julia is sometimes still preferable. Julia does quite well though, so this only seems necessary for things like core NLA libraries, for example I would not try rewriting ARPACK in Julia for a performance gain. The gain in flexibility for writing these in Julia is absolutely huge though, and I definitely recommend it. The Grassmann.jl library has many examples of what a language like Julia makes possible, or DifferentialEquations.jl.

Xcelerate 6y ago
I love Julia and have used it for years, but this is really my only complaint with my language. I used to do HPC work and ended up spending a lot of time tracking down where tiny bits of memory were being allocated in for loops. A lot of my routines ended up having “blocks” of memory to pass in as arguments, which as you mention, moves the language a bit closer to C or Fortran IMO.
- snicker7 6y ago
  Minimizing allocations (and deallocations) and optimizing memory layout are just things you have to do for HPC in any language. Though I think that it is easier in Julia than in C/Fortran.
etik 6y ago
There is a recent effort [1] to provide low-level support for faster operations by transforming user code to take advantage of a compiler's instruction set, memory packing, etc. This is being expanded upon to essentially provide a Julia-native BLAS. Some of the benchmarks are even competitive with or beat Intel MKL (calibrate that statement appropriately to your level of trust in benchmarks). I wouldn't count out a Julia ARPACK implementation just yet. [1] LoopVectorization: https://github.com/chriselrod/LoopVectorization.jl https://github.com/chriselrod/LoopVectorization.jl Announcement post and discussion: https://discourse.julialang.org/t/ann-loopvectorization/32843 https://discourse.julialang.org/t/ann-loopvectorization/3284...