Could you suggest a good binary operator that is not commonly used for matrices in Julia Base or packages? I need to use it for some common matrix operations in my code. I don’t want to step on something else.

what’s the operation?

I don’t think there’s much problem with conflicts for everything beyond standard arithmetic. I personally like `⊕`

(`\oplus`

) or '⊗` (`

\otimes`).

If there is a conflict Julia just informs and asks to make the statement precise, e.g `MyLibrary.:⊕`

multiplication of a homogeneous matrix with an inhomogeneous vector, e.g.,

hmul(A,x) = A[1:end,1:end-1]*x.+A[:,end]

This is somewhat like a affine (translation) operation.

Some possibilities

`\otimesrhrim`

: `⨵`

`\oplusrhrim`

: `⨮`

since they have a sense of directionality to them…

Yes, homogeneous translation would be one application of it. Thanks for the tip!

is the precedence the same as * for those operators?

Seems to be faster to do

```
hmul2(A,x) = A * [x; 1]
```

Also, check out StaticArrays for this, it will be *much* faster, if you are working on e.g. 3-4D vectors.

```
julia> using StaticArrays
julia> A = @SArray rand(4,4);
julia> x = @SArray rand(3);
julia> @btime hmul($A, $x)
110.808 ns (3 allocations: 352 bytes)
4-element SizedVector{4, Float64, Vector{Float64}} with indices SOneTo(4):
1.281701572863999
1.45284156801398
0.8569634443644434
1.3118854951752263
julia> @btime hmul2($A, $x)
2.900 ns (0 allocations: 0 bytes)
4-element SVector{4, Float64} with indices SOneTo(4):
1.281701572863999
1.45284156801398
0.8569634443644434
1.3118854951752263
```

Wow, thanks for the tip. Yes, I am using 3-4D vectors and StaticArrays. I’m impressed that the augmentation does not cost an allocation. I see that it is constructed on the stack as an SVector. Very cool!

You don’t normally get any allocations as long as the sizes are statically known. It is easy for the compiler to see that `[x; 1]`

will have type `SVector{4, Float64}`

, though you could also use `StaticArrays.push`

with no difference.

If you want to get rid of the last element of the output vector (which is common), you can use `StaticArrays.pop`

. This can be quite useful, since

```
julia> typeof(x[1:end-1])
Vector{Float64} (alias for Array{Float64, 1})
julia> typeof(pop(x))
SVector{2, Float64} (alias for SArray{Tuple{2}, Float64, 1, 2})
```

Yes, type stability for these operations is tricky. Thanks for the tip. I have a more general question.

Suppose I have a vector of static vectors, e.g., something like,

```
100-element reinterpret(SVector{2, Float64}, ::Vector{Float64}
```

but I’d like to be able to do

```
y = A*x
```

What is the Julian and most performant way to accomplish this? E.g., this works, but seems a bit awkward,

```
y = Ref(A).*x
100-element Vector{SVector{2, Float64}}
```

Now suppose I want the output y to be a matrix. Tullio is an option, e.g.,

```
@tullio y[i,k] := A[i,j]*x[k][j]
2×100 Matrix{Float64}
```

Is there a solution using reinterpret?

This is the usual way.

That being said, if you are worried about the performance of this operation then you may be making a more fundamental mistake — the Matlab/Python style of “vectorized” programming leaves a lot of performance on the table by doing a bunch of separate “vectorized” loops that each perform one small operation (producing a bunch of temporary arrays) rather than combining the calculations into a single pass over the data whenever possible. Presumably, in your real application, you aren’t *just* doing `Ref(A) .* x`

, and you should try to exploit that by writing a *single* loop that does all of the operations you need on each `x[k]`

at once.

I suppose that loop fusion should eliminate temporaries for element-wise matrix operations. (matrix multiplication clearly cant be fused). Thanks for the response. I find some tension going between vector of vectors representations and matrices (or more generally tensors).

Just using a loop should be a good solution here, but you could also have a look at GitHub - JuliaArrays/HybridArrays.jl: Arrays with both statically and dynamically sized axes in Julia which lets you mix static and dynamic dimensions, e.g. the number of rows is static, while the number of columns is dynamic. I believe it builds on top of StaticArrays.

```
julia> As = @SMatrix rand(4,4);
julia> Xs = [@SVector rand(3) for _ in 1:10];
julia> reinterpret(reshape, Float64, hmul2.(Ref(As), Xs))
4×10 reinterpret(reshape, Float64, ::Vector{SVector{4, Float64}}) with eltype Float64:
1.68023 1.41434 1.33775 1.3449 2.34028 0.655653 2.2069 2.09013 1.45916 1.42397
2.12982 1.87084 1.80543 1.79198 2.79192 1.03538 2.60938 2.55969 1.67285 1.79207
1.11249 1.06255 1.02429 1.0432 1.2869 0.845418 1.24341 1.22446 1.01764 1.02749
1.09398 0.907248 0.901433 0.866234 1.47011 0.479292 1.39198 1.32744 0.951519 0.951254
```

The reference operator precedence table is here: