# Julia is very slow compared to python/numpy. What am I doing wrong here?

**URL:** <https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598>\
**Category:** Performance\
**Tags:** benchmark, matrices\
**Created:** [March 28, 2022, 9:21am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598 "2022-03-28T09:21:12Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [March 28, 2022, 9:21am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/1 "2022-03-28T09:21:13Z")

</div>

Hi guys,  
I wrote a recursive function which construct a matrix in python using Qutip module:  
(I just used the tensor function and Pauli matrices definitions from Qutip)

```julia
import numpy as np
from qutip import *

def psitest(K):
    if K > 1 and K % 1 == 0:
        return tensor(psitest(K-1), - sigmaz())
    elif K == 1:
        return (1/np.sqrt(2)) * sigmax()

```

then for example for K=16 I have

```julia
Time=0.017930984

```

I wrote exactly the same function in Julia:

```julia
using LinearAlgebra
function ψtest(K)
    if K > 1 && K % 1 == 0
        return kron(ψtest(K-1), [-1. 0.;0. 1.])
    elseif K == 1
        return 1/√(2) * [0. 1.;1. 0.]
    end
end

```

where for the same K=16 gives

```julia
39.406 s (65 allocations: 42.67 GiB)

```

for K =25 actually python gives it in about 2 second where in Julia I can not even test it.  
I know I should not naively transform a code in python to Julia, but I am wondering what is here that makes this huge difference.  
Thank you all for your help

---

<div class="post-metadata">

**Author:** ![bkamins](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bkamins/32/208538_2.png) [@bkamins](https://discourse.julialang.org/u/bkamins)\
**Post date:** [March 28, 2022, 9:28am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/2 "2022-03-28T09:28:02Z")

</div>

I assume Qutip uses sparse arrays. You can do the same in Julia (like you - I have not tried to optimize code - just added `sparse`):

```julia
julia> using SparseArrays

julia> function ψtest(K)
           if K > 1 && K % 1 == 0
               return kron(ψtest(K-1), sparse([-1. 0.;0. 1.]))
           elseif K == 1
               return 1/√(2) * sparse([0. 1.;1. 0.])
           end
       end
ψtest (generic function with 1 method)

julia> @time ψtest(16);
  0.002063 seconds (182 allocations: 3.012 MiB)

julia> @time ψtest(25);
  1.169190 seconds (290 allocations: 1.500 GiB, 17.61% gc time)

```

with dense matrix `K=25` for sure will not fit in your RAM.

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [March 28, 2022, 9:42am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/3 "2022-03-28T09:42:49Z")

</div>

Wow, thank you very much 🙏. @bkamins.  
I appreciate if you also tell me other kind of optimizations I can do in this kind of computations.  
My problem is finally the SYK model in quantum mechanics.

---

<div class="post-metadata">

**Author:** ![bkamins](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bkamins/32/208538_2.png) [@bkamins](https://discourse.julialang.org/u/bkamins)\
**Post date:** [March 28, 2022, 9:52am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/4 "2022-03-28T09:52:59Z")

</div>

I am not a linear algebra expert, but assuming you want to keep the recursive definition then I think this is mostly OK. There might be a way to reduce the amount of allocations by pre-allocating the final matrix and using `kron!` on views, but it would need to be benchmarked.

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [March 28, 2022, 10:51am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/5 "2022-03-28T10:51:24Z")

</div>

I should check these. Thanks @bkamins again.

---

<div class="post-metadata">

**Author:** ![pitsianis](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pitsianis/32/26588_2.png) [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Post date:** [March 28, 2022, 4:59pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/6 "2022-03-28T16:59:20Z")

</div>

> [@Vahid\_Hosseinzadeh](#):
>
> I appreciate if you also tell me other kind of optimizations I can do in this kind of computations.

You are building a linear operator that you will eventually be applying to one or more vectors to transform them. If so, it will be much more effective to actually compute the result of applying the operator onto the data vector instead of the operator itself. That is, write a function that computes `ψtest(k)(v)` for any vector `v`.

For a general methodology to translate recursive operators that involve Kronecker products into code, please allow me to advertise my own work “The Kronecker Product in Approximation and Fast Transform Generation” [abstract only](https://dl.acm.org/doi/book/10.5555/266698) or the [full PhD pdf](https://users.cs.duke.edu/~nikos/reprints/T-001-KronProdApprox.pdf).

With Julia’s symbolic manipulation, term rewriting and introspection capabilities, I have been itching to redo the above in a single environment instead of the combination of conditional term rewriting in Mathematica, and the recurrence unrolling, expression simplification and C/FORTRAN/MATLAB code generation in Maple that I used 25 years ago…

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [March 30, 2022, 7:37pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/7 "2022-03-30T19:37:36Z")

</div>

Thanks @pitsianis for the explanations. I will look into your work for sure.  
My problem is coding for many body quantum systems like SYK model. Here for example the matrix finally is the Hamiltonian where I am interested in its eigenvalues. These \psi's are Majorana fermions which satisfy Clifford Algebra. Actually I do not know what I am really to compute 😄. For now I am just re-deriving what people have obtained.

---

<div class="post-metadata">

**Author:** ![pitsianis](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pitsianis/32/26588_2.png) [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Post date:** [March 31, 2022, 6:59am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/8 "2022-03-31T06:59:10Z")

</div>

> [@Vahid\_Hosseinzadeh](#):
>
> I am interested in its eigenvalues

The same advice holds for eigenvalues, since `eig(kron(A,B)) = kron(eig(A),eig(B))`

The right-hand-side expression is much cheaper to compute than the left.

In your specific recurrence example, it is even easier.

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [April 3, 2022, 8:51pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/9 "2022-04-03T20:51:58Z")

</div>

Thanks @pitsianis, this formula `eig(kron(A,B)) = kron(eig(A),eig(B))` is really cool and I should use it definitely. However in this part of my problem at the end I have H = \psi\_1 \psi\_2 \psi\_3 ... where \psi\_i = \sigma\_1 \otimes \sigma\_2 .... and I want the spectrum of H. Do you have any suggestion in this case? Beside that do you know any cheatsheet like references that I can find these kind of magic formulas?!  
thank you very much

---

<div class="post-metadata">

**Author:** ![pitsianis](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pitsianis/32/26588_2.png) [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Post date:** [April 7, 2022, 1:05pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/10 "2022-04-07T13:05:23Z")

</div>

> [@Vahid\_Hosseinzadeh](#):
>
> … do you know any cheatsheet like references that I can find these kind of magic formulas?

Start from pages 4-8 in the [PhD Thesis](https://users.cs.duke.edu/~nikos/reprints/T-001-KronProdApprox.pdf) that I pointed above. A good reference for me at the time was the book by Willi-Hans Steeb. “Kronecker Product of Matrices and Applications”, Wissenschaftsverlag, 1991.

I also recommend: Charles F. Van Loan, [The ubiquitous Kronecker product](https://doi.org/10.1016/S0377-0427(00)00393-9), Journal of Computational and Applied Mathematics, Volume 123, Issues 1–2, 2000, Pages 85-100

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [April 7, 2022, 1:29pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/11 "2022-04-07T13:29:32Z")

</div>

Thanks for the references. I start with your thesis.

---

<div class="post-metadata">

**Author:** ![MarcMush](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/marcmush/32/18006_2.png) [@MarcMush](https://discourse.julialang.org/u/MarcMush)\
**Post date:** [April 7, 2022, 2:34pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/12 "2022-04-07T14:34:30Z")

</div>

also, your function is not type-stable: it will return `nothing` when K is less than 1 or non integer

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [April 7, 2022, 2:38pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/13 "2022-04-07T14:38:39Z")

</div>

Thank you. yes, I should take care of that.

---

<div class="post-metadata">

**Author:** ![jlapeyre](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jlapeyre/32/4514_2.png) [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)\
**Post date:** [April 7, 2022, 2:56pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/14 "2022-04-07T14:56:59Z")

</div>

In case you are not aware [QuantumOptics.jl](https://github.com/qojulia/QuantumOptics.jl) has significant overlap with QuTip.

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [April 7, 2022, 3:05pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/15 "2022-04-07T15:05:08Z")

</div>

Yeah, I just add it and trying to learn it. However I am not quit sure that it is suitable for my kind of problem as I want some rare used algebras like Clifford algebra or others. Maybe it has it. Thanks for the suggestion though.  
Sorry a quick question: is numpy.complex128 is the same as Complex{Float64} in Julia. I think they are same but not sure. I am quit surprised that how much Julia is easy and fast by the way. I wrote the code naively in Julia and is 4-5x faster even without using any special package (I thought first that this might be due to precision).

---

<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:** [April 7, 2022, 3:06pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/16 "2022-04-07T15:06:11Z")

</div>

They are the same.

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [April 7, 2022, 3:07pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/17 "2022-04-07T15:07:05Z")

</div>

That is a relief. THanks

---

<div class="post-metadata">

**Author:** ![jlapeyre](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jlapeyre/32/4514_2.png) [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)\
**Post date:** [April 7, 2022, 3:54pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/18 "2022-04-07T15:54:20Z")

</div>

```julia
julia> ComplexF64
ComplexF64 (alias for Complex{Float64})

```

IIRC Julia called this `Complex128` , but it was changed before 1.0

---

<div class="post-metadata">

**Author:** ![alastair-marshall](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/alastair-marshall/32/21418_2.png) [@alastair-marshall](https://discourse.julialang.org/u/alastair-marshall)\
**Post date:** [April 8, 2022, 12:23pm UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/19 "2022-04-08T12:23:18Z")

</div>

You might want to check something like this, I don’t know if it’s good for what you want to do

[https://github.com/Krastanov/QuantumClifford.jl](https://github.com/Krastanov/QuantumClifford.jl)

---

<div class="post-metadata">

**Author:** ![Vahid\_Hosseinzadeh](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vahid_hosseinzadeh/32/33293_2.png) [@Vahid\_Hosseinzadeh](https://discourse.julialang.org/u/Vahid_Hosseinzadeh)\
**Post date:** [April 9, 2022, 7:59am UTC](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598/20 "2022-04-09T07:59:08Z")

</div>

Thanks for sharing this. It is not for my current job. But it is interesting. Maybe in the future.

[Next page](https://discourse.julialang.org/t/julia-is-very-slow-compared-to-python-numpy-what-am-i-doing-wrong-here/78598.md?page=2)