## Overview

This public domain single header file library (Sobol.h) implements basic functions for sampling Sobol sequences. It implements two types of generators:

- Progressive, which is the standard sampling method. Specifically a generator can be queried for any number of samples. The library supports 1-3 dimensional generators.
- Stratified generators. These generate a sequence one dimension lower, require specifying the number of samples which will be generated and the remaining dimension is equispaced over the interval.

## Documentation

The library consists of a data table and mostly inline functions. Any functions which are not defined as inline are noted below. As such one file must define a macro (`#define SOBOL_IMPLEMENTATION`

) prior to including the header.

#### Data structures

Progressive generators contains an implied position element $ i $ followed by a state element per dimension.

Stratified generators start with the same data layout as progressive generators of one dimension lower followed by the element $r$ which is $\frac{1}{number\ of\ samples}$ .

#### Initalizing

Initializing a generator sets its position to zero and caller provides the initial value for the state of each dimension. This is logically equivalent to XORing this value to the given dimension result at each query. (SEE: implementation details)

Standard generators

Progressive generators additionally need the number of samples which will be generated.

#### Position in sequence manipulation

Returns the current position in the sequence of the generator. This accepts any type of generator.

Moves the current position by the relative value *offset*. Only standard generators are directly supported. These functions are not defined as inline.

#### Reading sequence members

Each dimension returns results on the unit interval: $ \left[0,\thinspace1\right) $ by default. (SEE: implementation details)

The 1D version returns the next member of the sequence as either a single or double and increments the position.

nD versions place the $n$ dimensional elements sequentially at $d$, in either as singles or doubles, and increments the position.

#### Exposed internals

The state updating functions. (SEE: implementation details)

#### Optional features

The following funtions are only available with `#define SOBOL_EXTRA`

. These additionally require a single precision square root function be defined as `sqrtf`

. None of these functions are defined as inline.

Uniform sampling of (open) unit disk ($\mathbb{D}^1$),half-disk ($x\geq 0$) and quater-disk ($x,y\geq 0$)are performed using a rejection method. This has no impact on the well-space properties. The average number of iterations is $\frac{4}{\pi} \approx 1.27324$. They both take a 2D generator $s$ and place the point in $p$. The returned value is $p \cdot p$.

Uniform sampling of (open) unit 3D sphere ($\mathbb{S}^2$) and half-sphere ($z\geq 0$) are generated using Marsaglia’s transform^{1} of the unit disk.

## Implementation details and configuration options

#### State update

The basic construction is the classic bit-reversed gray codes^{2}. The state in each explicit dimension is a 32-bit integer which at each update has some number of bits flipped based on the trailing one count of the position. Instead of storing and updating the explicit position, the library instead stores its bit-not and adding one becomes subtracting one ($ \sim x = -x-1 $). This allows the state update to directly use a trailing zero count (without an extra bit-not op). The input to the trailing count is only zero when an update is wrapping the position from -1 to 0, so intel’s `BSF` (bit scan forward) can be used and is in fact the libraries default.

#### Configuration

if a macro `#define SOBOL_BIAS`

is defined (default is undefined) then single precision queries do not perform a bit-shift to convert the integer to floating point. This causes single results to be on the interval $ \left[0,\thinspace1\right] $, double results are uneffected.

The library needs to be able to perform a trailing zero count. Currently this is only set-up for x86-a-likes.

`#define SOBOL_NTZ(X) number_of_trailing_zeroes_function(X)`