Create a custom Random.Sampler or Matrix sampler via Distribution.Sampleable

I am trying to define `rand` methods to “randomize” points `x` represented in a matrix `(n, s)` (`n` is the number of point and `s` is the dimension).
I could directly define

``````function rand(RM::RandomizationMethod)
x = RM.x
return shift(x)
end
``````

where `shift(x)` can be a random translation e.g. `shift(x) = x+rand(size(x))` and `RM` is a structure containing information on the original set of point to randomize (here directly the points but it can be their bit representation in some base and other properties).
All my different randomization methods are gathered in a type `RandomizationMethod`.

My questions:

• I have read about `Random.Sampler` but I am not sure to understand what they actually do (especially for nontrivial cases) and why one need them.
In my case are they necessary?
From what I understood, if I manage to define the appropriate `Sampler` then `rand` and all associate functions `rand!` or `rand(RM, 5)` (for 5 randomization so a `(n,s,5)` array) or `rand(RM, 5,2)` (for a `(n,s,5,2)` array) should work almost directly.
In particular, I could have a `Sampler` for one randomization and a slightly different if I want to generate multiple ones more efficiently.

• I saw Matrix Sampler from the `Distribution.jl` package. It seems to be doing sort of what I want with the `Sampleable` type.
Is it more appropriate for my problem?

• I am a bit lost on how to properly use that for matrix sampling: especially the part to specialize my code for multiple realization.
Here is a sort of MWE of my issue

``````using Random
using Distributions: Matrixvariate, Continuous
import Distributions: Sampleable
import Random: _rand!, rand!

abstract type SamplerShifting{F<:Matrixvariate,S<:Continuous} <: Sampleable{F,S} end
Base.size(sampler::SamplerShifting) = size(sampler.points) # the size of each matrix sample
struct Shifter{F<:AbstractMatrix{<:Real}} <: SamplerShifting{Matrixvariate,Continuous}
points::F
#... potentially other things
end

function shift!(rng, sampler, x)
for i in eachindex(x)
x[i] = sampler.points[i] + rand(rng)
end
end

function _rand!(rng::AbstractRNG, sampler::SamplerShifting, A::DenseMatrix{T}) where {T<:Real}
# ... generate a single matrix sample to A
#... potentially some initialization
shift!(rng, sampler, A)
return A #? is it not bad to return the value here for performance ?
end

function rand!(rng::AbstractRNG, s::Sampleable{Matrixvariate}, A::AbstractMatrix)
size(A) == size(s) ||
throw(DimensionMismatch("Output size inconsistent with sample length."))
_rand!(rng, s, A)
end

m = 5
n = 2^m
s = 10
x = zeros(n, s)
sp = Shifter(x)

rand(sp, 5, 2) #! this work! Returns a Matrix of my matrix samples

#* Now I want to generate mutliple samples to do all the initialization only once
#* so I try the following (which do not work)
function _rand!(rng::AbstractRNG, sampler::SamplerShifting, A::Array{DenseMatrix{T}}) where {T<:Real}
# ... generate a multiple matrix sample to A
#... potentially some initialization
for i in eachindex(A)
shift!(rng, sampler, view(A, i)) #? this in place replacment do not work proprely
end
return A
end

function rand!(rng::AbstractRNG, s::Sampleable{Matrixvariate}, A::Array{DenseMatrix{T}}) where {T<:Real}
size(A[firstindex(A)]) == size(s) ||
throw(DimensionMismatch("Output size inconsistent with sample length."))
_rand!(rng, s, A)
end

rand(sp, 5, 2) #!! Run but this does NOT use the new method for multiple call
A = [x, x]
rand!(sp, A) #! Run but produce the same realization
``````