Is there an elegant/optimized way to take the outer product of 3 sparse vectors?
Problem 1: Triple Outer Product
In numpy I can just chain calls to np.outer(np.outer(v,v),v) but I cannot do the same in Julia.
I see this example but that’s a bit gross and does not take advantage of the sparsity which leads to question 2.
Problem 2: Taking advantage of sparsity
How can I do this operation and take advantage of the fact that these are sparse vectors and as a result only a few of the multiplications actually need to be done. I can take the outer product of two sparse vectors fine but the taking the outer product of the sparse matrix and sparse vector does not seem to work.
Any help would be greatly appreciated!
Do you want a 3d sparse array as the result? The SparseArrays package does not support this, but there are other packages that do: Sparse 3D arrays - #7 by jovansam
Yeah basically just arbitrary dimension sparse arrays (4th is the most I’d need). There are a decent amount of threads asking for this functionality, is there any reason this could not or should not be added to SparseArrays in the Julia Base?
Like you mention there are a few packages out there, but none are well maintained. I will still try them though as I’d rather not have to implement this. Or maybe I’ll just use MATLAB.
Multidimensional sparse arrays are best stored in totally different data structures than sparse matrices and vectors, so it’s more natural to put them in a separate package IMO.
What’s wrong with e.g. SparseArrayKit.jl? It seems to be maintained (though not under heavy development), and indeed the author immediately responded to your query about a triple outer product.
Oh nothing seems wrong.
This is more my personal aspirations for packages in Julia. I know Julia is still relatively new, but one day I would like to see a lot of these smaller packages wrapped into consistently maintained and documented packages instead of having the sprawl of undocumented and dead packages we have right now. Some kind of happy medium between Base containing everything (e.g. C++ stdlib) and what we have now.
Guess I should mark a solution:
This can be done using SparseArrayKit.jl and TensorOperations.jl via:
N = 10
u = SparseArray{Float64}(undef,(N,))
@tensor sparse_triple_outer_product[i,j,k] := u[i] * u[j] * u[k]