# Usage of lazy vector instead of function

**URL:** <https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766>\
**Category:** Performance\
**Created:** [November 4, 2022, 9:11am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766 "2022-11-04T09:11:59Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![Salmon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/salmon/32/22968_2.png) [@Salmon](https://discourse.julialang.org/u/Salmon)\
**Post date:** [November 4, 2022, 9:11am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/1 "2022-11-04T09:11:59Z")

</div>

Hi everyone,  
recently I stumbled across a nice benchmark between different languages and was having a look at the Julia implementation:

> <https://github.com/niklas-heer/speed-comparison/blob/master/src/leibniz.jl>

Now in general, I understand why there is a use case for lazy arrays, but is this really one?  
As I understand it, this code should essentially be lowered to the following:

```julia
function f(rounds)
    pi = 1.0

    @simd for i in 2:(rounds + 2)
        x = (-1)^iseven(i)
        pi += x / (2 * i - 1)
    end

    return pi*4
end

```

This vector was introduced in a commit on purpose so there is probably a reason for it.  
Is `@simd` doing more magic if it is seeing a vector instead?  
I would be interested if there can be a performance gain by using a lazy array instead of essentially calling a function.

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 4, 2022, 9:28am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/2 "2022-11-04T09:28:08Z")

</div>

I do not know about lazy stuff, but I am pretty sure that `@simd` or `@turbo` would largely prefer the loop written as

```julia
function f2(rounds)
    pi = 1.0

	ub = floor(rounds/2)+1

    @simd for i in 1:ub
        pi += -2 / (16i*i-1)
    end
	if iseven(rounds)
		pi -= 1 / (4ub+3)
	end
    return pi*4
end

```

Where I took out the factor 4 and dealt with the parity of `rounds` separately. Less things inside the loop is always better: here there is only one addition, two multiplicatins and one division, while in your two versions, there was fetching values and/or powers.

If you are allowed to break a little the API, then i would add one to `ub` and remove the `if` block completely.

Edit: I might have the signs wrong, I’m on my phone 🙂

---

<div class="post-metadata">

