Julia vs Matlab - building a Jacobian matrix

New to Julia
21

I see, very clear once again! I also guess from this that in the case of nested loops:

@inbounds for i in 1:10
              for j in 1:10
                  ... 
              end
           end

the function @inbounds do not apply for all loops and I should call it again when needed.
I did no know about the begin environment. Thank you for that as well!

1 Like
22

Another performance tip when using nested for loops:
Make sure the inner loop traverses the columns of your matrix (a.k.a. the first index).
This can make a dramatic difference, especially when your matrices are large and when the logic using them is simple.

See this example:

using BenchmarkTools

function f1(n)
    res = zeros(n, n)
    for i in 1:n, j in 1:n
        res[i, j] = i + j
    end
    res
end

function f2(n)
    res = zeros(n, n)
    for j in 1:n, i in 1:n
        res[i, j] = i + j
    end
    res
end

@btime f1(2000);
  58.354 ms (2 allocations: 30.52 MiB)

@btime f2(2000);
  12.590 ms (2 allocations: 30.52 MiB)
1 Like
23

I did not know about this! Thank you very much for the tip!

24

Just to make a point:

julia> sum_of_diff_3(x) = last(x) - first(x)
julia> @btime sum_of_diff_3($(rand(10_000)));
  1.849 ns (0 allocations: 0 bytes)

(this computes the same thing).

So the sum of diffs is a toy example, I just wanted to underscore that if you really want to speed things up, you look for a better algorithm. In many cases that will trump all tricks and optimizations. See if you are calculating things that you can do on paper, or maybe skip entirely.

4 Likes
25

Hi @rdeits
Your code here using btime says 0 allocations. If I use @time I get

@time sum_of_diff_2(rand(100))
  0.000006 seconds (6 allocations: 1.047 KiB)

Why the apparent difference?

26

You’re measuring the time and allocation of the construction of the rand(100) vector. In my example I used $ to interpolate that vector into the benchmarked expression: https://github.com/JuliaCI/BenchmarkTools.jl/blob/master/doc/manual.md#interpolating-values-into-benchmark-expressions

1 Like
27

Matlab also uses column-major ordering, so mtsch’s tip actually applies to both languages :slight_smile: (you’ll note that summing a matrix column in Matlab is about 5 times faster than summing an equivalent length row)

28

A note on @inbounds, I want to emphasize @rdeits suggestion that you not implement @inbounds in your code until you’re certain that it is correct. I cannot count the number of hours I have spent tracking down a stupid non-deterministic bug because I used @inbounds without knowing that in some cases it would index out of bounds and just grab the next thing available in memory. These are the hardest types of bugs to catch because their behavior can change each time you run the code depending on what is in memory “near” the thing you’re indexing into. Also, I’ve noticed that on V0.7 the benefits to @inbounds appear to be lower than in previous versions of Julia. I have no formal proof of this other than my own small amount of observational data.

4 Likes
29

There was a very nice PR that fixed this (https://github.com/JuliaLang/julia/pull/27386).

Note that you can turn off your @inbounds by starting julia with --check-bounds=yes in case you suspect that you have made a mistake somewhere.

4 Likes
30

It is very hard to say that, some times the numerical job depends on the third parts code (plotting etc)

31

Hey and welcome to the forum, please try not to “necropost” (post in very old threads) unless there is a strong reason to. It is often better to start a new thread than necroposting.

1 Like