# Very slow execution time in comparison even to Python

**URL:** <https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223>\
**Category:** Performance\
**Created:** [October 29, 2020, 6:25am UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223 "2020-10-29T06:25:39Z")\
**Posts on this page:** 16\
**Page:** 2

<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:** [October 29, 2020, 10:20am UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/21 "2020-10-29T10:20:28Z")

</div>

> [@moeddel](#):
>
> `rand(SMatrix{3,POINTS_PER_BATCH,Float64})`

I would be a bit careful about using `SMatrix` here. `SMatrix` is for cases where the compiler can be sure of the sizes, and the sizes should be small. In this particular case `POINTS_PER_BATCH` is small and fixed, but at least be aware of this.

---

<div class="post-metadata">

**Author:** ![moeddel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moeddel/32/18641_2.png) [@moeddel](https://discourse.julialang.org/u/moeddel)\
**Post date:** [October 29, 2020, 10:38am UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/22 "2020-10-29T10:38:33Z")

</div>

> I would be a bit careful about using `SMatrix` here. `SMatrix` is for cases where the compiler can be sure of the sizes, and the sizes should be small. In this particular case `POINTS_PER_BATCH` is small and fixed, but at least be aware of this.

For this example everything is known, since `POINTS_PER_BATCH` is known at compile time, but in general this can lead to problems. I am aware of this.

How fast does this run, if you use the above mentioned within your code?

---

<div class="post-metadata">

**Author:** ![GunnarFarneback](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gunnarfarneback/32/1827_2.png) [@GunnarFarneback](https://discourse.julialang.org/u/GunnarFarneback)\
**Post date:** [October 29, 2020, 10:47am UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/23 "2020-10-29T10:47:39Z")

</div>

> [@moeddel](#):
>
> Since we are only interested in one part of the SVD, i.e. U, we can calculate the SVD of MM’ like shown [here](https://www.rdocumentation.org/packages/corpcor/versions/1.6.9/topics/fast.svd). The SVDs of MM’ and M share U.

It’s even better than that. Since `MM'` is symmetric and positive definite, the SVD is equivalent to the eigen factorization, and StaticArrays has a specialized `eigen` for 3x3 matrices. Only problem is that `svd` and `eigen` doesn’t use the same ordering of singular/eigen values and vectors, so you need to take that into account.

---

<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:** [October 29, 2020, 12:18pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/25 "2020-10-29T12:18:17Z")

</div>

After playing around with this a bit, I think it should be possible to get the runtime down to approximately 1 second + printing (on my old laptop) for a batch size of 10^6. Then multithreading comes on top of that.

That’s compared to an original runtime of 18 minutes, including printing and compilation.

---

<div class="post-metadata">

**Author:** ![moeddel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moeddel/32/18641_2.png) [@moeddel](https://discourse.julialang.org/u/moeddel)\
**Post date:** [October 29, 2020, 12:31pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/26 "2020-10-29T12:31:11Z")

</div>

> [@GunnarFarneback](#):
>
> It’s even better than that. Since `MM'` is symmetric and positive definite, the SVD is equivalent to the eigen factorization, and StaticArrays has a specialized `eigen` for 3x3 matrices. Only problem is that `svd` and `eigen` doesn’t use the same ordering of singular/eigen values and vectors, so you need to take that into account.

That does indeed work. It seems the eigenvalues and corresponding vectors are ordered in reverse of what the SVD does. So here my updated code. The speedup does not so much depend on `POINTS_PER_BATCH`. For `POINTS_PER_BATCH` of `8` and `256` the “fast SVD trick” does provide a speed up of about 2.

```julia
using LinearAlgebra, Statistics, Random, StaticArrays, BenchmarkTools

struct Plane
    point::SVector{3,Float64}
    normal::SVector{3,Float64}
end

function mean_center(points::SMatrix)
    centroid = mean(points,dims=2)
    points_centered = points .- centroid
    centroid, points_centered
end

function best_fit_plane(points::SMatrix)
    centroid, centered_points = mean_center(points)

	# slowest
	#u = svd(centered_points).U[:,end]

	# better
	#m = centered_points*transpose(centered_points)
	#u = svd(m).U[:,end]

	# best
	#m = Symmetric(centered_points*transpose(centered_points))
	#u = eigen(m).vectors[:,1]

    Plane(centroid, u)
end

function main()
	MIN = -2000
	MAX = 2000
	BATCHES_COUNT = 1250
	POINTS_PER_BATCH = 2^3 # or 2^8
	BATCHES_COUNT = div(2^16,POINTS_PER_BATCH)
	
	generateBatch() = MIN .+ (MAX - MIN)*rand(SMatrix{3,POINTS_PER_BATCH,Float64})

	out = Vector{Plane}(undef,BATCHES_COUNT) # initialize output
	for i in eachindex(out)
		b = generateBatch()
		out[i] = best_fit_plane(b)
	end
	return out
end

@benchmark main()

```

We start off here for `POINTS_PER_BATCH = 8` using the slowest approach

```julia
BenchmarkTools.Trial: 
  memory estimate: 62.38 MiB
  allocs estimate: 434178
  --------------
  minimum time: 129.477 ms (2.97% GC)
  median time: 139.976 ms (2.76% GC)
  mean time: 146.015 ms (3.71% GC)
  maximum time: 214.297 ms (1.98% GC)
  --------------
  samples: 35
  evals/sample: 1

```

The slow approach is slightly better

```julia
BenchmarkTools.Trial: 
  memory estimate: 58.63 MiB
  allocs estimate: 450562
  --------------
  minimum time: 111.597 ms (3.00% GC)
  median time: 116.901 ms (3.71% GC)
  mean time: 122.655 ms (4.18% GC)
  maximum time: 201.346 ms (4.26% GC)
  --------------
  samples: 41
  evals/sample: 1

```

The fast approach yields the best results

```julia
BenchmarkTools.Trial: 
  memory estimate: 36.88 MiB
  allocs estimate: 393218
  --------------
  minimum time: 57.961 ms (0.00% GC)
  median time: 66.490 ms (6.86% GC)
  mean time: 67.735 ms (6.30% GC)
  maximum time: 87.808 ms (7.00% GC)
  --------------
  samples: 74
  evals/sample: 1

```

We can get better results for larger `POINTS_PER_BATCH = 256` using the slowest approach

```julia
BenchmarkTools.Trial: 
  memory estimate: 25.48 MiB
  allocs estimate: 13825
  --------------
  minimum time: 14.025 ms (0.00% GC)
  median time: 17.141 ms (11.42% GC)
  mean time: 16.768 ms (7.38% GC)
  maximum time: 25.514 ms (9.19% GC)
  --------------
  samples: 299
  evals/sample: 1

```

The slow approach is slightly better

```julia
BenchmarkTools.Trial: 
  memory estimate: 19.49 MiB
  allocs estimate: 14337
  --------------
  minimum time: 8.721 ms (0.00% GC)
  median time: 11.221 ms (0.00% GC)
  mean time: 11.340 ms (8.76% GC)
  maximum time: 17.238 ms (14.03% GC)
  --------------
  samples: 441
  evals/sample: 1

```

The fast approach yields the best results

```julia
BenchmarkTools.Trial: 
  memory estimate: 18.81 MiB
  allocs estimate: 12545
  --------------
  minimum time: 6.666 ms (0.00% GC)
  median time: 8.965 ms (0.00% GC)
  mean time: 8.979 ms (10.28% GC)
  maximum time: 14.986 ms (26.97% GC)
  --------------
  samples: 557
  evals/sample: 1

```

> [@DNF](#):
>
> After playing around with this a bit, I think it should be possible to get the runtime down to approximately 1 second + printing (on my old laptop) for a batch size of 10^6. Then multithreading comes on top of that.

Around 1 sec on my machine for a batch size 2^{20} and `POINTS_PER_BATCH = 8`. No parallelization either.

---

<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:** [October 29, 2020, 1:22pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/27 "2020-10-29T13:22:31Z")

</div>

> [@moeddel](#):
>
> ```julia
> function mean_center(points::SMatrix)
> centroid = mean(points,dims=2)
> points_centered = points .- centroid
> centroid, points_centered
> end
> 
> ```

It looks to me like `mean(points; dims=Val(2))` improves performance slightly (actually, it dramatically speeds up `mean_center`, but the overall performance isn’t so much affected.)

> [@moeddel](#):
>
> `u = eigen(m).vectors[:,1]`

I think it should be the last, not the first eigenvector. Also, I think it’s better to not access the fields like that. You can either use

```julia
_, V = eigen(m)
u = V[:, end]

```

or

```julia
u = eigvecs(m)[:, end]

```

---

<div class="post-metadata">

**Author:** ![mschauer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mschauer/32/13946_2.png) [@mschauer](https://discourse.julialang.org/u/mschauer)\
**Post date:** [October 29, 2020, 1:27pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/28 "2020-10-29T13:27:33Z")

</div>

_Cunningham’s Law_ states “the best way to get the right answer on the internet is not to ask a question; it’s to post the wrong answer.” I think this thread suggests a new corollary 😉

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [October 29, 2020, 1:35pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/29 "2020-10-29T13:35:36Z")

</div>

This is great from the first time unoptimized code posted here has been optimized to death. The main surprising part is that so far this thread have gotten to the typical 1000x speedup 😉

---

<div class="post-metadata">

**Author:** ![moeddel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moeddel/32/18641_2.png) [@moeddel](https://discourse.julialang.org/u/moeddel)\
**Post date:** [October 29, 2020, 1:58pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/30 "2020-10-29T13:58:44Z")

</div>

> [@DNF](#):
>
> I think it should be the last, not the first eigenvector. Also, I think it’s better to not access the fields like that. You can either use
> 
> ```julia
> _, V = eigen(m)
> u = V[:, end]
> 
> ```
> 
> or
> 
> ```julia
> u = eigvecs(m)[:, end]
> 
> ```

I have tried this out and it seems `eigen` and `svd` have a different order. `svd` puts the largest eigenvalue first, whereas `eigen` of a symmetric 3\times 3 matrix seems to put the largest eigenvalue last.

```julia
using Test, StaticArrays
a = rand(SMatrix{3,8,Float64})
U, = svd(a)
λs, Ue = eigen(Symmetric(a*transpose(a)))
@test U[:,end] ≈ -Ue[:,1]

```

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [October 29, 2020, 2:05pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/31 "2020-10-29T14:05:10Z")

</div>

Sharing a version “slower than Python” on julia discourse is best way to get the fastest version of your algorithm ever written?

---

<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:** [October 29, 2020, 2:17pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/32 "2020-10-29T14:17:55Z")

</div>

> [@moeddel](#):
>
> it seems `eigen` and `svd` have a different order.

I know, that’s why I said it. But I managed to swap the two cases in my mind 🤦‍♂️

---

<div class="post-metadata">

**Author:** ![Yury\_Kotlyarov](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/yury_kotlyarov/32/19020_2.png) [@Yury\_Kotlyarov](https://discourse.julialang.org/u/Yury_Kotlyarov)\
**Post date:** [October 29, 2020, 2:21pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/33 "2020-10-29T14:21:55Z")

</div>

yeah, sounds like "an idiomatic way to get a lot of support ". really unexpected attention to our small issue…

thanks to all of you, guys! now we are carefully checking all your suggestions. Of course, the most time-consuming call is `println` to the terminal.

btw now code is 10x faster than haskell version, will try to do the same trick on the haskell discourse

---

<div class="post-metadata">

**Author:** ![moeddel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moeddel/32/18641_2.png) [@moeddel](https://discourse.julialang.org/u/moeddel)\
**Post date:** [October 29, 2020, 2:29pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/34 "2020-10-29T14:29:56Z")

</div>

Could you share the specifics, on how you did the benchmark in the end?

---

<div class="post-metadata">

**Author:** ![Yury\_Kotlyarov](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/yury_kotlyarov/32/19020_2.png) [@Yury\_Kotlyarov](https://discourse.julialang.org/u/Yury_Kotlyarov)\
**Post date:** [October 29, 2020, 2:31pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/35 "2020-10-29T14:31:09Z")

</div>

Sure, will post tomorrow. First of all we are going to rewrite all three source bases to exclude random and writing to terminal. Next thing to we are going to replace Haskell’s `String` with `Text`, because `String` in Haskell is linked list of `Char`s

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [October 29, 2020, 3:55pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/36 "2020-10-29T15:55:38Z")

</div>

Ahhh - recursive discourse sniping.

“Very slow exectution time in comparison to Julia”

---

<div class="post-metadata">

**Author:** ![moeddel](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/moeddel/32/18641_2.png) [@moeddel](https://discourse.julialang.org/u/moeddel)\
**Post date:** [October 29, 2020, 9:38pm UTC](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223/37 "2020-10-29T21:38:53Z")

</div>

While splitting up the procedure into random data generation, I remembered that each call to `rand` has an overhead as stated in the documentation, so I replaced it with a single call. Now we are yet again faster.

```julia
using LinearAlgebra, Statistics, Random, StaticArrays, BenchmarkTools

struct Plane{M,T}
    point::SVector{M,T}
    normal::SVector{M,T}
end

function mean_center(points::SMatrix{M,N,T}) where {M,N,T}
    centroid = mean(points,dims=2)
    points_centered = points .- centroid
    centroid, points_centered
end

function best_fit_plane(points::SMatrix{M,N,T}) where {M,N,T}
    centroid, centered_points = mean_center(points)
	m = Symmetric(centered_points*transpose(centered_points))
	λs,U = eigen(m)
    Plane{M,T}(centroid, U[:,1])
end

function generateRandomBatchData()
	MIN = -2000
	MAX = 2000
	POINTS_PER_BATCH = 2^3
	BATCHES_COUNT = div(800000,POINTS_PER_BATCH)
	DIM = 3
	ELEMENT_TYPE = Float64
	
	randomBatchData = MIN .+ (MAX - MIN)*rand(SMatrix{DIM,POINTS_PER_BATCH,ELEMENT_TYPE},BATCHES_COUNT)
	return randomBatchData
end

function coreComputation(batchData)
	return map(best_fit_plane,batchData)
end

```

The data generation is now now significantly faster than the main computation just by eliminating the repeated calls to `rand`.

```julia
@benchmark generateRandomBatchData()
BenchmarkTools.Trial: 
  memory estimate: 80.10 MiB
  allocs estimate: 399500
  --------------
  minimum time: 23.918 ms (0.00% GC)
  median time: 46.357 ms (0.00% GC)
  mean time: 45.681 ms (16.30% GC)
  maximum time: 94.225 ms (25.72% GC)
  --------------
  samples: 110
  evals/sample: 1

```

The main computation is now well below 1s without multi-threading

```julia
batchData = generateRandomBatchData()
@benchmark coreComputation($batchData)
BenchmarkTools.Trial: 
  memory estimate: 215.15 MiB
  allocs estimate: 1900002
  --------------
  minimum time: 326.180 ms (7.11% GC)
  median time: 355.630 ms (6.67% GC)
  mean time: 357.212 ms (8.27% GC)
  maximum time: 418.457 ms (5.71% GC)
  --------------
  samples: 14
  evals/sample: 1

```

So at best something around 350ms on my laptop, which is not the fastest.

[Previous page](https://discourse.julialang.org/t/very-slow-execution-time-in-comparison-even-to-python/49223.md?page=1)