Transpose operator

General Usage
1

Hello everyone,

I am trying to use a transpose operator over a vector in order to perform am element-wise addition.
For example, I want to add a column vector a = [a1;a2;a3] to a row vector b = [b1,b2] I should get matrix
M = a+b = [a1+b1, a1+b2; a2+b1, a2+b2; a3+b1, a3+b2].
In MATLAB it is equivalent (if both vectors are row vectors) M = a.'+b

I am trying to get the same in Julia but here is the problem, there is no .’ operator in Julia starting from 1.0 version. There is the transpose operator which does not work in broadcasting mode. The adjoint operator is not Valid for me because I work with complex numbers.

a = Vector{ComplexF64}([1+3im,2])
b = Vector{ComplexF64}([0,0,0])
Z = zeros(ComplexF64,3,2)
G = zeros(ComplexF64,3,2)


@. Z = b + a'           # Works but takes the complex conjugate
@. Z = b + transpose(a) # DOES NOT WORK!!!! The error is " DimensionMismatch("array could not be broadcast to match destination") "
Z = b .+ transpose(a)   # Works but not efficient
@. Z = b + conj(a')

The third case Z = b .+ transpose(a) is not efficient because it makes 2 loops first one for addition b .+ transpose(a), than it runs the second loop one for the assignment of b .+ transpose(a) to Z. While the other 3 cases do it within one loop.
So which is the fastest way?
And why transpose doesn’t within Broadcasting?

Thank you in advance

2 
Z .= b .+ transpose(a)
3 Likes
3

To expand on @DNF s answer: The @. macro is broadcasting every function call, so the function transpose is not applied to the whole array, but is broadcasted and applied to every entry of a on its own and therefore does nothing. In your case, you don’t even need to preallocate Z at all. If you just do Z = b .+ transpose(a) without initializing Z first, it should be just as, if not a bit more efficient than initializing Z with zeros and modifying it in place.

4 Likes
4

Using broadcasted assignment is good when you already have a predefined array Z, but no point otherwise.

3 Likes
5

Thank you, that works!
It’s actually as fast as @. Z = b + conj(a') .

6

Thank you for your comment.
Now I understand why @. Z = b + transpose(a).This is equivalent to transpose.(a) which does nothing because each element of the vector “a” is just a number, so there is nothing to transpose!|
I better understand now how the macro @. operator works.

This is also true, without prealocation it is just 5 times slower. But in my case it is very critical because I am running an optimization. Thus 5 times speedup means to me waiting 1 day instead of 5 :slight_smile:

7

You can also do @. Z = b + $transpose(a) to prevent the transpose call from being “dotted”.

2 Likes
8

Thank you. It is good to know that operator $ can prevent element-wise operation.
I checked it out. That works but a bit slower than what @DNF suggests.

9

I think they should be completely identical, the difference is probably an issue with the benchmarking.

2 Likes
10

BTW, does anyone know why the ' operator apparently isn’t broadcasted? Are there any other examples of operators that are ignored by the @. macro?

11

Maybe I am missing something, but I don’t think ' is treated specially by @. (aka __dot__):

julia> v = fill(1+2im, 3)
3-element Array{Complex{Int64},1}:
 1 + 2im
 1 + 2im
 1 + 2im

julia> @. (v + 3)'
1×3 Adjoint{Complex{Int64},Array{Complex{Int64},1}}:
 4-2im  4-2im  4-2im

Also, it should be noted that transpose is recursive, and for the use case here it may be better style to use permutedims.

12

' is a special case, because it isn’t actually a regular function call, but instead gets its own AST representation. It is then lowered to a function call at the lowering stage, but because @. works on expressions, not on lowered IR, it doesn’t see ' as a function call and therefore doesn’t add any dots. I actually proposed making ' just a regular function call in RFC: parse `a'` as call expression by simeonschaub · Pull Request #33683 · JuliaLang/julia · GitHub, but unfortunately, that’s a bit breaking.
Another example of such an operator would be <::

julia> @. [Int, UInt, Float64] <: Integer
ERROR: TypeError: in <:, expected Type, got Array{DataType,1}
Stacktrace:
 [1] top-level scope at REPL[4]:1
1 Like