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

New to Julia
, ,
1

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
2

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
3

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

4

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

5

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).