**Author:** ![mcabbott](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcabbott/32/6603_2.png) [@mcabbott](https://discourse.julialang.org/u/mcabbott)\
**Post date:** [November 4, 2022, 10:00am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/3 "2022-11-04T10:00:39Z")

</div>

There is no difference, in fact the version here came before the one with the lazy array. The only question was whether computing `x` inside the loop was legal under the implicit rules of the game; IIRC no other solutions did this, all updated `x = -1.0 * x` in the loop, but the R solution did pre-compute an array of `x` values. This Julia version was 2x faster than C etc, because it allows for SIMD.

Afterwards, many others were updated to compute `x` inside the loop, catching up with Julia again. Julia could be too, to exactly this code.

And some were updated to re-order the sum. In  
[this PR](https://github.com/niklas-heer/speed-comparison/pull/57) @giordano proposes just `f(n) = sum(i -> 4/(4i-2n-3), 1:n)`.

So far the prize for the most devilish answer (not mine!) goes to this even faster version:

```julia
function f(rounds)
     рi = 1.0
     @simd for i in 2:(rounds + 2)
         x = (-1)^iseven(i)
         рi += x / (2 * i - 1)
     end
     рi *= 4
     return float(pi)
 end

```

@lrnv this seems to give the wrong answer, `f(10^6) - pi ≈ 1.367f0`. And is also slower.

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 4, 2022, 10:07am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/4 "2022-11-04T10:07:19Z")

</div>

Yes I saw, my bad. I corrected it with a correct answer, whihc is indeed much slower.

---

<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:** [November 4, 2022, 10:13am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/5 "2022-11-04T10:13:37Z")

</div>

> [@mcabbott](#):
>
> So far the prize for the most devilish answer (not mine!) goes to this even faster version:
> 
> ```julia
> function f(rounds)
> рi = 1.0
> @simd for i in 2:(rounds + 2)
> x = (-1)^iseven(i)
> рi += x / (2 * i - 1)
> end
> рi *= 4
> return float(pi)
> end
> 
> ```

For the curious, this is the [benchmark](https://github.com/giordano/speed-comparison/commit/ca3f9936206242bd36fa1ee4eaa27a075e5940b6#commitcomment-88008220):

 ![image](https://global.discourse-cdn.com/julialang/original/3X/d/9/d90069a109d11bdaa101bea7a007cd43c56bcb86.png)

---

<div class="post-metadata">

**Author:** ![Salmon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/salmon/32/22968_2.png) [@Salmon](https://discourse.julialang.org/u/Salmon)\
**Post date:** [November 4, 2022, 10:45am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/6 "2022-11-04T10:45:34Z")

</div>

Fantastic, not only does it run faster than all others, but it also has O(1) complexity and is by far more accurate than the others even for a small number of rounds. Truly, a good algorithm can make all the difference 😛 .

Jokes aside, trying out a bit, I think a faster way is not to rely on `iseven`:

```julia
function f(rounds)
    pi = 1.0

    @simd for i in 2:2:(rounds + 2)
        pi += -1. / (2 * i - 1)
        pi += 1. / (2 * i + 1)
    end

    return pi*4
end

```

But I presume this is against the rules, since I am doing some reordering and dont end up with the (numerically) same result. Anyway, it was fun 🙂

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 4, 2022, 11:20am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/7 "2022-11-04T11:20:16Z")

</div>

> [@Salmon](#):
>
> Jokes aside, trying out a bit, I think a faster way is not to rely on `iseven`:

On my machine, the runtime is exactly the same. Note also that simd IS reordering stuffs anyway

---

<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:** [November 4, 2022, 12:10pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/8 "2022-11-04T12:10:09Z")

</div>

> [@Salmon](#):
>
> But I presume this is against the rules

It isn’t like there are well defined rules and the maintainer of the repo has been pretty open to manual loop unrolling so far.

---

<div class="post-metadata">

**Author:** ![Salmon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/salmon/32/22968_2.png) [@Salmon](https://discourse.julialang.org/u/Salmon)\
**Post date:** [November 4, 2022, 12:11pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/9 "2022-11-04T12:11:03Z")

</div>

Interesting. Ive only tried it on my laptop (for which simd doesnt seem to be so great, admittedly) and here its quite a bit faster.

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [November 4, 2022, 5:13pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/10 "2022-11-04T17:13:12Z")

</div>

This might be sorta not fair, but combining two divisions into one cuts time in half.

```julia
function f2(rounds)
    p = 1.0
    @simd for i in 2:2:(rounds + 2)
        p += -2.0 / (4*i^2 -1)
    end
    return p*4
end

```

Is it possible for Julia to understand this tranformation natively?

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 4, 2022, 5:49pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/11 "2022-11-04T17:49:05Z")

</div>

With one less multiplication per iter :

```julia
function f3(rounds)
    p = 2.0
    @simd for i in 2:2:(rounds + 2)
        p += 1 / (0.25-i*i)
    end
    return 2p
end

```

Although it’s negligeable compared to the division cost indeed.

EDIT : The @simd seems to be by default, so that :

```julia
f4(n) = 4-2*sum(i -> 1/(i*i-0.25), 2:2:(n+2))

```

is actually enough…

---

<div class="post-metadata">

**Author:** ![Salmon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/salmon/32/22968_2.png) [@Salmon](https://discourse.julialang.org/u/Salmon)\
**Post date:** [November 6, 2022, 11:34am UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/12 "2022-11-06T11:34:26Z")

</div>

> [@lrnv](#):
>
> EDIT : The @simd seems to be by default, so that :
> 
> ```julia
> f4(n) = 4-2*sum(i -> 1/(i*i-0.25), 2:2:(n+2))
> 
> ```
> 
> is actually enough…

This is pretty beatiful, especially when you compare this nice one-liner to the equally fast avx c++ code which looks quite cumbersome 🙂  
would be quite cool if this would be the julia code in the benchmark I think. While you can optimize the hell out of your code if you need performance, It is especially cool that most of the time you already hit the mark by simply writing the code that you might think of first.

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 6, 2022, 1:07pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/13 "2022-11-06T13:07:30Z")

</div>

> [@Salmon](#):
>
> It is especially cool that most of the time you already hit the mark by simply writing the code that you might think of first.

Yeah, except telling people that “this is the code that we think of first” is clearly a lie in this case 😉

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [November 6, 2022, 2:28pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/14 "2022-11-06T14:28:07Z")

</div>

Also there might be a little correctness problem. Contracting two terms for efficiency needs a correction on even or odd inputs (but it is an O(1) correction).

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 6, 2022, 2:58pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/15 "2022-11-06T14:58:24Z")

</div>

> [@Dan](#):
>
> little correctness problem. Contracting two terms for efficiency needs a correction on even or odd inputs (but it is an O(1) correction).

Well,

```julia
julia> [f.(1:10) f4.(1:10)]
10×2 Matrix{Float64}:
 3.46667 3.46667
 2.89524 3.33968
 3.33968 3.33968
 2.97605 3.28374
 3.28374 3.28374
 3.01707 3.25237
 3.25237 3.25237
 3.04184 3.23232
 3.23232 3.23232
 3.0584 3.2184

```

And thus the value is the same when n is even, but is the next one then n is odd, so we are in fact _more precise_ without the correction. It is not that we omit one iteration, but we do one more.

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [November 6, 2022, 5:27pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/16 "2022-11-06T17:27:13Z")

</div>

If combining 2 terms leads to a big saving… why not keep going?

```julia
using Underscores

g(n,t=1) = @_ 4*sum(ifelse(isodd(_),-1,1)*(1/(2*_*t+1)), 0:(n+1))
g2(n,t=1) = 4*(@_ sum((2/(4*(_*t+1)^2-1)), 0:2:(n+1)))
g4(n,t=1) = 4*(sum(x -> (b = (x+2)*t; a = 4*b^2; 4*((a+3)/((a+3)^2-(8*b)^2))), 0:4:(n+1)))

```

we have:

```julia
julia> g(6)
3.017071817071818

julia> g(6) ≈ g2(6) ≈ g4(6)
true

```

But the real Julia fun begins with the second parameter:

```julia
julia> @_ map(g(_,1//1),1:6)
6-element Vector{Rational{Int64}}:
     52//15
    304//105
   1052//315
  10312//3465
 147916//45045
 135904//45045

julia> @_ map(g2(_,1//1),1:6)
6-element Vector{Rational{Int64}}:
    304//105
    304//105
  10312//3465
  10312//3465
 135904//45045
 135904//45045

```

and also:

```julia
julia> using Symbolics

julia> @variables k
1-element Vector{Num}:
 k

julia> g(4,k)
4.0 + -4 / (1 + 2k) + 4 / (1 + 4k) + -4 / (1 + 6k) + 4 / (1 + 8k) + -4 / (1 + 10k)

julia> g2(6,k)
2.6666666666666665 + 8 / (4((1 + 2k)^2) - 1) + 8 / (4((1 + 4k)^2) - 1) + 8 / (4((1 + 6k)^2) - 1)

julia> g4(6,k)
2.895238095238095 + (48 + 64((2 + 4k)^2)) / ((3 + 4((2 + 4k)^2))^2 - ((16 + 32k)^2))

```

and `g4` is faster:

```julia
julia> begin
       @btime g(1000)
       @btime g2(1000)
       @btime g4(1000)
       end
  1.287 μs (0 allocations: 0 bytes)
  989.600 ns (0 allocations: 0 bytes)
  765.291 ns (0 allocations: 0 bytes)
3.1405966379005616

```

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 6, 2022, 5:38pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/17 "2022-11-06T17:38:46Z")

</div>

Again, removing a few multiplications give a ~20% speedup (`g5`), and removing the `t` gives a ~20% more (`g6`):

```julia
using Underscores

g(n,t=1) = @_ 4*sum(ifelse(isodd(_),-1,1)*(1/(2*_*t+1)), 0:(n+1))
g2(n,t=1) = 4*(@_ sum((2/(4*(_*t+1)^2-1)), 0:2:(n+1)))
g4(n,t=1) = 4*(sum(x -> (b = (x+2)*t; a = 4*b^2; 4*((a+3)/((a+3)^2-(8*b)^2))), 0:4:(n+1)))
g5(n,t=1) = 16*(sum(x -> (b = 4((x+2)*t)^2; a = b+3; a/(a^2-16b)), 0:4:(n+1)))

g6(n) = 16*(sum(x -> (b = 4(x+2)^2; a = b+3; a/(a^2-16b)), 0:4:(n+1)))

@assert g(6) ≈ g2(6) ≈ g4(6) ≈ g5(6) ≈ g6(6)

using BenchmarkTools
@btime g(10000)
@btime g2(10000)
@btime g4(10000)
@btime g5(10000)
@btime g6(10000)

```

outputs :

```julia
  4.429 μs (0 allocations: 0 bytes)
  2.444 μs (0 allocations: 0 bytes)
  1.930 μs (0 allocations: 0 bytes)
  1.690 μs (0 allocations: 0 bytes)
  1.360 μs (0 allocations: 0 bytes)

```

---

<div class="post-metadata">

**Author:** ![mcabbott](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcabbott/32/6603_2.png) [@mcabbott](https://discourse.julialang.org/u/mcabbott)\
**Post date:** [November 6, 2022, 5:47pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/18 "2022-11-06T17:47:45Z")

</div>

These are fast but lose accuracy at large `n`:

```julia
julia> g2(10^6) - pi # similar to original
-9.999980012942444e-7

julia> g6(10^6) - pi # similar for g4, g5. Much worse at 10^9
-2.1787090125053705

```

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [November 6, 2022, 5:57pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/19 "2022-11-06T17:57:38Z")

</div>

The reason for accuracy loss, is Int64 overflow (there is no numerical methods reason for accuracy loss).

And with:

```julia
g6f(n) = 16*(sum(x -> (b = 4(x+2.0)^2; a = b+3; a/(a^2-16b)), 0:4:(n+1)))

```

(note `+2.0` to get Float64)

```julia
julia> g6f(10^6) - pi
-9.999960020046217e-7

julia> g(10^6) - pi
-9.999979990737984e-7

```

and it is faster in Float64,

```julia
julia> @btime g6f(10000)
  3.881 μs (0 allocations: 0 bytes)
3.1414926935740497

julia> @btime g6(10000)
  6.365 μs (0 allocations: 0 bytes)
3.1414926935740497

```

---

<div class="post-metadata">

**Author:** ![mcabbott](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcabbott/32/6603_2.png) [@mcabbott](https://discourse.julialang.org/u/mcabbott)\
**Post date:** [November 6, 2022, 6:07pm UTC](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766/20 "2022-11-06T18:07:47Z")

</div>

That must be it. For me the fix slows it down again, unfortunately (M1 mac):

```julia
julia> @btime f_in(10^6) - pi # top message
  335.333 μs (0 allocations: 0 bytes)
-9.999979067032427e-7

julia> @btime f4(10^6) - pi # from above, g2 similar
  230.916 μs (0 allocations: 0 bytes)
9.999970012053438e-7

julia> @btime g6(10^6) - pi # with overflow
  134.291 μs (0 allocations: 0 bytes)
-2.1787090125053705

julia> @btime g6f(10^6) - pi # version with x+2.0
  194.208 μs (0 allocations: 0 bytes)
-9.999960020046217e-7

```

I believe the StepRange `0:4:(n+1)` may cause issues for static compilation (due to error cases, zero step?).

Edit: solving the StepRange also speeds this up for me:

```julia
julia> @btime g7(10^6) - pi # version from below
  185.083 μs (0 allocations: 0 bytes)
-1.0000000019161348e-6

julia> g8(n) = 16*(sum(x -> (b = float(8x-4)^2; a = b+3; a/(a^2-16b)), 1:n÷4));

julia> @btime g8(10^6) - pi # with UnitRange not StepRange
  129.416 μs (0 allocations: 0 bytes)
-1.0000000019161348e-6

```

[Next page](https://discourse.julialang.org/t/usage-of-lazy-vector-instead-of-function/89766.md?page=2)