# Variation in computation time for execution of the exact same code

**URL:** <https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333>\
**Category:** Performance\
**Created:** [March 22, 2020, 5:15am UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333 "2020-03-22T05:15:42Z")\
**Posts on this page:** 16\
**Page:** 2

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [March 24, 2020, 2:35pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/21 "2020-03-24T14:35:08Z")

</div>

[https://github.com/JuliaCI/BenchmarkTools.jl/blob/master/doc/manual.md#interpolating-values-into-benchmark-expressions](https://github.com/JuliaCI/BenchmarkTools.jl/blob/master/doc/manual.md#interpolating-values-into-benchmark-expressions)

---

<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:** [March 24, 2020, 3:40pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/22 "2020-03-24T15:40:59Z")

</div>

> [@jleman](#):
>
> Regarding the extra dimension of U: Later in the code (in a line I deleted from this abbreviated block, but shown below). I copy the values of U to another matrix that has to be in three dimensions to play well with other calculations.

But `stepSource` returns a 2d array, not a 3d one. What I’m saying is that _inside_ `stepSource` you create a 3d array, and then create a 2d slice which you return. The fact that it started with a 3d array has no effect outside the function, except that it causes an extra memory allocation.

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [March 24, 2020, 5:24pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/23 "2020-03-24T17:24:44Z")

</div>

That makes sense. My recollection is that that using slices forced me into extra dimensions due to later calculations with other 3D elements. However, I am reworking the code with the StaticArray format you suggested. As an example here is what I have attempted for the three-phase sine source. Each sub-vector stores three voltages (one for each phase of a three-phase voltage source) at a given time step. This seems to output the right numerical data. One item that puzzles me is illustrated in the print statements. The first two make sense. In the third I am trying a quick way to get all the first elements of each vector for later plotting purposes. I could do this with another variable and loop, but is there an easier way?

```julia
using Plots
using LinearAlgebra
using StaticArrays

function threePhaseSine(θ)
    Δt = 2.0e-3; ω = 2.0*π*60.0 #
    Aphase(nt) = sqrt(2/3)*cos.(ω * Δt * nt .+ θ) #nt is time step
    Bphase(nt) = sqrt(2/3)*cos.(ω * Δt * nt .- 2*π/3 .+ θ)
    Cphase(nt) = sqrt(2/3)*cos.(ω * Δt * nt .+ 2*π/3 .+ θ)
    Vsin = [MVector{3,Float64}([Aphase(nt), Bphase(nt), Cphase(nt)]) for nt in 1:10]

    return Vsin 
end

Vabc = threePhaseSine(0)

display(Vabc)
println("First vector of Vabc = ",Vabc[1])
println("First element of vector 1 of Vabc = ",Vabc[1][1])
println("First element of all vectors of Vabc = ",Vabc[:][1]) #this is unexpected

plot(Vabc[:][1]) #attempt to plot the voltage for Phase A as a function of time (doesn't work)

```

---

<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:** [March 24, 2020, 6:20pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/24 "2020-03-24T18:20:00Z")

</div>

`Vabs[:]` is just a copy of `Vabs`, so `Vabs[:][1]` is the same as writing `Vabs[1]`, except you make an extra copy.

You probably want `first.(Vabs)`.

> [@jleman](#):
>
> ```julia
> Vsin = [MVector{3,Float64}([Aphase(nt), Bphase(nt), Cphase(nt)])
> 
> ```

Drop the square brackets around the phase functions here, this allocates a normal array and then converts it to an MVector.

In fact, it’s better to drop the `{3,Float64}` since the former is automatically understood from the number of inputs, and the element type is inferred, and will then allow also other types than `Float64`, so:

```julia
Vsin = [MVector(Aphase(nt), Bphase(nt), Cphase(nt)) for nt in 1:10]

```

Also consider whether you can use an SVector instead.

BTW: The dots you put inside the `Aphase` etc. functions don’t really make any sense, since `threePhaseSine` only works for scalar input. You can, however, broadcast the _outer_ function over vector input, as `threePhaseSine.(θ)`.

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [March 24, 2020, 6:59pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/25 "2020-03-24T18:59:55Z")

</div>

Roger on the dots. I see there is no `second.(Vabs)` (i.e. for accessing the second or third elements of each vector). Is there a more general way to access?

Likewise, if I am using MVector or Mmatrix elements (because they have to be updated), how do I index certain blocks and update in place? For example in the code below, Vout is a vector of vector of vectors (?? multiple conductors, space, and time). I want to write a small array of Voltages (trappedChargeVoltage) to each conductor (all nc), each spatial step (all NDZ), but only the first time step (NDT = 1).

```julia
using StaticArrays
using LinearAlgebra

nc = 3
NDZ = 5
NDT = 10

VV = [MVector{nc,Float64}(zeros(nc)) for _ in 1:NDZ] #initialize vector of vectors to store update voltages at eachs spatial node
Vout = [[MVector{nc,Float64}(zeros(nc)) for _ in 1:NDZ] for _ in 1:NDT] #initialize array to store all output variables
#Vout contains voltages for conductors (nc) at at each spatial step (NDZ), and each time step (NDT)

trappedChargeV = [SVector{3,Float64}(10.0,0.0,-10.0)]
VV .= trappedChargeV
first.(Vout) .= trappedChargeV
display(first.(Vout))

```

---

<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:** [March 24, 2020, 7:22pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/26 "2020-03-24T19:22:00Z")

</div>

> [@jleman](#):
>
> I see there is no `second.(Vabs)` (i.e. for accessing the second or third elements of each vector. is there a more general way access?

There is `first` and `last`, otherwise use `getindex.(Vabs, n)` to get the nth index.

> [@jleman](#):
>
> Likewise, if I am using MVector or Mmatrix elements (because they have to be updated)

You can update with an `SVector` too, but then you have to replace the entire subvector, for example:

```julia
VV = [zero(SVector{nc, Float64}) for _ in 1:NDZ]
# or 
VV = Vector{SVector{nc, Float64}}(undef, NDZ)

trappedChargeV = SVector(10.0,0.0,-10.0)
VV .= (trappedChargeV,)

```

Note that again you create unnecessary temporary arrays with `MVector{nc,Float64}(zeros(nc))` since `zeros(nc)` creates an ordinary array which is converted to a StaticArray. Use unstead `zero(MVector{3,Float64})`, etc.

But here you may as well make it easier on yourself by just writing

```julia
trappedChargeV = SVector(10.0,0.0,-10.0) # <- no need for {3, Float64}
VV = fill(trappedChargeV, NDZ)

```

This will not work:

```julia
first.(Vout) .= trappedChargeV

```

The left side creates a temporary array and assigns to that, and it disappears without touching `Vout`.

I would do this:

```julia
for v in Vout
    v[1] = trappedChargeV
end
# or
setindex!.(Vout, (trappedChargeV,), 1)

```

Notice that this works for `SVector`s.

---

<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:** [March 24, 2020, 8:36pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/27 "2020-03-24T20:36:59Z")

</div>

Note that keeping `nc = 3` as a variable is a performance issue.

If you can, write them as literal values. Compare the following three functions. `foo` is yours, `bar` is the same but with literal instead of `nc`, and `baz` uses `SVector`. Also note that the first two don’t do what you (probably) want.

```julia
using BenchmarkTools, StaticArrays

function foo()
    nc = 3
    NDZ = 5
    NDT = 10

    VV = [MVector{nc,Float64}(zeros(nc)) for _ in 1:NDZ]
    
    trappedChargeV = [SVector{3,Float64}(10.0,0.0,-10.0)]
    VV .= trappedChargeV
    
    Vout = [[MVector{nc,Float64}(zeros(nc)) for _ in 1:NDZ] for _ in 1:NDT]
    first.(Vout) .= trappedChargeV
     
    return Vout, VV
end

function bar()
    NDZ = 5
    NDT = 10

    VV = [MVector{3,Float64}(zeros(3)) for _ in 1:NDZ]
    
    trappedChargeV = [SVector{3,Float64}(10.0,0.0,-10.0)]
    VV .= trappedChargeV
    
    Vout = [[MVector{3,Float64}(zeros(3)) for _ in 1:NDZ] for _ in 1:NDT]
    first.(Vout) .= trappedChargeV
     
    return Vout, VV
end

function baz()
    NDZ = 5
    NDT = 10

    trappedChargeV = SVector(10.0, 0.0, -10.0)
    VV = fill(trappedChargeV, NDZ)

    Vout = [[zero(SVector{3, Float64}) for _ in 1:NDZ] for _ in 1:NDT]
    setindex!.(Vout, (trappedChargeV,), 1)
     
    return Vout, VV
end

```

Benchmarks:

```julia
julia> @btime foo();
  74.461 μs (665 allocations: 45.28 KiB)

julia> @btime bar();
  3.046 μs (140 allocations: 10.03 KiB)

julia> @btime baz();
  642.392 ns (14 allocations: 2.58 KiB)

```

---

<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:** [March 24, 2020, 11:00pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/28 "2020-03-24T23:00:50Z")

</div>

It’s surprising to me that constant propagation didn’t handle the `foo` case… 😕

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [March 25, 2020, 8:05pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/29 "2020-03-25T20:05:51Z")

</div>

Is there a good way to make something like this work, or is it best just to go to MVector.

```julia
vec1 = SVector{3}(rand(3))
for i = 1:4    
    vec1 .+= SVector{3}(rand(3))
end

```

---

<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:** [March 25, 2020, 8:09pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/30 "2020-03-25T20:09:56Z")

</div>

Just do `vec1 += `.

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [March 25, 2020, 8:17pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/31 "2020-03-25T20:17:49Z")

</div>

It updates in-place, correct?

---

<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:** [March 25, 2020, 8:19pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/32 "2020-03-25T20:19:47Z")

</div>

Yeah, there’s no need to update in-place. You can think of `StaticArray`s as _values_ instead of containers.

> [@jleman](#):
>
> It updates in-place, correct?

No, it creates a new instance. But it will be super-fast, since it’s a stack-allocated value type.

Also, instead of

```julia
SVector{3}(rand(3))

```

write

```julia
rand(SVector{3})

```

then you will avoid creating an unnecessary temporary array.

Observe:

```julia
julia> function foo()
           vec1 = rand(MVector{3})
           for i = 1:4    
               vec1 .+= rand(MVector{3})
           end
           return vec1
       end

julia> function bar()
           vec1 = rand(SVector{3})
           for i = 1:4    
               vec1 += rand(SVector{3}) # <- no dot!
           end
           return vec1
       end

julia> @btime foo();
  88.801 ns (5 allocations: 160 bytes)

julia> @btime bar();
  74.856 ns (0 allocations: 0 bytes) # <- no allocations!

```

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [March 26, 2020, 6:52am UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/33 "2020-03-26T06:52:40Z")

</div>

Almost done implementing the staticarray approach and will post some figures showing the performance gains. Two more items that are stumping me.

1. I need to perform an elementwise absolute value operation on a result in the form of that shown below. I think I understand why what I have won’t work. But is there a straightforward way to do this?

```julia
vec1 = [[SVector{3,Float64}(-i,-i,-i) for i in 1:4] for _ in 1:10]
abs.(vec1)

```

1. I need to save the vec1 static-array output to a file. I normally save matrix results out as a csv file using writedlm, but that seems tedious as it keeps all the brackets. Any suggestions for a convenient way to store out the data? The intent would be to read the data back in for later post processing and plotting. Or is there a convenient way to convert the data back to a non-static array form where manipulating contents is a bit easier.

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [March 26, 2020, 8:21am UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/34 "2020-03-26T08:21:53Z")

</div>

> [@DNF](#):
>
> then you will avoid creating an unnecessary temporary array.

May be worth pointing out that with a couple modifications, `MArray`s don’t allocate either:

```julia
julia> @inline inlinerand(::Type{<:MVector{3}}) = MVector(rand(),rand(),rand())
inlinerand (generic function with 1 method)

julia> function foo2()
           vec1 = inlinerand(MVector{3})
           for i = 1:4
               vec1 .+= inlinerand(MVector{3})
           end
           return SArray(vec1)
       end
foo2 (generic function with 1 method)

julia> @btime foo2()
  47.758 ns (0 allocations: 0 bytes)
3-element SArray{Tuple{3},Float64,1,3} with indices SOneTo(3):
 2.4573027826314195
 2.1526819466919456
 2.4786087715660896

```

---

<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:** [March 26, 2020, 8:41am UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/35 "2020-03-26T08:41:17Z")

</div>

> [@jleman](#):
>
> ```julia
> vec1 = [[SVector{3,Float64}(-i,-i,-i) for i in 1:4] for _ in 1:10]
> 
> ```

Yes, this doesn’t work since `vec1` is a vector of vectors of vectors. It’s a bit awkward to perform broadcasting to arbitrarily deep levels, but if you dig around on this forum you might find some nice solutions. Here’s one, I don’t think it’s the best one:

```julia
myabs(x) = abs(x)
myabs(x::AbstractArray) = myabs.(x)

```

Then `myabs(vec1)` should work.

BTW FYI, you can make a matrix of `SVector`s, if you want:

```julia
julia> vec1 = [SVector{3,Float64}(-i,-i,-i) for i in 1:4, j in 1:10]
4×10 Array{SArray{Tuple{3},Float64,1,3},2}:
 [-1.0, -1.0, -1.0] [-1.0, -1.0, -1.0] [-1.0, -1.0, -1.0] … [-1.0, -1.0, -1.0] [-1.0, -1.0, -1.0]
 [-2.0, -2.0, -2.0] [-2.0, -2.0, -2.0] [-2.0, -2.0, -2.0] [-2.0, -2.0, -2.0] [-2.0, -2.0, -2.0]
 [-3.0, -3.0, -3.0] [-3.0, -3.0, -3.0] [-3.0, -3.0, -3.0] [-3.0, -3.0, -3.0] [-3.0, -3.0, -3.0]
 [-4.0, -4.0, -4.0] [-4.0, -4.0, -4.0] [-4.0, -4.0, -4.0] [-4.0, -4.0, -4.0] [-4.0, -4.0, -4.0]

```

Or maybe you could define a recursive helper function, like this:

```julia
rdot(f, x) = f(x) # recursive dot
rdot(f, x::AbstractArray) = rdot.(f, x)

```

and then call

```julia
vec2 = rdot(abs, vec1)

```

This would work with other functions beside `abs` as well.

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [March 27, 2020, 6:57pm UTC](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333/36 "2020-03-27T18:57:38Z")

</div>

**Results update**

It took me a couple of days to re-write the old allocation heavy code per the static array suggestions above, but got it working and here are the results:

Old Code (allocation heavy due to use of base array and slicing approach):

```julia
455.657741 seconds (585.01 M allocations: 73.176 GiB, 88.12% gc time)

```

Revised Code (making use of static vector/matricies (StaticArrays.jl))

```julia
23.023161 seconds (71.79 M allocations: 3.910 GiB, 66.59% gc time)

```

Between the above improvements and those discussed on a related thread that’s about a 1200X improvement to my original ugly code. Thanks to all who contributed comments and suggestions!

> [@Example C vs Julia vs Python](https://discourse.julialang.org/t/example-c-vs-julia-vs-python/33120):
>
> I have been using Julia for about a year. I am an EE, not a programmer. In one of my graduate EE courses I requested to use Julia instead of C for an assignment and got the task of comparing performance of C, Julia, and Python. As a test case I did a 3D finite difference time domain simulation of electromagnetic wave propagation. I compared timing and got the following results. More detail on the approach is below. Also the code is posted below. Any suggestions for how to speed things up …

[Previous page](https://discourse.julialang.org/t/variation-in-computation-time-for-execution-of-the-exact-same-code/36333.md?page=1)