# How to make this piece of code faster? optimising it or using multiple threads?

**URL:** <https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590>\
**Category:** New to Julia\
**Created:** [November 1, 2019, 3:17pm UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590 "2019-11-01T15:17:44Z")\
**Posts on this page:** 20\
**Page:** 2

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 10, 2019, 4:53am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/21 "2019-11-10T04:53:59Z")

</div>

I tried what you suggested. I get the following result

```julia
 6.819 s (121276468 allocations: 3.79 GiB)

```

---

<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:** [November 10, 2019, 8:25am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/22 "2019-11-10T08:25:15Z")

</div>

> [@ennvvy](#):
>
> cedg=setdiff([1:1:N;],inneighbors(A.G,i))

Allocating this `[1:1:N]` outside of the loop pretty much can’t be more costly overall… Its a fundamental performance of julia and C and basically all programming languages. Allocating memory is slower than accessing memory.

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 11, 2019, 5:12am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/23 "2019-11-11T05:12:46Z")

</div>

I tried using `filter!(x->!(x in inneighbors(A.G,i)),cedg)`, but it doesn’t seem to improve significantly.

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 11, 2019, 6:41am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/24 "2019-11-11T06:41:06Z")

</div>

@carstenbauer,@Karajan,@Raf:

1. I realise that most of the computations yield results that are of type Float64. Is there a way to ensure that by default all the inputs and outputs are in either `Int32` or `Float32` format? Would this help save time?
2. Is `Int64[]` or `Float64[]` the best way to initialise an array whose size varies in every iteration?

---

<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:** [November 11, 2019, 9:06am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/25 "2019-11-11T09:06:02Z")

</div>

1. I doubt that is where you can gain performance. That’s a highly specialized optimization, that most likely will not help at all in this case.

2. `Int64[]` is best for initializing integer arrays, and `Float64[]` is best for initializing float arrays 😉

Perhaps you can update your posted code? You have replied that you have changed several things, but I cannot see that you are updating the code itself.

---

<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:** [November 11, 2019, 9:16am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/26 "2019-11-11T09:16:22Z")

</div>

