CODEWITHSUNDEEP
julia
RM-
RM-

Reputation: 1008

Vector matrix element wise multiplication (by rows) in Julia, efficiently

I want to multiply each row of the matrix A by the same vector v. For example

A =[1.0 3.0; 1.0 1.0]
v = [1.0, 2.0]

I would like to output

 [1.0 6.0; 1.0 2.0]

So far I am doing:

(v.*A')'

But I doubt this is computationally efficient given that I transpose a matrix twice.

As a note this question was answered for Matlab (https://uk.mathworks.com/matlabcentral/answers/243307-vector-matrix-multiplication-row-wise)

Upvotes: 5

Views: 4344

Answers (2)

carstenbauer
carstenbauer

Reputation: 10127

You have at least the following options:

  • (v.*A')' (OP's suggestion)
  • v'.*A (shortest way)
  • mapslices(row->v.*row, A, 2)
  • Manual implementation from @AborAmmar's post (fastest way)

i.e.

function tt(v, A)
    r = similar(A)
    @inbounds for j = 1:size(A,2) 
        @simd for i = 1:size(A,1) 
            r[i,j] = v[j] * A[i,j] # fixed a typo here!
        end
    end 
    r
end

Speed comparison (in ascending order)

julia> @btime tt($v,$A); # @AborAmmar's answer
  38.826 ns (1 allocation: 112 bytes)

julia> @btime ($v)'.*$A;               
  49.920 ns (1 allocation: 112 bytes)

julia> @btime (($v).*($A)')';           
  123.797 ns (3 allocations: 336 bytes)                           

julia> @btime mapslices(row->($v).*row, $A, 2);
  25.174 μs (106 allocations: 5.16 KiB)

EDIT: Took more careful and faster manual implementation from @AborAmmar's post (replaced old one) and updated speed comparison.

Upvotes: 5

AboAmmar
AboAmmar

Reputation: 5559

What is really awesome about Julia is that you can create your own functions with the same performance as highly optimized builtin functions. Here is another implementation with the same performance (or slightly better) as (v' .* A):

function tt1(v, A)
    v'.*A
end

function tt(v, A)
    r = similar(A)
    @inbounds for j = 1:size(A,2) 
        @simd for i = 1:size(A,1) 
            r[i,j] = v[j] * A[i,j]
        end
    end 
    r
end   

const n = 10^4
const A = rand(n,n)
const v = rand(n)

@time tt1(v, A);
@time tt(v, A);

  0.283537 seconds (6 allocations: 762.940 MiB, 19.70% gc time)
  0.235699 seconds (6 allocations: 762.940 MiB, 0.44% gc time)

Upvotes: 3

Related Questions