Skip to content

Commit

Permalink
reapply formatter
Browse files Browse the repository at this point in the history
  • Loading branch information
ArnoStrouwen committed Feb 22, 2024
1 parent 54d803e commit 47f54c7
Show file tree
Hide file tree
Showing 18 changed files with 414 additions and 397 deletions.
3 changes: 2 additions & 1 deletion .JuliaFormatter.toml
Original file line number Diff line number Diff line change
@@ -1,2 +1,3 @@
style = "sciml"
format_markdown = true
format_markdown = true
format_docstrings = true
2 changes: 1 addition & 1 deletion README.md
Original file line number Diff line number Diff line change
Expand Up @@ -197,7 +197,7 @@ names = [
"Shift",
"Digital Shift",
"Matousek Scramble",
"Owen Scramble",
"Owen Scramble"
]
p = [plot(thickness_scaling = 1.5, aspect_ratio = :equal) for i in sequences]
for (i, x) in enumerate(sequences)
Expand Down
2 changes: 1 addition & 1 deletion docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@ makedocs(sitename = "QuasiMonteCarlo.jl",
warnonly = [:missing_docs],
linkcheck_ignore = [
"https://www.sciencedirect.com/science/article/pii/S0010465509003087",
"https://artowen.su.domains/mc/qmcstuff.pdf",
"https://artowen.su.domains/mc/qmcstuff.pdf"
],
format = Documenter.HTML(assets = ["assets/favicon.ico"],
canonical = "https://docs.sciml.ai/QuasiMonteCarlo/stable/"),
Expand Down
2 changes: 1 addition & 1 deletion docs/pages.jl
Original file line number Diff line number Diff line change
Expand Up @@ -5,5 +5,5 @@ pages = [
"samplers.md",
"randomization.md",
"design_matrix.md",
"types.md",
"types.md"
]
3 changes: 2 additions & 1 deletion docs/src/design_matrix.md
Original file line number Diff line number Diff line change
Expand Up @@ -48,7 +48,8 @@ Using `std(μ)` then gives you the estimated variance of your RQMC prediction.

```@example 2
# Or using `generate_design_matrices`. Note that this is less memory efficient since it allocate space for 5 large big matrices.
μ = [mean(f(c) for c in eachcol(X)) for X in QuasiMonteCarlo.generate_design_matrices(N,
μ = [mean(f(c) for c in eachcol(X))
for X in QuasiMonteCarlo.generate_design_matrices(N,
d,
FaureSample(R = OwenScramble(base = b, pad = pad)),
num_mats)]
Expand Down
2 changes: 1 addition & 1 deletion docs/src/randomization.md
Original file line number Diff line number Diff line change
Expand Up @@ -143,7 +143,7 @@ names = [
"Shift",
"Digital Shift",
"Matousek Scramble",
"Owen Scramble",
"Owen Scramble"
]
p = [plot(thickness_scaling = 1.5, aspect_ratio = :equal) for i in sequences]
for (i, x) in enumerate(sequences)
Expand Down
18 changes: 9 additions & 9 deletions ext/QuasiMonteCarloDistributionsExt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -22,9 +22,9 @@ Return a point set from a distribution `D`:
- `ub` is the upper bound. Its dimension must match `length(lb)`.
"""
function QuasiMonteCarlo.sample(n::Integer,
d::Integer,
D::Distributions.Sampleable,
T = eltype(D))
d::Integer,
D::Distributions.Sampleable,
T = eltype(D))
@assert n>0 QuasiMonteCarlo.ZERO_SAMPLES_MESSAGE
x = [[rand(D) for j in 1:d] for i in 1:n]
return reduce(hcat, x)
Expand All @@ -51,8 +51,8 @@ Return a QMC point set where:
In the first method the type of the point set is specified by `T` while in the second method the output type is inferred from the bound types.
"""
function QuasiMonteCarlo.sample(n::Integer, lb::T, ub::T,
S::D) where {T <: Union{Base.AbstractVecOrTuple, Number},
D <: Distributions.Sampleable}
S::D) where {T <: Union{Base.AbstractVecOrTuple, Number},
D <: Distributions.Sampleable}
QuasiMonteCarlo._check_sequence(lb, ub, n)
lb = float.(lb)
ub = float.(ub)
Expand All @@ -61,10 +61,10 @@ function QuasiMonteCarlo.sample(n::Integer, lb::T, ub::T,
end

function QuasiMonteCarlo.DesignMatrix(N,
d,
D::Distributions.Sampleable,
num_mats,
T = Float64)
d,
D::Distributions.Sampleable,
num_mats,
T = Float64)
X = QuasiMonteCarlo.initialize(N, d, D, T)
return QuasiMonteCarlo.DistributionDesignMat(X, D, num_mats)
end
Expand Down
6 changes: 3 additions & 3 deletions src/Faure.jl
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ end

# raise a Pascal matrix to a given power.
@fastmath function pascal_power!(result::UpperTriangular, pascal::UpperTriangular,
power::Integer, base::Integer)
power::Integer, base::Integer)
@inbounds @simd for idx in eachindex(pascal)
i, j = Tuple(idx)
i j ? (result[idx] = powermod(power, j - i, base) * pascal[idx]) : result[idx]
Expand Down Expand Up @@ -45,7 +45,7 @@ Base.@kwdef struct FaureSample <: DeterministicSamplingAlgorithm
end

@fastmath function sample(n::Integer, dimension::Integer, S::FaureSample, T = Float64;
skipchecks = false)
skipchecks = false)
base = nextprime(dimension)
n_digits = ceil(Int, log(base, n))
power = n_digits - 1
Expand All @@ -60,7 +60,7 @@ end
end

@fastmath @views function _faure_samples(n_samples::I, n_digits::I, dimension::I,
::Type{F} = Float64) where {I <: Integer, F}
::Type{F} = Float64) where {I <: Integer, F}
base = nextprime(dimension)
inv_base = inv(base)

Expand Down
2 changes: 1 addition & 1 deletion src/Kronecker.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,7 +36,7 @@ function sample(n::Integer, d::Integer, S::KroneckerSample{Missing}, T = Float64
end

function sample(n::Integer, d::Integer, k::KroneckerSample{V},
T) where {V <: AbstractVector}
T) where {V <: AbstractVector}
@assert eltype(V)==T "Sample type must match generator type."
return randomize(sample(n, d, k), k.R)
end
Expand Down
2 changes: 1 addition & 1 deletion src/Lattices.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@ Base.@kwdef @concrete struct GridSample <: DeterministicSamplingAlgorithm
end

function sample(n::Integer, d::Integer, S::GridSample, T = Float64)
samples = rand.(range.(zeros(T, d), ones(T,d); length = n+1), Ref(n))
samples = rand.(range.(zeros(T, d), ones(T, d); length = n + 1), Ref(n))

Check warning on line 15 in src/Lattices.jl

View check run for this annotation

Codecov / codecov/patch

src/Lattices.jl#L15

Added line #L15 was not covered by tests
randomize(mapreduce(permutedims, vcat, samples), S.R)
end

Expand Down
40 changes: 20 additions & 20 deletions src/QuasiMonteCarlo.jl
Original file line number Diff line number Diff line change
Expand Up @@ -97,24 +97,24 @@ include("RandomizedQuasiMonteCarlo/scrambling_base_b.jl")
include("RandomizedQuasiMonteCarlo/iterators.jl")

export SamplingAlgorithm,
GridSample,
SobolSample,
LatinHypercubeSample,
RandomizedHaltonSample,
LatticeRuleSample,
RandomSample,
HaltonSample,
VanDerCorputSample,
GoldenSample,
KroneckerSample,
FaureSample,
randomize,
RandomizationMethod,
NoRand,
Shift,
ScrambleMethod,
OwenScramble,
MatousekScramble,
DigitalShift,
DesignMatrix
GridSample,
SobolSample,
LatinHypercubeSample,
RandomizedHaltonSample,
LatticeRuleSample,
RandomSample,
HaltonSample,
VanDerCorputSample,
GoldenSample,
KroneckerSample,
FaureSample,
randomize,
RandomizationMethod,
NoRand,
Shift,
ScrambleMethod,
OwenScramble,
MatousekScramble,
DigitalShift,
DesignMatrix
end # module
110 changes: 62 additions & 48 deletions src/RandomizedQuasiMonteCarlo/iterators.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@ end

"""
OwenDesignMat{T<:Real, I<:Integer, F<:Integer} <: AbstractDesignMatrix
Create an Owen scrambling iterator for doing multiple i.i.d. [`OwenScrambling`](@ref) randomization.
One can use the common[`DesignMatrix`](@ref) interface to create the iterator.
"""
Expand All @@ -22,6 +23,7 @@ end

"""
ScrambleDesignMat{T<:Real, I<:Integer} <: AbstractDesignMatrix
Create a scrambling iterator (Digital Shift or Matousek depending on the `R` field) for doing multiple i.i.d. [`DigitalShift`](@ref) or [`MatousekScrambling`](@ref) randomization.
One can use the common[`DesignMatrix`](@ref) interface to create the iterator.
"""
Expand All @@ -35,6 +37,7 @@ end

"""
ShiftDesignMat{T<:Real} <: AbstractDesignMatrix
Create a Shift iterator for doing multiple i.i.d [`Shift`](@ref) randomization.
One can use the common[`DesignMatrix`](@ref) interface to create the iterator.
"""
Expand All @@ -46,6 +49,7 @@ end

"""
DistributionDesignMat{T<:Real} <: AbstractDesignMatrix
Create an iterator for multiple distribution randomization. The distribution is chosen with the field `D`.
This is equivalent to using `rand!(D, X)` for some matrix `X`.
One can use the common[`DesignMatrix`](@ref) interface to create the iterator.
Expand All @@ -58,6 +62,7 @@ end

"""
RandomDesignMat{T<:Real} <: AbstractDesignMatrix
Create an iterator for multiple uniform randomization. This it similar to [`DistributionDesignMat`](@ref) with the field `D = Uniform()`
One can use the common[`DesignMatrix`](@ref) interface to create the iterator.
"""
Expand All @@ -76,44 +81,47 @@ Base.eltype(::Type{RandomDesignMat{T}}) where {T} = Matrix{T} # TODO one could c
```julia
DesignMatrix(n, d, sample_method::DeterministicSamplingAlgorithm, num_mats, T = Float64)
```
Create an iterator for doing multiple i.i.d. randomization over QMC sequences where
- `num_mats` is the length of the iterator
- `n` is the number of points to sample.
- `d` is the dimensionality of the point set in `[0, 1)ᵈ`,
- `sample_method` is the quasi-Monte Carlo sampling strategy used to create a deterministic point set `out`.
- `T` is the `eltype` of the point sets. For some QMC methods (Faure, Sobol) this can be `Rational`
It is now compatible with all scrambling methods and shifting. One can also use it with `Distributions.Sampleable` or `RandomSample`.
- `num_mats` is the length of the iterator
- `n` is the number of points to sample.
- `d` is the dimensionality of the point set in `[0, 1)ᵈ`,
- `sample_method` is the quasi-Monte Carlo sampling strategy used to create a deterministic point set `out`.
- `T` is the `eltype` of the point sets. For some QMC methods (Faure, Sobol) this can be `Rational`
It is now compatible with all scrambling methods and shifting. One can also use it with `Distributions.Sampleable` or `RandomSample`.
"""
function DesignMatrix(N, d, S::DeterministicSamplingAlgorithm, num_mats::Integer, T::DataType = Float64)
function DesignMatrix(
N, d, S::DeterministicSamplingAlgorithm, num_mats::Integer, T::DataType = Float64)
return DesignMatrix(N, d, S, S.R, num_mats, T)
end

function DesignMatrix(N,
d,
S::DeterministicSamplingAlgorithm,
R::OwenScramble,
num_mats::Integer,
T::DataType = Float64)
d,
S::DeterministicSamplingAlgorithm,
R::OwenScramble,
num_mats::Integer,
T::DataType = Float64)
X, random_bits, bits, indices = initialize(N, d, S, R, T)
return OwenDesignMat(X, random_bits, bits, indices, R, num_mats)
end

function DesignMatrix(N,
d,
S::DeterministicSamplingAlgorithm,
R::ScrambleMethod,
num_mats::Integer,
T = Float64)
d,
S::DeterministicSamplingAlgorithm,
R::ScrambleMethod,
num_mats::Integer,
T = Float64)
X, random_bits, bits = initialize(N, d, S, R, T)
return ScrambleDesignMat(X, random_bits, bits, R, num_mats)
end

function DesignMatrix(N,
d,
S::DeterministicSamplingAlgorithm,
R::Shift,
num_mats,
T::DataType = Float64)
d,
S::DeterministicSamplingAlgorithm,
R::Shift,
num_mats,
T::DataType = Float64)
X = initialize(N, d, S, R, T)
return ShiftDesignMat(X, R, num_mats)
end
Expand Down Expand Up @@ -144,10 +152,10 @@ function initialize(n, d, sampler, R::OwenScramble, T::DataType = Float64)
end

function scramble!(random_points::AbstractArray{T},
random_bits,
bits,
indices,
R::OwenScramble) where {T <: Real}
random_bits,
bits,
indices,
R::OwenScramble) where {T <: Real}
randomize_bits!(random_bits, bits, indices, R)
for i in CartesianIndices(random_points)
random_points[i] = bits2unif(T, @view(random_bits[:, i]), R.base)
Expand All @@ -170,9 +178,9 @@ function initialize(n, d, sampler, R::ScrambleMethod, T::DataType = Float64)
end

function scramble!(random_points::AbstractArray{T},
random_bits,
bits,
R::ScrambleMethod) where {T <: Real}
random_bits,
bits,
R::ScrambleMethod) where {T <: Real}
randomize_bits!(random_bits, bits, R)
for i in CartesianIndices(random_points)
random_points[i] = bits2unif(T, @view(random_bits[:, i]), R.base)
Expand Down Expand Up @@ -202,31 +210,35 @@ end
"""
```julia
generate_design_matrices(n, d, sample_method::DeterministicSamplingAlgorithm,
num_mats, T = Float64)
num_mats, T = Float64)
generate_design_matrices(n, d, sample_method::RandomSamplingAlgorithm,
num_mats, T = Float64)
num_mats, T = Float64)
generate_design_matrices(n, lb, ub, sample_method,
num_mats = 2)
num_mats = 2)
```
Create `num_mats` matrices each containing a QMC point set, where:
- `n` is the number of points to sample.
- `d` is the dimensionality of the point set in `[0, 1)ᵈ`,
- `sample_method` is the quasi-Monte Carlo sampling strategy used to create a deterministic point set `out`.
- `T` is the `eltype` of the point sets. For some QMC methods (Faure, Sobol) this can be `Rational`
If the bound `lb` and `ub` are specified instead of `d`, the samples will be transformed into the box `[lb, ub]`.
- `n` is the number of points to sample.
- `d` is the dimensionality of the point set in `[0, 1)ᵈ`,
- `sample_method` is the quasi-Monte Carlo sampling strategy used to create a deterministic point set `out`.
- `T` is the `eltype` of the point sets. For some QMC methods (Faure, Sobol) this can be `Rational`
If the bound `lb` and `ub` are specified instead of `d`, the samples will be transformed into the box `[lb, ub]`.
"""
function generate_design_matrices(n, d, sampler::DeterministicSamplingAlgorithm, num_mats::Integer,
T::DataType = Float64)
function generate_design_matrices(
n, d, sampler::DeterministicSamplingAlgorithm, num_mats::Integer,
T::DataType = Float64)
return generate_design_matrices(n, d, sampler, sampler.R, num_mats, T)
end

function generate_design_matrices(n, d, sampler::RandomSamplingAlgorithm, num_mats::Integer,
T::DataType = Float64)
function generate_design_matrices(

Check warning on line 234 in src/RandomizedQuasiMonteCarlo/iterators.jl

View check run for this annotation

Codecov / codecov/patch

src/RandomizedQuasiMonteCarlo/iterators.jl#L234

Added line #L234 was not covered by tests
n, d, sampler::RandomSamplingAlgorithm, num_mats::Integer,
T::DataType = Float64)
return [sample(n, d, sampler, T) for j in 1:num_mats]
end

function generate_design_matrices(n, lb, ub, sampler,
num_mats = 2)
num_mats = 2)
if n <= 0
throw(ZeroSamplesError())
end
Expand All @@ -246,22 +258,24 @@ end

"""
generate_design_matrices(n, d, sampler, R::NoRand, num_mats, T = Float64)
`R = NoRand()` produces `num_mats` matrices each containing a different deterministic point set in `[0, 1)ᵈ`.
Note that this is an ad hoc way to produce i.i.d sequence as it creates a deterministic point in dimension `d × num_mats` and split it in `num_mats` point set of dimension `d`.
This does not have any QMC garantuees.
"""
function generate_design_matrices(n, d, sampler, R::NoRand, num_mats::Integer, T::DataType = Float64)
function generate_design_matrices(

Check warning on line 266 in src/RandomizedQuasiMonteCarlo/iterators.jl

View check run for this annotation

Codecov / codecov/patch

src/RandomizedQuasiMonteCarlo/iterators.jl#L266

Added line #L266 was not covered by tests
n, d, sampler, R::NoRand, num_mats::Integer, T::DataType = Float64)
out = sample(n, num_mats * d, sampler, T)
@warn "The `generate_design_matrices(n, d, sampler, R = NoRand(), num_mats)` method does not produces true and independent QMC matrices, see [this doc warning](https://docs.sciml.ai/QuasiMonteCarlo/stable/design_matrix/) for more context.
Prefer using randomization methods such as `R = Shift()`, `R = MatousekScrambling()`, etc., see [documentation](https://docs.sciml.ai/QuasiMonteCarlo/stable/randomization/)"
return [out[(j * d + 1):((j + 1) * d), :] for j in 0:(num_mats - 1)]
end

function generate_design_matrices(n,
d,
sampler,
R::RandomizationMethod,
num_mats::Integer,
T::DataType = Float64)
d,
sampler,
R::RandomizationMethod,
num_mats::Integer,
T::DataType = Float64)
return collect(DesignMatrix(n, d, sampler, R, num_mats, T))
end
Loading

0 comments on commit 47f54c7

Please sign in to comment.