Dot operator, 1. with array-indexing expression on the left-hand side and 2. with Matrix*Vector-multiplication

Hallo,
I would like to understand how the dot operator works. There are two specific problems:

at Performance Tips · The Julia Language (headline: Consider using views for slicing)
I learned that "…on the left-hand side of an assignment, where array[1:5, :] = ... assigns in-place to that portion of array". As I understand it, no unnecessary copies are stored in memory.

Now I read here:
Functions · The Julia Language that if I use the dot-operator .= :
"If the left-hand side is an array-indexing expression, e.g. X[begin+1:end] .= sin.(Y) , then it translates to broadcast! on a view, e.g. broadcast!(sin, view(X, firstindex(X)+1:lastindex(X)), Y) , so that the left-hand side is updated in-place. "

Now I wonder why the dot near the equal sign in X[begin+1:end] .= sin.(Y) is necessary. Shouldn’t it be updated in-place anyway? What is the advantage of the View in this case?

If I have an expression like: Θ[:,k+1] = K*Θ[:,k]+B
and K is a 10x10 matrix , Θ[:,k] is a 10x1 vector and B is another 10x1 vector.
Can I take advantage of the dot vector in this case?
I can’t use . * because I don’t want a broadcast. So I can’t use @.
Because of the array-indexing expression on the left-hand side the .= seems unnecessary (it’s in-place anyway). And .+ without fusion doesn’t seem to make sense.
Did I get that right?

It would be great if someone could explain this to me :flushed:

The difference between = and .= is whether you allocate an array for the answer before copying it to the view.

Is it wrong that an array-indexing expression at the left-hand side assignes in place by default?

This allocates first a slice Θ[:,k], then the product K*Θ[:,k], then for the sum + B before writing into Θ. Better would be Θ[:,k+1] .= K*view(Θ,:,k) .+ B, the .+ is useful precisely because it can be fused with the .=. Even better would be

Θ[:,k+1] .= B
@views mul!(Θ[:,k+1], K, Θ[:,k], true, true)
2 Likes

In case it wasn’t clear, you absolutely do need the dotted assignment operator here. Without the dot, the values are assigned ‘in-place’ in the variable X. But before that happens, the expression sin.(Y) on the right hand side will create a temporary array. If you want to avoid that temporary array you must dot the entire expression, including the assignment operator.

Without the dot, this

X[begin+1:end] = sin.(Y)

is equivalent to this

temp = sin.(Y)  # create temporary array
X[begin+1:end] = temp  # assign in-place to X
5 Likes

Thanks a lot macabbott for that really helpful answer!
But I’m still confused about
Θ[:,k+1] .= B
I’m not familiar with views. If I don’t assign a copy of B but a view, wouldn’t a subsequent change in Θ then also have to change B? I’ve tried it and seen that it doesn’t happen, but I still don’t understand how this works.

Maybe what you miss is that Θ is only one pointer deep. It owns a continuous block of memory. Or, it’s a mutable container of just numbers, not boxes which contain numbers. We can’t link part of it to B's memory, all we can do is copy the numbers one by one from B to each overwrite one in Θ. Which is what Θ[:,k+1] = B and Θ[:,k+1] .= B and Θ[:,k+1] .= view(B,:) all do, although the precise chain of functions involved will differ.

That makes it a lot clearer to me. Thank you!
Now I ask myself why I didn’t think of that :blush:

Sorry for being so persistent.
But I did

@time begin
       Θ[:,k+1] .= B
       @views mul!(Θ[:,k+1], K, Θ[:,k], true, true)
 end

and got 0.000030 seconds and 14 allocations: 400 bytes.

For
@time Θ[:,k+1] .= K*Θ[:,k] .+ B

I got 0.000029 seconds and 8 allocations: 560 bytes

These are more allocations than expected!
Is for efficiency only the number of allocated bytes critical or does every allocation cost time?

Don’t use the @time macro for micro benchmarks, it is not suitable. Install the BenchmarkTools package, read the manual, and then use that. The allocation estimates you get from @time do not tell you what you are looking for.

1 Like

To be honest, I wonder why @time is in Base. It seems like 99% of the time BenchmarkTools should be used instead. Perhaps it would be better if @time was removed, and then you would have to make a conscious choice when timing code, instead of relying on the built-in default.

2 Likes

Do you mean that?

