How to write generic low-rank updates?

General Usage
1

I have to make computations with a low-rank update of a matrix A+uv.

Is there a way to modify the statement

x = A*w + u * (v*w)

so that it works generically for all reasonable types? u,v could be a vector and its adjoint, or two matrices of compatible size. In Matlab it would work out of the box because “everything has two indices” there, but in Julia, as it is written, for instance it fails for

A = randn(5, 5)
u = randn(5)
v = abs.(randn(5)')
w = randn(5)

because v is an Array{Float64,2}, so the product in parentheses returns an Array{Float64,1} and the product outside parentheses fails with

ERROR: MethodError: no method matching *(::Array{Float64,1}, ::Array{Float64,1})

Ideally, x should be a Vector if w is a Vector, and a Matrix if w is a Matrix.

2

The problem here is that abs.(randn(5)') doesn’t have the same type as randn(5)' and returns a matrix instead of an adjoint. I guess that’s a problem with broadcasting which should preserve the input type here?

3

Note that map(abs, randn(5)') preserves the Adjoint-ness.

1 Like
4

Thanks to both! I am not expert enough to tell if this is a problem with broadcasting that should be fixed at the language level or if it’s just me doing it wrong. But in any case switching to map worked for my use case.