Why for-loop of vector multiplication is slower that dot product?

General Usage
1

I just want to know why the .* below is faster that for-loop product?

using BenchmarkTools, LinearAlgebra
S1 = rand(10000);
S2 = rand(10000);
S = zeros(10000);
S3 = zeros(10000);
n =size(S1,1)
function mm!(S, S1, S2, n)
  for i in 1:n
    S[i] = S1[i]*S2[i]
  end  
end
julia> @btime $S .= $S1.*$S2;
  1.910 μs (0 allocations: 0 bytes)
julia> @btime mm!($S3, $S1, $S2, $n)
  2.089 μs (0 allocations: 0 bytes)
2

Your code seems to have some clear bugs in it that cause you to use global variables:

function mm!(S, S1, S3, n)
  for i in 1:n
    S[i] = S1[i]*S2[i]
  end  
end

You take in S, S1 and S3 as arguments, but operate on S, S1 and S2.

2 Likes
3

Sorry it was typo and updated it! Thanks for the correction.

4

add that @inbounds and you get the same performance, because the dot operation has that guaranteed.

That said, you can write that loop more cleanly as:

julia> function mm!(S, S1, S2)
         for i in eachindex(S, S1, S2)
            @inbounds S[i] = S1[i]*S2[i]
         end  
       end
mm! (generic function with 2 methods)

julia> @btime mm!($S3, $S1, $S2)
  1.706 μs (0 allocations: 0 bytes)

that eachindex(S, S1, S2) also guarantees that you are inbounds (it will error if the lengths of the arrays differ) - and, I bet, in the near future you won’t need, in this case, the @inbounds flag there to recover the performance of the dot operation.

4 Likes
5

Thank you very much for your comments. Why eachindex and @inbounds give better performance? Is it related to the compiler or the memory?

Where can I find the function of dot operation?

6
3 Likes
7

@inbounds disables bounds checking (which is disabled already on the dot operation, because the input guarantees that it is inbounds). By not having to check bounds, the compiler might be able (it is in this case) to use SIMD computations.

The eachindex does not improve performance here, but it will probably do in future Julia versions, by automatically guaranteeing the @inbounds inside the loop.

It is some broadcast operation, but that machinery is pretty complicated.

1 Like
8

Thank you. I got the

9

Also, it is not safe to use @inbounds with for i in 1:n, both because the arrays may not have 1-based indexing, and also because it is easy for a user-provided n to be wrong. It is not a good idea to pass the length of an array as a separate input argument, that can easily, quickly and safely be read directly from the array itself.

With eachindex you don’t need to worry about this.

2 Likes
10

How is about the index of matrix. For example:

n =4;
m=5;
s =zeros(n,m)
for i in 1:n
for j in 1:m
println(s[i,j]]
end
end
11

It’s similar. You can use eachindex with matrices too:

for i in eachindex(s)
    println(s[i])
end

or use axes:

for j in axes(s, 2)
    for i in axes(s, 1)  # first index should be inner loop, if possible
        println(s[i,j]]
    end
end
2 Likes
12

Got it. Thank you very much!