# Any suggestions for speeding this simple scalar function up?

**URL:** <https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748>\
**Category:** Performance\
**Tags:** performance\
**Created:** [August 20, 2021, 8:19pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748 "2021-08-20T20:19:06Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 20, 2021, 8:19pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/1 "2021-08-20T20:19:06Z")

</div>

Hi!

I’ve been playing around calculating this expression:

```julia
P = b * ((rho/rho0)^gamma - 1)

```

Where values for `b` are shown how to be calculated in the code (see bottom spoiler). Running the different variants I have of the code, the minimum times I get:

```julia
Minimum Time EOS1:
Trial(62.999 μs)
Minimum Time EOS2:
Trial(68.600 μs)
Minimum Time EOS3:
Trial(69.401 μs)
Minimum Time EOS1AVX:
Trial(61.799 μs)
Minimum Time EOS1SIMILAR:
Trial(63.399 μs)

```

Where `EOS1` which is written in the most “naive” way, seems to perform the best. For fun and curiosity, would anyone know how to further speed this up?

Thanks for reading.

> **Code**
>
> ```julia
> using BenchmarkTools
> using LoopVectorization
> rho = rand(990:rand():1010,10^3)
> 
> function EOS1(rho;c0=100,gamma=7,rho0=1000)
> b=(c0^2*rho0)/gamma
> 
> P = b*((rho/rho0).^gamma .- 1)
> end
> 
> function EOS1SIMILAR(rho;c0=100,gamma=7,rho0=1000)
> b=(c0^2*rho0)/gamma
> 
> P = similar(rho)
> 
> P .= b*((rho/rho0).^gamma .- 1)
> end
> 
> function EOS2(rho;c0=100,gamma=7,rho0=1000)
> b=(c0^2*rho0)/gamma
> 
> rho0_denominator = 1/(rho0^gamma);
> 
> P = b*(rho.^gamma * rho0_denominator .- 1)
> end
> 
> function EOS3(rho;c0=100,gamma=7,rho0=1000)
> b=(c0^2*rho0)/gamma
> 
> rho0_denominator = 1/(rho0^gamma);
> 
> P = rho.^gamma * rho0_denominator * b .- b
> end
> 
> bEOS1 = @benchmark EOS1($rho)
> println("Minimum Time EOS1:\n",bEOS1)
> bEOS2 = @benchmark EOS2($rho)
> println("Minimum Time EOS2:\n",bEOS2)
> bEOS3 = @benchmark EOS3($rho)
> println("Minimum Time EOS3:\n",bEOS3)
> 
> bEOS1AVX = @benchmark @avx EOS1($rho)
> println("Minimum Time EOS1AVX:\n",bEOS1AVX)
> bEOS1SIMILAR = @benchmark EOS1SIMILAR($rho)
> println("Minimum Time EOS1SIMILAR:\n",bEOS1SIMILAR)
> 
> ```

---

<div class="post-metadata">

**Author:** ![Sukera](https://avatars.discourse-cdn.com/v4/letter/s/ce7236/32.png) [@Sukera](https://discourse.julialang.org/u/Sukera)\
**Post date:** [August 20, 2021, 8:24pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/2 "2021-08-20T20:24:40Z")

</div>

Use more dots. In all functions, only `rho` seems to be a vector, while all other objects are scalars. Dots in the same expression are only fused into a single loop if all operations that are touched are also dotted. If they’re seperated, the broadcast will have to materialize an array first by allocating and creating more loops than necessary.

Another way would be to use `@.`.

E.g.

```julia
function EOS4(rho; c0=100, gamma=7, rho0=1000)
    b=(c0^2*rho0)/gamma

    rho0_denominator = 1/(rho0^gamma);

    return rho .^ gamma .* rho0_denominator .* b .- b
end

```

---

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 20, 2021, 8:31pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/3 "2021-08-20T20:31:05Z")

</div>

Thanks for the suggestion.

When I try on my system with your function:

```julia
bEOS4 = @benchmark EOS4($rho)
BenchmarkTools.Trial: 
  memory estimate: 7.94 KiB
  allocs estimate: 1
  --------------
  minimum time: 67.399 μs (0.00% GC)
  median time: 69.799 μs (0.00% GC)
  mean time: 70.641 μs (0.00% GC)
  maximum time: 130.999 μs (0.00% GC)
  --------------
  samples: 10000
  evals/sample: 1

```

Which is still ~10% slower than `EOS1`. I tried adding more “.” in `EOS1` too, but did not affect performance.

Kind regards

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [August 20, 2021, 9:04pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/4 "2021-08-20T21:04:48Z")

</div>

You seem to have tried some `avx` variant, but is it there? Here `@turbo` makes a _huge_ difference:

```julia
julia> using LoopVectorization

julia> function EOS1(rho;c0=100,gamma=7,rho0=1000)
           b=(c0^2*rho0)/gamma
           P = b*((rho/rho0).^gamma .- 1)
       end
EOS1 (generic function with 1 method)

julia> function EOS1_turbo(rho;c0=100,gamma=7,rho0=1000)
           b=(c0^2*rho0)/gamma
           P = @turbo b*((rho/rho0).^gamma .- 1)
       end
EOS1_turbo (generic function with 1 method)

julia> EOS1(rho) ≈ EOS1_turbo(rho)
true

julia> @btime EOS1($rho);
  68.043 μs (3 allocations: 23.81 KiB)

julia> @btime EOS1_turbo($rho);
  1.809 μs (3 allocations: 23.81 KiB)

```

If this will be called from within a hot loop, you probably want to preallocate `b`.

edit:

If you use `@.` to have better chances of not forgetting some loop fusion, that gets even better:

```julia
julia> function EOS1_turbo(rho;c0=100,gamma=7,rho0=1000)
           b=(c0^2*rho0)/gamma
           P = @turbo @. b*((rho/rho0)^gamma - 1)
       end
EOS1_turbo (generic function with 1 method)

julia> @btime EOS1_turbo($rho);
  701.667 ns (1 allocation: 7.94 KiB)

```

(the reason one would notice that that was needed is that there were 3 allocations in the previous version, where only `P` should be a new allocation there. Thus, some `.` is missing and an intermediate vector is being generated in the previous version, which should be: `P = @turbo b .* ((rho ./ rho0).^gamma .- 1)`.

---

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 20, 2021, 9:17pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/5 "2021-08-20T21:17:20Z")

</div>

> [@lmiq](#):
>
> ```julia
> function EOS1_turbo(rho;c0=100,gamma=7,rho0=1000)
> b=(c0^2*rho0)/gamma
> P = @turbo @. b*((rho/rho0)^gamma - 1)
> end
> 
> ```

Wow, you are indeed right, it has a huge impact! I also agree on using `@.`, a lot less thinking of putting dots and gives better looking equation.

But some questions:

What is a hot loop?

A simple explanation of what @turbo does to make this much faster?

If you have time and want to.

Thanks!

---

<div class="post-metadata">

**Author:** ![stillyslalom](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stillyslalom/32/45687_2.png) [@stillyslalom](https://discourse.julialang.org/u/stillyslalom)\
**Post date:** [August 20, 2021, 9:20pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/6 "2021-08-20T21:20:49Z")

</div>

`@turbo` has a lot of tricks up its sleeve, but faster integer exponentiation makes the biggest difference:

```julia
julia> pow(x, n) = x^n
pow (generic function with 1 method)

julia> @btime Base.map(x -> pow(x, 7), x) setup=(x=rand(1000));
  30.900 μs (1 allocation: 7.94 KiB)

julia> @btime LoopVectorization.vmap(x -> pow(x, 7), x) setup=(x=rand(1000));
  897.872 ns (1 allocation: 7.94 KiB)

```

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [August 20, 2021, 9:22pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/7 "2021-08-20T21:22:42Z")

</div>

> [@Ahmed\_Salih](#):
>
> What is a hot loop?

If you have to call that same function thousand of times inside a loop, you will like to avoid creating a new `P` every time. Thus, create a variant of the function that receives a preallocated P and updates it.

> [@Ahmed\_Salih](#):
>
> A simple explanation of what @turbo does to make this much faster?

`@turbo` Loop-vectorizes the operation, meaning that uses the capacity of the processors to perform more than one operation at once (search for SIMD or loop vectorization in google). When the operations are amenable to that, it is very good (sometimes that happens automatically, even without the `@turbo` flag, and sometimes you get a good speedup with `@simd`).

Now, how exactly in this specific case the speedup is that great, I think some expert should kick in. Normally I would expect something of the order of 4x to 8x speedup.

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [August 20, 2021, 9:23pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/8 "2021-08-20T21:23:37Z")

</div>

> [@stillyslalom](#):
>
> faster integer exponentiation makes the biggest difference

Ah, well noticed. Then you should probably take a look at: [GitHub - JuliaMath/FastPow.jl: optimal addition-chain exponentiation for Julia](https://github.com/JuliaMath/FastPow.jl)

(for some reason I don’t know, it is not affecting the result here:)

```julia
julia> function EOS1_fastpow(rho;c0=100,gamma=7,rho0=1000)
           b=(c0^2*rho0)/gamma
           @fastpow P = @. b*((rho/rho0)^gamma - 1)
       end
EOS1_fastpow (generic function with 1 method)

julia> @btime EOS1_fastpow($rho);
  67.417 μs (1 allocation: 7.94 KiB)

```

(tried moving the macros to see if the result changes, with no success)

---

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 20, 2021, 9:37pm UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/9 "2021-08-20T21:37:40Z")

</div>

Thanks for the observation @stillyslalom !

I also tried using FastPow and pow(x,n), but unfortunately the first one I could not get to work as you mentioned Leandro and the other one made it slower with @turbo @.

Thanks for the explanation of hot loop too !

Kind regards

---

<div class="post-metadata">

**Author:** ![CRquantum](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/crquantum/32/27824_2.png) [@CRquantum](https://discourse.julialang.org/u/CRquantum)\
**Post date:** [August 21, 2021, 8:18am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/10 "2021-08-21T08:18:59Z")

</div>

Sorry for interrupting a stupid question, if I do not use @ turbo trick for those .+ .- .\* ./ .^2 etc stuff, will the speed be slow down a lot?  
if I do not use these tricks, by default, what is Julia’s optimization level? like equivalent to Fortran’s O2 or something?

---

<div class="post-metadata">

**Author:** ![Sukera](https://avatars.discourse-cdn.com/v4/letter/s/ce7236/32.png) [@Sukera](https://discourse.julialang.org/u/Sukera)\
**Post date:** [August 21, 2021, 9:12am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/11 "2021-08-21T09:12:25Z")

</div>

Julia and LLVM both do loop unrolling and SIMD vectorization to a certain extent, but `@turbo` by default takes more knowledge about the current CPU into account. It’s also much more aggressive and (with wrong use) can be unsafe. How much `@turbo` will help in that regard is not really possible to pin down to a specific number (it depends on the loop in question) - it’s just that through its aggressive nature, you can usually see a large speedup compared to base julia and its LLVM passes (as far as I know, mostly because `@turbo` requires independent loop iterations, whereas both LLVM and julia would have to prove that on their own and can’t just assume it).

While julia itself does have “optimization levels” as a command line flag, they pretty much don’t do anything nowadays as far as I’m aware.

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [August 21, 2021, 9:13am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/12 "2021-08-21T09:13:07Z")

</div>

You can shave off a bit of extra time (5%-15%) by removing one division:

```julia
function EOS2_turbo(ρ; c0=100, γ=7, ρ0=1000)
    b1 = (c0^2 * ρ0) / γ
    b2 = 1 / ρ0^γ
    P = @turbo b2 .* ρ.^γ .- b1
    return P
end

```

```julia
jl> @btime EOS1_turbo($ρ);
  638.608 ns (1 allocation: 7.94 KiB)

jl> @btime EOS2_turbo($ρ);
  601.775 ns (1 allocation: 7.94 KiB)

```

(Sorry about changing your symbols, but I can no longer stand seeing variable names like `gamma` 😉 )

**Edit:** There’s an error in `EOS2_turbo`, but I can’t find what it is. The values are extremely wrong, though.

**Update:** My suggestion is apparently numerically abhorrent:

```julia
jl> ρ = 990.1
990.1

jl> ρ^γ / ρ0^γ
240.65219205640526

jl> (ρ/ρ0)^γ
0.9327245837531053

```

So _don’t use it_!

---

<div class="post-metadata">

**Author:** ![giordano](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/giordano/32/2166_2.png) [@giordano](https://discourse.julialang.org/u/giordano)\
**Post date:** [August 21, 2021, 9:42am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/13 "2021-08-21T09:42:07Z")

</div>

What if gamma and rho0 are floats instead of integers that can overflow?

---

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 21, 2021, 9:44am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/14 "2021-08-21T09:44:11Z")

</div>

Haha no problem - I also like to use unicode symbols, just not for posting questions on Discourse generally 😅

Interesting speed up, but weird behaviour. It seems like the division and ^ takes the same precendence, so it goes from left to right:

```julia
julia> 990.1^7 / 1000^7
240.65219205640526

julia> (990.1/1000)^7
0.9327245837531053

julia> (990.1^7) / (1000^7)
240.65219205640526

```

So the error is that it is calculating:

(ρ^γ / ρ0)^γ

But that gives:

```julia
(990.1^7/1000)^7
6.14149305499489e125

```

So I am a bit my confused self now too 🙂

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [August 21, 2021, 9:45am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/15 "2021-08-21T09:45:50Z")

</div>

> [@CRquantum](#):
>
> Sorry for interrupting a stupid question, if I do not use @ turbo trick for those .+ .- .\* ./ .^2 etc stuff, will the speed be slow down a lot?  
> if I do not use these tricks, by default, what is Julia’s optimization level? like equivalent to Fortran’s O2 or something?

In this particular case I would bet that turbo is much more agressive than any Fortran compiler. You can try it. AFAIK you can get more agressive loop optimizations with --march=native in Fortran. (you have asked that in the Fortran forum, by the way)

---

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 21, 2021, 9:56am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/16 "2021-08-21T09:56:08Z")

</div>

rho0 can be a float value, I have changed it to 1000.0

gamma can too in theory, but most often in my cases it would just be “7”. I found using “7” instead of “7.0” sped up the power calculation.

Kind regards

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [August 21, 2021, 10:32am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/17 "2021-08-21T10:32:56Z")

</div>

It’s not about precedence, but a numerical issue. It works ok up to 6, and fails for 7.

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [August 21, 2021, 10:35am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/18 "2021-08-21T10:35:49Z")

</div>

Probably integer overflow for 1000^7

---

<div class="post-metadata">

**Author:** ![Sukera](https://avatars.discourse-cdn.com/v4/letter/s/ce7236/32.png) [@Sukera](https://discourse.julialang.org/u/Sukera)\
**Post date:** [August 21, 2021, 10:39am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/20 "2021-08-21T10:39:32Z")

</div>

Indeed:

```julia
julia> Base.Checked.checked_mul(1000^6, 1000)                             
ERROR: OverflowError: 1000000000000000000 * 1000 overflowed for type Int64
Stacktrace:                                                               
 [1] throw_overflowerr_binaryop(op::Symbol, x::Int64, y::Int64)           
   @ Base.Checked ./checked.jl:154                                        
 [2] checked_mul(x::Int64, y::Int64)                                      
   @ Base.Checked ./checked.jl:288                                        
 [3] top-level scope                                                      
   @ REPL[14]:1                                                           

```

---

<div class="post-metadata">

**Author:** ![Ahmed\_Salih](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ahmed_salih/32/206579_2.png) [@Ahmed\_Salih](https://discourse.julialang.org/u/Ahmed_Salih)\
**Post date:** [August 21, 2021, 10:41am UTC](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748/21 "2021-08-21T10:41:04Z")

</div>

> [@Sukera](#):
>
> `Base.Checked.checked_mul(1000^6, 1000)`

And @DNF ,does that mean I should always use floats for this? I.e. 1000.0^7?

Or how to interpret this result.

EDIT: Actually it is not as such a problem for my equation, since it just does (rho/rho0)^7, so no overflow. But it means I cannot get the speed up from calculating 1/rho0^7 first I suppose.

Kind regards

[Next page](https://discourse.julialang.org/t/any-suggestions-for-speeding-this-simple-scalar-function-up/66748.md?page=2)