Converting a complex vector to a sparse diagonal array in Julia

Viewed 485

Trying to replicate a calculation from Matlab in Julia but am having trouble converting a single column complex array into a sparse diagonalized array for matrix multiplication.

Here is the Matlab code I am trying to replicate in Julia:

x*diag(sparse(y)) 

where x is of size 60,600000 and is of type: double, and y is size 600000,1 and is of type: complex double.

2 Answers

You can use Diagonal for that:

using LinearAlgebra
x=rand(6,60)
y=rand(Complex{Float64},60,1)
x*Diagonal(view(y,:))

I have used view(y,:) to convert y to a Vector - this is here a dimension drop operator you can also use shorter form vec(y) instead. Depending on what you want to do you might explicitly say that you want first column by view(y,:,1).

Note that Diagonal is just a sparse representation of a matrix.

julia> Diagonal(1:4)
4×4 Diagonal{Int64,UnitRange{Int64}}:
 1  ⋅  ⋅  ⋅
 ⋅  2  ⋅  ⋅
 ⋅  ⋅  3  ⋅
 ⋅  ⋅  ⋅  4

Another option that might cover more use case scenarios are BandedMatrices:

using BandedMatrices
x*BandedMatrix(0=>view(y,:))

Note that BandedMatrix uses set of pairs for bands, where band 0 is actually the diagonal.

I guess you don't mean it like this, but one can also interpret the question in the way that y is a sparse vector in the Julia sense, and you want to construct a sparse diagonal matrix out of it. In that case you can do the following:

julia> y = sprand(10, 0.2)
10-element SparseVector{Float64,Int64} with 2 stored entries:
  [4 ]  =  0.389682
  [5 ]  =  0.232429

julia> I, V = findnz(y)
([4, 5], [0.3896822408908356, 0.2324294021548845])

julia> sparse(I, I, V)
5×5 SparseMatrixCSC{Float64,Int64} with 2 stored entries:
  [4, 4]  =  0.389682
  [5, 5]  =  0.232429

Unfortunately, spdiagm does not preserve structural zeros for a sparse input:

julia> spdiagm(0 => y)
10×10 SparseMatrixCSC{Float64,Int64} with 10 stored entries:
  [1 ,  1]  =  0.0
  [2 ,  2]  =  0.0
  [3 ,  3]  =  0.0
  [4 ,  4]  =  0.389682
  [5 ,  5]  =  0.232429
  [6 ,  6]  =  0.0
  [7 ,  7]  =  0.0
  [8 ,  8]  =  0.0
  [9 ,  9]  =  0.0
  [10, 10]  =  0.0

I don't know whether this is intentional, but I filed an issue about this behaviour.

Related