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?

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

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

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.

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)