> [@ennvvy](#):
>
> Is there a way to ensure that by default all the inputs and outputs are in either `Int32` or `Float32` format?

Yes.

> [@ennvvy](#):
>
> Would this help save time?

No.

Well, that’s not quite true. Potentially it can but until you have managed to remove all allocations within your loops it would at best be of marginal benefit.

The key at this point is to understand your algorithm. What do you actually need to compute and what do you actually need to store? And even more importantly, what do you have to recompute within each loop and what can be precomputed outside the loop?

To take a concrete example, that has been highlighted as expensive in the profiling, what do you need the setdiff computation **for**? You clearly don’t need it long term since it’s discarded fairly quickly. What it’s used for is to iterate over the nodes in the graph that are not inneighbors to the current node. Maybe you can just loop over all nodes and query the graph in some other way whether there’s an edge between the nodes in question? Depending on how the graph is stored this may or may not be more efficient. Or you can use a temporary array, that you allocate only once, to mark the inneighbors and refer to that in a loop over all nodes, like in this modification (I also de-aliased the uses of `edg` but otherwise left things as they were):

```julia
function main(Cmin::Float64,Cmax::Float64,N::Int64,n::Int64,outdeg::Int64,Pᵣ::Float64,Pᵢ::Float64,P::Int64,μ::Float64,σ::Float64)
    A=Arr(Array{Float64}(undef,N),Array{Float64}(undef,N),watts_strogatz(N,outdeg,Pᵣ,is_directed=true),Float64[],Any[])
    rand!(Normal(μ,σ),A.θ)
    round.(rand!(Uniform(Cmin,Cmax),A.C),digits=2)
    A.W=1.0.-A.C
    A.W/=sum(A.W)
    current_neighbors = falses(N)
    @inbounds for t in 1:n
        push!(A.θₜ,A.θ)
        upagents=Int64[]
        upagents=sample([1:1:N;],Weights(A.W),P,replace=false)
        @inbounds for i in upagents
            current_neighbors .= false
            for j in inneighbors(A.G, i)
                current_neighbors[j] = true
            end
            current_neighbors[i] = true
            
            edg = Int[]
            for j = 1:N
                if !current_neighbors[j] && rand() < Pᵢ
                    push!(edg, j)
                end
            end

            cedg=exp.(-abs.(A.θₜ[t][i].-A.θₜ[t][edg])/(1.0-A.C[i]))*(1.0-A.C[i])
            cedg/=sum(cedg)
            s=0.0
            c=0.0
            @inbounds for k in eachindex(cedg)
                s+=cedg[k]*sin(A.θₜ[t][edg[k]])
                c+=cedg[k]*cos(A.θₜ[t][edg[k]])
            end
            s+=A.C[i]*sin(A.θₜ[t][i])
            c+=A.C[i]*cos(A.θₜ[t][i])
            A.θ[i] = atan(s,c)
        end
        edg2=round.(A.θ,digits=2)
        if(all(y->y.==edg2[1],edg2)==true)
            break
        end
    end
    return A.θₜ
end

```

The next thing you need, and maybe you already have this, is a way to verify that the implementation is doing the right thing. My changes speed things up but maybe I made a mistake and it computes the wrong result faster. That you really need to know.

For more gains think about this. Do you _really_ need to store `edg` and `cedg` at all or can you make the necessary computations on the fly?

---

<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:** [November 11, 2019, 2:22pm UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/27 "2019-11-11T14:22:08Z")

</div>

> [@DNF](#):
>
> `Int64[]` is best for initializing integer arrays

you mean `Int[]` so it doesn’t necessarily fail on Win32

---

<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:** [November 11, 2019, 2:35pm UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/28 "2019-11-11T14:35:20Z")

</div>

I considered that, but concluded that @ennvvy specifically wanted `Int64`.

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 18, 2019, 11:09am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/29 "2019-11-18T11:09:14Z")

</div>

@carstenbauer,@Karajan,@Raf,@DNF: Thank you so much for the leads. I have now managed to get the following performance for the same parametric values as earlier:

```julia
using Distributions,Distributed,StatsBase,LinearAlgebra,LightGraphs,BenchmarkTools,Random

mutable struct Arr
    t::Int64
    agents::Array{Int64,1}
    θ::Array{Float64,1}
    C::Array{Float64,1}
    G::SimpleDiGraph{Int64}
    W::Array{Float64,1}
end

function main(Cmin::Float64,Cmax::Float64,N::Int64,n::Int64,outdeg::Int64,Pᵣ::Float64,Pᵢ::Float64,P::Int64,μ::Float64,σ::Float64)
    A=Arr(0,[1:N;],Array{Float64}(undef,N),Array{Float64}(undef,N),watts_strogatz(N,outdeg,Pᵣ,is_directed=true),
    Array{Float64}(undef,N))
    rand!(Normal(μ,σ),A.θ)
    rand!(Uniform(Cmin,Cmax),A.C)
    A.W.=1.0.-A.C
    A.W./=sum(A.W)
    θₜ=Vector{Float64}[]
    upagents=Array{Int64}[]
    cedg=Array{Float64}(undef,N)
    while(iszero(1 .-isapprox.(A.θ[1],A.θ,atol=1e-2))==false)
        A.t+=1
        push!(θₜ,copy(A.θ))
        upagents=sample(A.agents,Weights(A.W),P,replace=false)
        @inbounds for i in upagents
            cedg.=0.0
            @inbounds for j in inneighbors(A.G, i) 
              if(rand()<Pᵢ && i!=j)
                cedg[j]=exp(-abs(θₜ[A.t][i]-θₜ[A.t][j])/(1.0-A.C[i]))
              end
            end
            cedg./=(sum(cedg)/(1.0-A.C[i]))
            cedg[i]=A.C[i]
            s=0.0
            c=0.0
            @inbounds for k in 1:N
              s+=cedg[k]*sin(θₜ[A.t][k])
              c+=cedg[k]*cos(θₜ[A.t][k])
            end
            A.θ[i]=atan(s,c)
        end
    end
    return A.t,A.θ[1],sum(A.C.*θₜ[1])/sum(A.C),mean(A.C)
end

```

```julia
Cmin=0.5;Cmax=0.9;N=1000;n=5000;outdeg=250;Pᵣ=0.0;Pᵢ=0.5;P=500;μ=pi/2;σ=pi/18;
@btime main(Cmin,Cmax,N,n,outdeg,Pᵣ,Pᵢ,P,μ,σ)

```

Yields:

```julia
2.674 s (14911 allocations: 9.13 MiB)

```

---

<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:** [November 18, 2019, 11:40am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/30 "2019-11-18T11:40:35Z")

</div>

There’s also a function `sincos` which is faster than sin and cos individually. Not sure if it’s a bottleneck here though.

---

<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:** [November 18, 2019, 2:58pm UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/31 "2019-11-18T14:58:33Z")

</div>

If `sin` and `cos` are bottlenecks there should be bigger gains from the observation that, except for the first iteration, their arguments are the result of an `atan` computation. It would take some extra bookkeeping and complexity to take advantage of it though.

But this looks like a very low-hanging fruit. No need to compute `sin` and `cos` just to multiply them by zero:

```julia
            @inbounds for k in 1:N
                if cedg[k] != 0
                    s+=cedg[k]*sin(θₜ[A.t][k])
                    c+=cedg[k]*cos(θₜ[A.t][k])
                end
            end

```

---

<div class="post-metadata">

**Author:** ![Vasily\_Pisarev](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vasily_pisarev/32/7929_2.png) [@Vasily\_Pisarev](https://discourse.julialang.org/u/Vasily_Pisarev)\
**Post date:** [November 18, 2019, 7:56pm UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/32 "2019-11-18T19:56:34Z")

</div>

```julia
1 .- isapprox.(A.θ[1],A.θ,atol=1e-2)

```

is clearly a source of allocations.  
I’d formulate the loop condition as

```julia
while any(abs(x - A.θ[1]) > 0.01 for x in A.θ)

```

Also,

```julia
upagents = sample(...)

```

is a candidate for optimization. It appears that you want

```julia
sample!(A.agents, Weights(A.W), P, upagents, replace=false)

```

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 19, 2019, 4:56am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/35 "2019-11-19T04:56:17Z")

</div>

@GunnarFarneback: Your suggestion significantly lowered the execution time to the order of `ms`  
@Vasily_Pisarev:Your suggestions did lower the allocations

The revised code:

```julia
using Distributions,Distributed,StatsBase,LinearAlgebra,LightGraphs,BenchmarkTools,Random,Plots

mutable struct Arr
    t::Int64
    agents::Array{Int64,1}
    θ::Array{Float64,1}
    C::Array{Float64,1}
    G::SimpleDiGraph{Int64}
    W::Array{Float64,1}
end

function main(Cmin::Float64,Cmax::Float64,N::Int64,n::Int64,outdeg::Int64,Pᵣ::Float64,Pᵢ::Float64,P::Int64,μ::Float64,σ::Float64)
    A=Arr(0,[1:N;],Array{Float64}(undef,N),Array{Float64}(undef,N),watts_strogatz(N,outdeg,Pᵣ,is_directed=true),
    Array{Float64}(undef,N))
    rand!(Normal(μ,σ),A.θ)
    rand!(Uniform(Cmin,Cmax),A.C)
    A.W.=1.0.-A.C
    A.W./=sum(A.W)
    θₜ=Vector{Float64}[]
    upagents=Array{Int64}(undef,round(Int,0.5*N))
    cedg=Array{Float64}(undef,N)
    while(any(abs(x - A.θ[1]) > 0.01 for x in A.θ))
        A.t+=1
        push!(θₜ,copy(A.θ))
        sample!(Random.GLOBAL_RNG,A.agents, Weights(A.W), upagents, replace=false)
        @inbounds for i in upagents
            cedg.=0.0
            @inbounds for j in inneighbors(A.G, i) 
              if(rand()<Pᵢ && i!=j)
                cedg[j]=exp(-abs(θₜ[A.t][i]-θₜ[A.t][j])/(1.0-A.C[i]))
              end
            end
            cedg./=(sum(cedg)/(1.0-A.C[i]))
            cedg[i]=A.C[i]
            s=0.0
            c=0.0
            @inbounds for k in 1:N
              if(cedg[k]!=0)
                s+=cedg[k]*sin(θₜ[A.t][k])
                c+=cedg[k]*cos(θₜ[A.t][k])
              end
            end
            A.θ[i]=atan(s,c)
        end
    end
    #return A.t,A.θ[1],sum(A.C.*θₜ[1])/sum(A.C),mean(A.C)
    return θₜ
end

```

```julia
Cmin=0.5;Cmax=0.9;N=1000;n=5000;outdeg=250;Pᵣ=0.0;Pᵢ=0.5;P=500;μ=pi/2;σ=pi/18;
@btime main(Cmin,Cmax,N,n,outdeg,Pᵣ,Pᵢ,P,μ,σ)
589.232 ms (14481 allocations: 6.67 MiB)

```

To be able to simulate for network sizes larger than 10,000, I am planning to tweak this code (the inner for loops) to be able to run on GPU. Is that the right approach? Using the above code, I get the following performance:

```julia
Cmin=0.5;Cmax=0.9;N=10000;n=5000;outdeg=2500;Pᵣ=0.0;Pᵢ=0.5;P=5000;μ=pi/2;σ=pi/18;
@btime main(Cmin,Cmax,N,n,outdeg,Pᵣ,Pᵢ,P,μ,σ)
68.707 s (220977 allocations: 663.53 MiB)

```

---

<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:** [November 19, 2019, 5:12am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/36 "2019-11-19T05:12:52Z")

</div>

With this modification

```julia
                if(cedg[k]!=0)
                    si,ci = sincos(θₜ[A.t][k])
                    s+=si
                    c+=ci
                end

```

I get it to run significantly faster

```julia
  452.775 ms (14568 allocations: 7.12 MiB) # Before
  54.714 ms (14106 allocations: 4.74 MiB) # After

```

It’s actually a bit strange that the runtime is cut in more than half by this. Maybe `sincos` allows further optimizations to be made, like simd:ing the loop or something.

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 19, 2019, 5:29am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/37 "2019-11-19T05:29:11Z")

</div>

Thats really surprising. When I tried using `sincos` there wasnt any obvious difference in the time and allocation. This is what I got:

```julia
672.970 ms (14515 allocations: 6.85 MiB)

```

---

<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:** [November 19, 2019, 5:42am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/38 "2019-11-19T05:42:36Z")

</div>

It seems your timing of 500+ ms is a bit higher than what I see, maybe you just have a different processor that does not have the same instruction set as the one I use?

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 19, 2019, 5:44am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/39 "2019-11-19T05:44:35Z")

