Creating random Static Vector by StaticVector.jl and Distributions.jl in a loop

Hi
I wand to create a random 2 vector in a loop and want it to be static.
Here is the code:

using StaticArrays
using Distributions
dist2dN = MvNormal([0., 0.],[10. 0;0 10.])
for _ in 1:5
    x = @SVector [i for i in rand(dist2dN,1)]
    println(x)
end

where the result is not what I expected

[0.5365971833689593, -3.37550125378134]
[0.5365971833689593, -3.37550125378134]
[0.5365971833689593, -3.37550125378134]
[0.5365971833689593, -3.37550125378134]
[0.5365971833689593, -3.37550125378134]

without @SVector macro of course it generates different random vectors (actually not vectors but 2in1 matrices) each time.
Appreciate any help.
bests

1 Like

FWIW, this is because:

julia> @macroexpand @SVector [i for i in rand(dist2dN,1)]
quote
    #= /Users/chriselrod/.julia/packages/StaticArrays/0bweZ/src/SVector.jl:66 =#
    var"##292" = (i->begin
                #= /Users/chriselrod/.julia/packages/StaticArrays/0bweZ/src/SVector.jl:62 =#
                i
            end)
    #= /Users/chriselrod/.julia/packages/StaticArrays/0bweZ/src/SVector.jl:67 =#
    SVector{2}((var"##292"(3.9955890431901353), var"##292"(-2.5563441673280245)))
end

Maybe the easiest way would be

using StaticArrays, LinearAlgebra
S = @SMatrix[10.0 0.0; 0.0 10.0];# I'm assuming not diagonal in general?
m = @SVector zeros(2);
L = cholesky(S).L
for _ in 1:5
    x = m + L* @SVector(randn(length(m)))
    println(x)
end

I get

julia> for _ in 1:5
           x = m + L* @SVector(randn(length(m)))
           println(x)
       end
[3.747661043744402, -3.966264996727186]
[4.2503672520999265, 5.887575258272834]
[1.2384273143564963, 1.7567526130166946]
[3.092399770360624, -1.6927998659077437]
[-6.3010084683219025, -1.1313133350009374]

You can of course hard-code optimizations like not adding a zero mean if the mean is known to be 0.

Might be worth looking into a PR to Distributions.jl so that this returns an SVector:

julia> dist2dN = MvNormal(@SVector([0., 0.]),10.0)
MvNormal{Float64, PDMats.ScalMat{Float64}, SVector{2, Float64}}(
dim: 2
μ: [0.0, 0.0]
Σ: [100.0 0.0; 0.0 100.0]
)


julia> isbits(dist2dN)
true

julia> rand(dist2dN)
2-element Vector{Float64}:
 8.587630259531462
 1.4579629555675246

Thanks @Elrod for the explanations. Yes indeed, when the vector is static and not mutable, it is not strange that is not changing :sweat_smile:. I’ve come up with this though,

for _ in 1:5 
		x = rand(dist2dN)
		xss = @SVector [x[1],x[2]] 
		println(xss)
		println(isbits(xss))
	end
[-0.7533134459789248, 2.296833259646034]
true
[-3.4243299602495116, -2.7404408816986447]
true
[-1.020142875580128, -4.253875133369519]
true
[4.083471250760358, 1.2642846753752661]
true
[0.3937145507723281, -4.209228775050337]
true

It will be very good if Distributions.jl can create static objects itself.
thanks

That allocates a vector on every iteration.
My suggestion does not.
Presumably, if you’re using StaticAtrays, performance was the point.

You are right. But for my problem it is just enough (or maybe not!).
let me explain what is my problem. This 2-vector suppose to be the initial condition of a differential equation where I am using DifferentialEquations.jl to solve it. So I am just creating this SVector for each orbit and it improves benchmark (not very significantly actually).