# Matrix multiplication - inconsistent behaviour

**URL:** <https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640>\
**Category:** General Usage\
**Tags:** question\
**Created:** [March 28, 2020, 2:41pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640 "2020-03-28T14:41:34Z")\
**Posts on this page:** 11\
**Page:** 3

<div class="post-metadata">

**Author:** ![kristoffer.carlsson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/kristoffer.carlsson/32/22_2.png) [@kristoffer.carlsson](https://discourse.julialang.org/u/kristoffer.carlsson)\
**Post date:** [May 15, 2020, 10:00am UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/41 "2020-05-15T10:00:33Z")

</div>

Yeah, for people interested in this, the way forward would be to create a package with an array type that does not drop dimensions and implement the required methods on that. Then everyone who wants that behaviour just uses that package.

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [May 15, 2020, 1:19pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/42 "2020-05-15T13:19:51Z")

</div>

I fear that’d just replace one pitfall with an endless stream of other pitfalls. Significantly better if we could locally shadow indexing like TwoBasedIndices used to do or to use a macro — those would actually work like the folks here want.

But I really encourage folks who are transitioning to Julia to give this some time. I used to think we might need some special syntax to very tersely and efficiently create `[i]` for this purpose, but this is the first big thread about this behavior in years.

---

<div class="post-metadata">

**Author:** ![Henrique\_Becker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/henrique_becker/32/15443_2.png) [@Henrique\_Becker](https://discourse.julialang.org/u/Henrique_Becker)\
**Post date:** [May 15, 2020, 1:27pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/43 "2020-05-15T13:27:01Z")

</div>

> [@mbauman](#):
>
> his is the first big thread about this behavior in years

_cough [my thread had 24 posts and it is from this year](https://discourse.julialang.org/t/eachrow-over-array-any-2-does-not-return-rows-why/34601) cough_

What? No, I said nothing.

---

<div class="post-metadata">

**Author:** ![Henrique\_Becker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/henrique_becker/32/15443_2.png) [@Henrique\_Becker](https://discourse.julialang.org/u/Henrique_Becker)\
**Post date:** [May 15, 2020, 1:57pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/44 "2020-05-15T13:57:25Z")

</div>

So, just for me to understand better, what you recommend is something more like [`UnsafeArrays.jl`](https://github.com/JuliaArrays/UnsafeArrays.jl) strategy? There you can do:

```julia
@uviews A begin
    # inside this block A is replaced by an view object to the outer scope A
    # so changes to this A reflect in the outer A, but making slices of A
    # inside this block do not copy nor allocate, they are just pointers.
end

```

In this case, inside the block `A` would be replaced too by a `view`-like object that just wraps and modifies original `A` but the wrapper type do not drop dimensions? And when the scope ends the `A` is not the wrapped object so it already have the dimensions dropped if it was the case.

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [May 15, 2020, 2:07pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/45 "2020-05-15T14:07:46Z")

</div>

No, it’d be more akin to the `@views` macro and would walk through your expression tree to replace all `getindex` calls with a `getindex_preserving_dimensions` function — just like how `@views` swaps in `view` for `getindex`.

There’s been a slow but steady effort to make Julia’s syntax more locally overridable — for example in 1.5 you’ll be able to lexically change what gets called upon `A'` by locally shadowing `var"'"`:

```nohighlight
julia> module Test
       const var"'" = x->x^2
       a = 2
       @show a'
       end
a' = 4
Main.Test

```

This isn’t (yet) done in the same manner for square brackets. One tricky consideration is the fact that square brackets are used for both indexing and array literals.

Thanks for the pointer to your previous thread — I can no longer keep up with all the discussions here!

---

<div class="post-metadata">

**Author:** ![Henrique\_Becker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/henrique_becker/32/15443_2.png) [@Henrique\_Becker](https://discourse.julialang.org/u/Henrique_Becker)\
**Post date:** [May 15, 2020, 2:26pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/46 "2020-05-15T14:26:56Z")

</div>

> [@mbauman](#):
>
> I can no longer keep up with all the discussions here!

Completely understandable. Also, I do not know if 24 posts is considered long, XD.

Julia 1.5 seems to be very interesting. I am waiting for the non-allocating views that will make `UnsafeArrays.jl` not necessary anymore.

I proposed the solution of `UnsafeArrays.jl` strategy just because I found it was a little more accessible to someone with less experience, as the wrapper type does the heavy lifting and wrapper types are a more commonly used/understood concept. The macro doing the heavy lifting has the good side effect of not needing a wrapper type, nor black non-allocating structs magic, basically no runtime overhead, just needs someone that has a little more experience (I, for one, have never written a macro, even if I make widespread use of them).

---

<div class="post-metadata">

**Author:** ![jne](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jne/32/6341_2.png) [@jne](https://discourse.julialang.org/u/jne)\
**Post date:** [May 15, 2020, 3:40pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/47 "2020-05-15T15:40:56Z")

</div>

@compleat, what you said in your 3rd post still catches my eye:

> [@compleat](#):
>
> If I have a situation where I have an indices ‘ind1’ and ‘ind2’
> 
> mat1[ind1,:]\*mat2[:,ind2]

> [@compleat](#):
>
> I don’t want to always have  
> to anticipate the special cases of one or both of them being of length 1.

I’d like to make sure, by specific example, that by now you’ve caught the little gem that the existing array behavior serves you well on this:

```julia-repl
julia> mat1=randn(3,3); mat2=randn(3,2);

julia> ind1=[1,3]; ind2=[1,2]; mat1[ind1,:]*mat2[:,ind2];

julia> ind1=[1,3]; ind2=[2]; mat1[ind1,:]*mat2[:,ind2];

julia> ind1=[1]; ind2=[2]; mat1[ind1,:]*mat2[:,ind2];

```

(Note each succeeded without error. i.e. a vector of indices is precisely needed, even if length 1.)

Once you get used to the wisdom of this, you might find it natural to replace `[1]` with `1:1`.

If unconvinced, I hope at least I’ve reached another reader of this thread coming from Matlab. Best of luck-

---

<div class="post-metadata">

**Author:** ![jne](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jne/32/6341_2.png) [@jne](https://discourse.julialang.org/u/jne)\
**Post date:** [May 15, 2020, 9:14pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/48 "2020-05-15T21:14:14Z")

</div>

Hey @mbauman, while hot on the grill maybe you could teach me as a user what I should consider trade-off (if any) between using `[i]` and `i:i`. Trivial, yeah.

It occurs to me `i:i` may consume double the itty bitty space and not even iterate any more efficiently. I don’t know. But I assume they should interact identically with other ranges and functions of vectors, etc, or at least aspirationally. Otherwise just style? I can see `[i]` might look less readable when inside square brackets itself, as is often the case, but OTOH `[i]` is appealing in not violating the Don’t Repeat Yourself principle, as in `[really_long_name ___or_worse___ expression]`.

I’m also striving out of curiosity to understand what you meant by

> [@mbauman](#):
>
> I used to think we might need some special syntax to very tersely and efficiently create `[i]`

---

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [May 15, 2020, 9:44pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/49 "2020-05-15T21:44:29Z")

</div>

The biggest distinction between `[i]` (a `Vector`) and `i:i` (a `UnitRange`) is simply that the `Vector` is mutable. This has huge implications in both how it behaves and the sorts of optimizations Julia can easily perform. As an example, writing `[1]` in two separate places will create two separate `Vector`s because they might change!

```nohighlight
julia> A = [1]; B = [1];

julia> A == B
true

julia> A === B
false

julia> push!(A, 2)
2-element Array{Int64,1}:
 1
 2

julia> A == B
false

```

Every time you create a `Vector`, Julia needs to be prepared for its value and/or its length to possibly change. More succinctly, every `Vector` you create has a unique identity. `UnitRanges`, on the other hand, are _immutable_. They cannot change their values, so there’s no behavior _at all_ that Julia needs to track. In fact, Julia cannot possibly tell the difference between two separate creations of `1:1`. It is entirely defined by the _value_ of its two endpoints.

So, ironically, while it “feels” like `1:1` _should_ take up twice as much space as `[1]`, it actually ends up taking _no space at all_ in many common situations as Julia doesn’t even need to construct it in the first place. Further, the `Vector` has a header that it uses to track the array’s length and such. And so you end up seeing something like this in practice:

```nohighlight
julia> @allocated [1]
96

julia> @allocated 1:1
0

```

I think `x:x` that is the most succinct way of currently writing an immutable one-vector, but you’re right it’s not DRY at all. In the future Julia might get smarter about identifying uses of array literals that _happen to not need_ the mutable behaviors and thus give them the same performance, but that’s a harder thing to identify than when it’s promised up front.

---

<div class="post-metadata">

**Author:** ![jne](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jne/32/6341_2.png) [@jne](https://discourse.julialang.org/u/jne)\
**Post date:** [May 15, 2020, 10:31pm UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/50 "2020-05-15T22:31:56Z")

</div>

Ah yes, I again overlooked the key mutability aspect. I was especially entertaining examples where the `[i]` is used without explicitly storing, e.g. `mat1[[1],:]*mat2[:,[2]]`, so it’s interesting to hear mutability then too has a strong overhead.

Yes, in terms of DRY I’m myself tempted to write `range(i, length=1)` instead of `i:i` when they arise, or of course `range(an_expression, length=1)` instead of `i:i` preceded by `i=an_expression`. Perhaps an infix will someday arise for `range(i, length=n)` since it may be argued to be as natural as `i:j`. (And in that case, spelled-out `range` usages call that infix as they do (:), which would seem to result in a low overhead infix operator.)

Julia getting “smart” about array literals sounds like Matlab; and I don’t mean that flatteringly…

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [May 16, 2020, 6:06am UTC](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640/51 "2020-05-16T06:06:08Z")

</div>

> [@jne](#):
>
> Perhaps an infix will someday arise for `range(i, length=n)` since it may be argued to be as natural as `i:j` .

If you use these constructs often, nothing prevents you from defining a function as a shorthand. Eg

```julia
range_from(i, n) = range(i; length = n)

```

or similar, or

```julia
singleton_range(i) = i:i

```

You can of course use a short function name, or possibly Unicode symbols, to make these succinct.

> [@jne](#):
>
> Julia getting “smart” about array literals sounds like Matlab; and I don’t mean that flatteringly…

I don’t understand what you mean here; the only array literals I see mentioned in this topic are `[]` for constructing vectors — nothing special is going on there.

[Previous page](https://discourse.julialang.org/t/matrix-multiplication-inconsistent-behaviour/36640.md?page=2)