</div>

I am currently running it on google colab. It has an Intel(R) Xeon(R) CPU @ 2.30GHz with a 12GB RAM

---

<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:** [November 19, 2019, 5:46am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/40 "2019-11-19T05:46:29Z")

</div>

Ah, then I would not be surprised if the full potential of the processor can not be used. Have you tried running the example locally, on a reasonably modern CPU?

---

<div class="post-metadata">

**Author:** ![ennvvy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ennvvy/32/7816_2.png) [@ennvvy](https://discourse.julialang.org/u/ennvvy)\
**Post date:** [November 19, 2019, 6:48am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/41 "2019-11-19T06:48:31Z")

</div>

I tried it on an 2.4 GHz Dual-Core Intel Core i5. It still shows similar performance.

```julia
667.504 ms (14488 allocations: 6.71 MiB)

```

Inside the `for loop` while adding `si` and `ci` values to `s` and `c` they are multiplied by weights `A.cedg[k]` in the original code. Could that make the difference?

---

<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:** [November 19, 2019, 6:55am UTC](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590/42 "2019-11-19T06:55:26Z")

</div>

> [@baggepinnen](#):
>
> It’s actually a bit strange that the runtime is cut in more than half by this. Maybe `sincos` allows further optimizations to be made, like simd:ing the loop or something.

It seems unlikely that that loop could get any automatic SIMD optimizations and there shouldn’t be any allocations to gain from your change. I would guess that the randomness has resulted in different iteration counts for your two runs.

Btw, it’s extremely hard to do a fair benchmarking with the randomness and the loop until convergence. Even with a fixed seed, `@btime` will vary the number of runs if you do something that significantly affects the runtime, which can give entirely different numbers.

[Previous page](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590.md?page=1)

[Next page](https://discourse.julialang.org/t/how-to-make-this-piece-of-code-faster-optimising-it-or-using-multiple-threads/30590.md?page=3)