Static matrix with size 400

New to Julia
1

Hi,

I want to write a function to generate random matrices with size 20x20. I need to use this function frequently but I don’t need to change the matrices generated by it. So the StaticArrays seem to be a good choice. But I further notice that it only works for arrays less than 100. Is there any way to overcome this?

Thanks

2

why?

julia> SMatrix{20,20}(rand(20,20))
20×20 SMatrix{20, 20, Float64, 400} with indices SOneTo(20)×SOneTo(20):

#or

julia> @SMatrix(rand(20,20))
20×20 SMatrix{20, 20, Float64, 400} with indices SOneTo(20)×SOneTo(20):

It is true that the performance benefit degrades as matrix gets larger but you should try it for your workload to determine if it actually speeds things up or not, I think 20x20 is still maybe possible to still see a benefit

1 Like
3

I found it in this discussion. Large matrix crashes StaticArrays. Maybe things changed. Thanks for this check. I will try it myself

1 Like
4

1000x1000 is quite different from 20x20

5

Compile and runtime performance may be bad.

3 Likes
6

What do you do with the matrices later? Allocating a single (standard) Matrix at the start of the procedure and then using rand! instead of rand might be efficient too

3 Likes
7

Hi

I basically want to generate some haar random unitary matrices. The code I used before is like

using StaticArrays
using LinearAlgebra
using StatsBase
using Random

function rand_haar2(::Val{n}) where n
    M = @SMatrix randn(ComplexF64,n,n) 
    q = qr(M).Q
    L = cispi.(2 .* @SVector(rand(Float64,n)))
    return q*diagm(L)
end

Then I want to frequently use the function for instance u = rand_haar2(Val{4}()) if I want random unitary matrices with size 4. This worked well before. But now I want to generate random unitary matrices with size 20. This seems to be not so efficient based on the discussion in the link I posted. I am using the version without static arrays.

function rand_haar2_new(::Val{n}) where n
    M = randn(ComplexF64,n,n) 
    q = qr(M).Q
    L = cispi.(2 .* (rand(Float64,n)))
    return q*diagm(L)
end

It makes my code a little faster. But I am thinking if there can be further improvement.

8

In the version without StaticArrays, you don’t need the n to be known statically, i.e. enclosed in a Val. This will only slow down compilation. Moreover, you will be better off reusing the memory you create, like so:

using LinearAlgebra
using Random: rand!, randn!

function rand_haar2!(A::AbstractMatrix, M::AbstractMatrix, L::AbstractVector)
    randn!(M)
    Q = qr!(M).Q
    rand!(L)
    copyto!(A, Q)
    transpose(A) .*= cispi.(2 .* L)  # double check that this does indeed multiply the columns of `A`
    return A
end

n = 20
M = randn(ComplexF64, n, n);
L = rand(Float64, n);
A = similar(M);

rand_haar2!(A, M, L)
3 Likes