# Julia 1.0, tight-binding benchmark and array slices

**URL:** <https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943>\
**Category:** Performance\
**Created:** [August 23, 2018, 2:34pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943 "2018-08-23T14:34:17Z")\
**Posts on this page:** 10\
**Page:** 1

<div class="post-metadata">

**Author:** ![jabl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jabl/32/5112_2.png) [@jabl](https://discourse.julialang.org/u/jabl)\
**Post date:** [August 23, 2018, 2:34pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/1 "2018-08-23T14:34:17Z")

</div>

Hi,

following the recent release of Julia 1.0 I updated a small benchmark tight-binding program that I have implemented in Fortran, C++ with Eigen, C++ with Armadillo, and python/numpy. Roughly, the Fortran and both C++ versions are equivalent both in terms of LOC and performance. The Julia and Numpy versions are roughly the same in terms of LOC, about half the LOC of the Fortran/C++ versions. The Numpy version, however, is very slow, roughly a factor of 50 slower than Fortran (excluding the part which is just a lapack call).

Now, previously in the Julia 0.4 timeframe, the Julia version was about half as fast as the Fortran/C++ versions. That version used the Devectorize package, which seems to have been unmaintained now for several years. I was unable to make it work with julia 0.6.x, not to mention 1.0. However, it seems that as of Julia 0.6 there is the “@.” macro which does roughly the same as the @devec macro from Devectorize(?). With @. for a few critical operations, Julia 1.0 is a factor of 1.7 slower than Fortran. Without @., about a factor of 2.1 slower.

However, if I rewrite those expression as manual loops, Julia is only a factor of 1.1 slower than Fortran, that is, more or less the same! Very impressive!

Although slightly disappointing that I had to resort to writing manual loops for performance. Is there some trick I’m missing? The expressions in question are all of the form

@. v[:] = atoms[bj,:] - atoms[bi,:]

which I rewrite as an explicit loop like:

for z = 1:3  
v[z] = atoms[bj,z] - atoms[bi,z]  
end

Does Julia create a copy as part of the slicing operation, or what makes the array syntax slow? The @time macro does report a lot of allocations due to this, whether it’s an actual copy, an array descriptor for the slice, or whatever. Allocations for one particular case:

- Unoptimized: 3.83 M
- Using @.: 2.55 M
- Explicit loops: 4

Is there anything that can be done here?

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [August 23, 2018, 2:36pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/2 "2018-08-23T14:36:29Z")

</div>

> [@jabl](#):
>
> The expressions in question are all of the form
> 
> @. v[:] = atoms[bj,:] - atoms[bi,:]

You want to do

```julia
@views @. v = atoms[bj,:] - atoms[bi,:]

```

But there is a performance regression on Julia v1.0 which does currently slow down broadcasting.

> <https://github.com/JuliaLang/julia/issues/28126>
>
> Here is a minimal working example.
> 
> \`\`\`julia
> julia\> using BenchmarkTools
> 
> j…ulia\> function foo(a::Vector{T}, b::Vector{T}, c::Vector{T}, d::Vector{T}, e::Vector{T}) where T
> @. a = b + 0.1 \* (0.2c + 0.3d + 0.4e)
> nothing
> end
> foo (generic function with 1 method)
> 
> julia\> function goo(a::Vector{T}, b::Vector{T}, c::Vector{T}, d::Vector{T}, e::Vector{T}) where T
> @assert length(a) == length(b) == length(c) == length(d) == length(e)
> @inbounds for i in eachindex(a)
> a\[i\] = b\[i\] + 0.1 \* (0.2c\[i\] + 0.3d\[i\] + 0.4e\[i\])
> end
> nothing
> end
> goo (generic function with 1 method)
> 
> julia\> a,b,c,d,e=(rand(1000) for i in 1:5)
> Base.Generator{UnitRange{Int64},getfield(Main, Symbol("##9#10"))}(getfield(Main, Symbol("##9#10"))(), 1:5)
> 
> julia\> @btime foo($a,$b,$c,$d,$e)
> 1.277 μs (0 allocations: 0 bytes)
> 
> julia\> @btime goo($a,$b,$c,$d,$e)
> 345.568 ns (0 allocations: 0 bytes)
> 
> julia\> versioninfo()
> Julia Version 0.7.0-beta2.12
> Commit a878341 (2018-07-15 15:57 UTC)
> Platform Info:
> OS: Linux (x86\_64-pc-linux-gnu)
> CPU: Intel(R) Core(TM) i7-6820HQ CPU @ 2.70GHz
> WORD\_SIZE: 64
> LIBM: libopenlibm
> LLVM: libLLVM-6.0.0 (ORCJIT, skylake)
> Environment:
> JULIA\_PKG3\_PRECOMPILE = 1
> \`\`\`

(Note: Don’t forget to `@inbounds` or `@inbounds @simd` that `for` loop in the benchmarks 🙂)

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [August 23, 2018, 2:37pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/3 "2018-08-23T14:37:06Z")

</div>

That sounds like a great benchmark.

> [@jabl](#):
>
> Does Julia create a copy as part of the slicing operation

Yes, that’s precisely what happens — slicing isn’t able to be “devectorized” like all other function calls. You can use `@views` alongside the `@.` macro to instead make slicing return a lazy view instead of a copy.

---

<div class="post-metadata">

**Author:** ![carstenbauer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/carstenbauer/32/4981_2.png) [@carstenbauer](https://discourse.julialang.org/u/carstenbauer)\
**Post date:** [August 23, 2018, 2:46pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/4 "2018-08-23T14:46:25Z")

</div>

Also, try `julia -O3` and `@simd`.

---

<div class="post-metadata">

**Author:** ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)\
**Post date:** [August 23, 2018, 3:40pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/5 "2018-08-23T15:40:23Z")

</div>

You are also using an access pattern optimal for row major storage, whereas Julia uses column major storage. Flip the array dimensions and you might see increased performance.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [August 23, 2018, 6:18pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/6 "2018-08-23T18:18:18Z")

</div>

If the second dimension is small and you always use it together, not only would it be better to transpose it, but using an Array of SVectors (from StaticArrays.jl) would likely help too.

---

<div class="post-metadata">

**Author:** ![jabl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jabl/32/5112_2.png) [@jabl](https://discourse.julialang.org/u/jabl)\
**Post date:** [August 23, 2018, 6:20pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/7 "2018-08-23T18:20:01Z")

</div>

Thanks for all the suggestions. A little more experimenting showed that for this particular case:

- @views helps a bit
- -O3 doesn’t seem to have any effect
- transposing the atoms array did help a little bit. Surprisingly little, but my largest atoms array is about 10 kB, so it all fits in L1 cache anyhow. I guess the biggest benefit might be to enable SIMD, but OTOH with only 3 elements the benefits of SIMD are, well, minute.

All in all, with transposed atoms array I got:

- Explicit loops + @inbounds @simd: 0.96 x Fortran
- Slices with @views @.: 1.34 x Fortran

Pretty nice!

---

<div class="post-metadata">

**Author:** ![kristoffer.carlsson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/kristoffer.carlsson/32/22_2.png) [@kristoffer.carlsson](https://discourse.julialang.org/u/kristoffer.carlsson)\
**Post date:** [August 23, 2018, 7:06pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/8 "2018-08-23T19:06:12Z")

</div>

Being able to see the code would be nice.

---

<div class="post-metadata">

**Author:** ![jabl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jabl/32/5112_2.png) [@jabl](https://discourse.julialang.org/u/jabl)\
**Post date:** [August 24, 2018, 6:19am UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/9 "2018-08-24T06:19:29Z")

</div>

I’d love to, but the code was originally a homework exercise for a course, and AFAIK they are still giving this course, so it’d be a bit bad style to give out the exercise answer. Come to think of it, I should ask if they are still using this same exercise, if not I guess there’s nothing to prevent releasing it.

---

<div class="post-metadata">

**Author:** ![jabl](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jabl/32/5112_2.png) [@jabl](https://discourse.julialang.org/u/jabl)\
**Post date:** [September 22, 2018, 7:02pm UTC](https://discourse.julialang.org/t/julia-1-0-tight-binding-benchmark-and-array-slices/13943/10 "2018-09-22T19:02:37Z")

</div>

Well, to wake up this semi-zombie thread, I asked and got permission for releasing the code, wooo! So here it is: [Janne Blomqvist / tb · GitLab](https://gitlab.com/jabl/tb)

Please let me know if you have issues running it, or any other feedback for that matter!