using LogExpFunctions, Random
function _sample_1(θ)
u = rand(eltype(θ), size(θ))
return ifelse.(u .* (1 .+ exp.(-2θ)) .≤ 1, Int8(1), Int8(-1))
end
function _sample_2(θ)
u = randexp(eltype(θ), size(θ))
return ifelse.(2θ .≥ u .- log1mexp.(u), Int8(1), Int8(-1))
end
Here are two ways of sampling spins.
and their benchmarks
So implement the first.