@benchmark Θ[:,k+1] .= K*Θ[:,k] .+ B

BenchmarkTools.Trial:
memory estimate: 560 bytes
allocs estimate: 9

minimum time: 1.324 μs (0.00% GC)
median time: 1.373 μs (0.00% GC)
mean time: 1.415 μs (0.00% GC)
maximum time: 5.119 μs (0.00% GC)

samples: 10000
evals/sample: 10

@benchmark  begin
              Θ[:,k+1] .= B
              @views mul!(Θ[:,k+1], K, Θ[:,k], true, true)
              end

BenchmarkTools.Trial:
memory estimate: 400 bytes
allocs estimate: 15

minimum time: 2.229 μs (0.00% GC)
median time: 2.259 μs (0.00% GC)
mean time: 2.316 μs (0.00% GC)
maximum time: 5.337 μs (0.00% GC)

samples: 10000
evals/sample: 9

If I do

@benchmark  Θ[:,k+1] .= B

I would expect no allocations but I get:

BenchmarkTools.Trial:
memory estimate: 176 bytes
allocs estimate: 7

minimum time: 941.893 ns (0.00% GC)
median time: 961.571 ns (0.00% GC)
mean time: 1.002 μs (0.00% GC)
maximum time: 8.047 μs (0.00% GC)

samples: 10000
evals/sample: 28

And if I do

@benchmark  Θ[:,k+1] = B

BenchmarkTools.Trial:
memory estimate: 16 bytes
allocs estimate: 1

minimum time: 196.114 ns (0.00% GC)
median time: 197.932 ns (0.00% GC)
mean time: 204.880 ns (0.45% GC)
maximum time: 9.356 μs (97.79% GC)

samples: 10000
evals/sample: 621

I’m really confused about it :roll_eyes:

You must read the manual, specifically the part about variable interpolation.

I did, but I’m not sure what to do. Both Θ[:,k+1] and B are globals in my little trial.

I put it inside a function like

function dotequal(A,B)
       A.=B
       end

@benchmark dotequal(Θ[:,k+1],B)

but get the same result for the function with and the one without the dot.

I think I will put this aside for now.
Nevertheless, thank you very much for your patience.

Look up variable interpolation using $ in the BenchmarkTools manual.

1 Like

They mean you should use the interpolation operator $ to get proper timings:

julia> @benchmark $Θ[:,$k+1] .= $B
BenchmarkTools.Trial:
  memory estimate:  0 bytes
  allocs estimate:  0
  --------------
  minimum time:     21.866 ns (0.00% GC)
  median time:      21.966 ns (0.00% GC)
  mean time:        22.064 ns (0.00% GC)
  maximum time:     51.956 ns (0.00% GC)
  --------------
  samples:          10000
  evals/sample:     997

julia> @benchmark $Θ[:,$k+1] = $B
BenchmarkTools.Trial:
  memory estimate:  0 bytes
  allocs estimate:  0
  --------------
  minimum time:     15.315 ns (0.00% GC)
  median time:      15.315 ns (0.00% GC)
  mean time:        15.476 ns (0.00% GC)
  maximum time:     126.126 ns (0.00% GC)
  --------------
  samples:          10000
  evals/sample:     999

julia> @benchmark $Θ[:,$k+1] .= $K * $Θ[:,$k] .+ $B
BenchmarkTools.Trial:
  memory estimate:  320 bytes
  allocs estimate:  2
  --------------
  minimum time:     160.842 ns (0.00% GC)
  median time:      167.347 ns (0.00% GC)
  mean time:        179.127 ns (1.53% GC)
  maximum time:     868.750 ns (75.86% GC)
  --------------
  samples:          10000
  evals/sample:     784

julia> @benchmark  begin
           $Θ[:,$k+1] .= $B
           @views mul!($Θ[:,$k+1], $K, $Θ[:,$k], true, true)
       end
BenchmarkTools.Trial:
  memory estimate:  0 bytes
  allocs estimate:  0
  --------------
  minimum time:     98.835 ns (0.00% GC)
  median time:      99.894 ns (0.00% GC)
  mean time:        104.734 ns (0.00% GC)
  maximum time:     208.369 ns (0.00% GC)
  --------------
  samples:          10000
  evals/sample:     944
1 Like

Thank you so much!
I I’ve tried that too, but forgot the $ next to k+1 and it didn’t work.