# Why does this Python code performs three times faster than Julia?

**URL:** <https://discourse.julialang.org/t/why-does-this-python-code-performs-three-times-faster-than-julia/104919>\
**Category:** Performance\
**Tags:** performance\
**Created:** [October 13, 2023, 6:25am UTC](https://discourse.julialang.org/t/why-does-this-python-code-performs-three-times-faster-than-julia/104919 "2023-10-13T06:25:56Z")\
**Posts on this page:** 3\
**Page:** 2

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [October 13, 2023, 6:22pm UTC](https://discourse.julialang.org/t/why-does-this-python-code-performs-three-times-faster-than-julia/104919/21 "2023-10-13T18:22:36Z")

</div>

> [@SantiagoOrtiz](#):
>
> I thought that type annotation for variables could improve the performance. That is not the case in Julia, I guess? What is usually the main reason one would annotate the type of a local variable, then?

There is a relevant section of the manual (though it talks about function arguments rather than local variables): [Argument-type declarations](https://docs.julialang.org/en/v1/manual/functions/#Argument-type-declarations).

In practice, explicit local-variable declarations are rarely used in Julia. Mainly to have finer-grained control over type-promotion/conversion rules (as an alternative to explicit type conversion at each assignment) or to control variable scope. Not for performance per se.

---

<div class="post-metadata">

**Author:** ![ufechner7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ufechner7/32/51363_2.png) [@ufechner7](https://discourse.julialang.org/u/ufechner7)\
**Post date:** [October 13, 2023, 6:26pm UTC](https://discourse.julialang.org/t/why-does-this-python-code-performs-three-times-faster-than-julia/104919/22 "2023-10-13T18:26:04Z")

</div>

> [@SantiagoOrtiz](#):
>
> What is usually the main reason one would annotate the type of a local variable

If you define a struct or an array, then it is worth to use concrete types. Otherwise I don’t do it (unless for defining functions that shall only work on a specific type). In other words: For a local variable, don’t do it.

---

<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:** [October 13, 2023, 11:59pm UTC](https://discourse.julialang.org/t/why-does-this-python-code-performs-three-times-faster-than-julia/104919/23 "2023-10-13T23:59:17Z")

</div>

Since you are playing with a MD code, and I may guess you speak Spanish, you can find some elemental implementations here, to play with:

> **[GitHub - m3g/FundamentosDMC.jl: Material para el curso Fundamentos de...](https://github.com/m3g/FundamentosDMC.jl)**
>
> Material para el curso Fundamentos de Mecánica Estadística y Simulaciones: Simulación Computacional Avanzada en Química, Bioquímica y Ciencias de Materiales. CELFI/Universidad de Buenos Aires - Git...

Or other single - file examples here:

> <https://github.com/m3g/2021_FortranCon/blob/main/celllistmap_vs_namd/simulate.jl>

ps: This is how I would/could write that code snippet:

```julia
@kwdef struct Params{T}
    σ::T = 1.0
    ϵ::T = 1.0
    m::T = 1.0
    Δt::T = 1.0
end

import FastPow: @fastpow 
force(r, ϵ, σ) = @fastpow -4*ϵ*(-12*σ^12/r^13 + 6*σ^6/r^7)

function md_model2(ʋᵢ::T, rᵢ::T, n; params = Params{T}()) where {T}
    (; σ, ϵ, m, Δt) = params
    rʋ = zeros(T, 2, n)
    for i in axes(rʋ, 2)
        ʋᵢ = ʋᵢ + Δt/2*(force(rᵢ, ϵ, σ)/m)
        rᵢ₊₁ = rᵢ + Δt*ʋᵢ
        ʋᵢ₊₁ = ʋᵢ + Δt/2*(force(rᵢ₊₁, ϵ, σ)/m)
        rʋ[1, i] = rᵢ = rᵢ₊₁
        rʋ[2, i] = ʋᵢ = ʋᵢ₊₁
    end
    return rʋ
end

```

With which you get:

```julia
# original
julia> @btime md_model($ʋ₀, $r₀, $n);
  1.864 s (50000004 allocations: 839.23 MiB)

# code above
julia> @btime md_model2($ʋ₀, $r₀, $n);
  137.313 ms (2 allocations: 76.29 MiB)

```

Notes:

- Having the parameters in a struct gives you flexibility (to change them), without creating performance problems, because then it is practical to pass the parameters to the functions.
- Julia is column-major, so better store the consecutive data (coordinates, for instance) in a `(2,N)` matrix than in a `(N,2)` matrix. For coordinates, you would normally use the `StaticArrays` package.
- The `@fastpow` is important for performance in this case (gives you a ~5x speedup). But you can avoid it if you want by splitting the powers in products of small powers.

With these small powers and some additional small optimizations you can get:

```julia
julia> @btime md_model3($ʋ₀, $r₀, $n);
  112.524 ms (2 allocations: 76.29 MiB)

```

> **code**
>
> ```julia
> @kwdef struct Params{T}
> σ::T = 1.0
> ϵ::T = 1.0
> m::T = 1.0
> Δt::T = 1.0
> end
> 
> # ϵ4 = 4*ϵ
> # σ6 = σ^6
> function force3(r, ϵ4, σ6)
> r6 = (r^3)^2
> r7 = r6 * r
> r13 = r6 * r7
> return -ϵ4*(-12*σ6^2/r13 + 6*σ6/r7)  
> end
> 
> function md_model3(ʋᵢ::T, rᵢ::T, n; params = Params{T}()) where {T}
> (; σ, ϵ, m, Δt) = params
> ϵ4 = 4ϵ
> σ6 = σ^6
> Δt2 = Δt/2
> invm = 1/m
> rʋ = zeros(T, 2, n)
> for i in axes(rʋ, 2)
> ʋᵢ = ʋᵢ + Δt2*invm*(force3(rᵢ, ϵ4, σ6))
> rᵢ₊₁ = rᵢ + Δt*ʋᵢ
> ʋᵢ₊₁ = ʋᵢ + Δt2*invm*(force3(rᵢ₊₁, ϵ4, σ6))
> rᵢ = rᵢ₊₁
> ʋᵢ = ʋᵢ₊₁
> rʋ[1, i] = rᵢ₊₁
> rʋ[2, i] = ʋᵢ₊₁
> end
> return rʋ
> end
> 
> ```

Forgot to mention: with the above code you can also run with `Float32` and get some additional speedup:

```julia
julia> @btime md_model3(1.f0, 1.f0, $n; params = Params{Float32}());
  92.385 ms (2 allocations: 38.15 MiB)

```

[Previous page](https://discourse.julialang.org/t/why-does-this-python-code-performs-three-times-faster-than-julia/104919.md?page=1)
