SSJ  3.3.1
Stochastic Simulation in Java
Public Member Functions | Static Public Member Functions | List of all members
RadicalInverse Class Reference

This class implements basic methods for working with radical inverses of integers in an arbitrary basis \(b\). More...

Public Member Functions

 RadicalInverse (int b, double x0)
 Initializes the base of this object to \(b\) and its first value of \(x\) to x0. More...
 
double nextRadicalInverse ()
 A fast method that incrementally computes the radical inverse \(x_{i+1}\) in base \(b\) from \(x_i\) = \(\psi_b(i)\), using addition with rigthward carry as described in [235] . More...
 

Static Public Member Functions

static int [] getPrimes (int n)
 Provides an elementary method for obtaining the first \(n\) prime numbers larger than 1. More...
 
static double radicalInverse (int b, long i)
 Computes the radical inverse of \(i\) in base \(b\). More...
 
static int radicalInverseInteger (int b, double x)
 Computes the radical inverse of \(x\) in base \(b\). More...
 
static long radicalInverseLong (int b, double x)
 Computes the radical inverse of \(x\) in base \(b\). More...
 
static double nextRadicalInverse (double invb, double x)
 A fast method that incrementally computes the radical inverse \(x_{i+1}\) in base \(b\) from \(x_i\) = x = \(\psi_b(i)\), using addition with rigthward carry. More...
 
static void reverseDigits (int k, int bdigits[], int idigits[])
 A fast method that incrementally computes the radical inverse \(x_{i+1}\) in base \(b\) from \(x_i =\) x, using the Struckmeier’s algorithm described in [227] . More...
 
static int integerRadicalInverse (int b, int i)
 Computes the integer radical inverse of \(i\) in base \(b\), equal to \(b^k \psi_b(i)\) if \(i\) has \(k\) \(b\)-ary digits. More...
 
static int nextRadicalInverseDigits (int b, int k, int idigits[])
 Given the \(k\) digits of the integer radical inverse of \(i\) in bdigits, in base \(b\), this method replaces them by the digits of the integer radical inverse of \(i+1\) and returns their number. More...
 
static void getFaureLemieuxPermutation (int coordinate, int[] pi)
 Computes the permutations as proposed in [58]  \(\sigma_b\) of the set \(\{0, …, b - 1\}\) and puts it in array pi. More...
 
static void getFaurePermutation (int b, int[] pi)
 Computes the Faure permutation [60]  \(\sigma_b\) of the set \(\{0, …, b - 1\}\) and puts it in array pi. More...
 
static double permutedRadicalInverse (int b, int[] pi, long i)
 Computes the radical inverse of \(i\) in base \(b\), where the digits are permuted using the permutation \(\pi\). More...
 

Detailed Description

This class implements basic methods for working with radical inverses of integers in an arbitrary basis \(b\).

These methods are used in classes that implement point sets and sequences based on the van der Corput sequence (the Hammersley nets and the Halton sequence, for example).

We recall that for a \(k\)-digit integer \(i\) whose digital \(b\)-ary expansion is

\[ i = a_0 + a_1 b + …+ a_{k-1} b^{k-1}, \]

the radical inverse in base \(b\) is

\[ \psi_b(i) = a_0 b^{-1} + a_1 b^{-2} + \cdots+ a_{k-1} b^{-k}. \]

The van der Corput sequence in base \(b\) is the sequence \(\psi_b(0), \psi_b(1), \psi_b(2), …\)

Note that \(\psi_b(i)\) cannot always be represented exactly as a floating-point number on the computer (e.g., if \(b\) is not a power of two). For an exact representation, one can use the integer

\[ b^k \psi_b(i) = a_{k-1} + \cdots+ a_1 b^{k-2} + a_0 b^{k-1}, \]

which we called the integer radical inverse representation. This representation is simply a mirror image of the digits of the usual \(b\)-ary representation of \(i\).

It is common practice to permute locally the values of the van der Corput sequence. One way of doing this is to apply a permutation to the digits of \(i\) before computing \(\psi_b(i)\). That is, for a permutation \(\pi\) of the digits \(\{0,…,b-1\}\),

\[ \psi_b(i) = \sum_{r=0}^{k-1} a_r b^{-r-1} \]

is replaced by

\[ \sum_{r=0}^{k-1} \pi(a_r) b^{-r-1}. \]

Applying such a permutation only changes the order in which the values of \(\psi_b(i)\) are enumerated. For every integer \(k\), the first \(b^k\) values that are enumerated remain the same (they are the values of \(\psi_b(i)\) for \(i=0,…,b^k-1\)), but they are enumerated in a different order. Often, different permutations \(\pi\) will be applied for different coordinates of a point set.

The permutation \(\pi\) can be deterministic or random. One (deterministic) possibility implemented here is the Faure permutation \(\sigma_b\) of \(\{0,…,b-1\}\) defined as follows [60] . For \(b=2\), take \(\sigma= I\), the identical permutation. For even \(b=2c > 2\), take

\begin{align} \sigma[i] & = 2\tau[i]\phantom{{} + 1} \qquad i = 0, 1, …, c-1 \\ \sigma[i+c] & = 2\tau[i] + 1 \qquad i = 0, 1, …, c-1 \end{align}

where \(\tau[i]\) is the Faure permutation for base \(c\). For odd \(b=2c+1\), take

\begin{align} \sigma[c] & = c \\ \sigma[i] & = \tau[i],\phantom{{} + 1} \qquad\mbox{ if } 0 \le\tau[i] < c \\ \sigma[i] & = \tau[i] + 1, \qquad\mbox{ if } c \le\tau[i] < 2c \end{align}

for \(0 \le i < c\), and take

\begin{align} \sigma[i] & = \tau[i-1],\phantom{{} + 1} \qquad\mbox{ if } 0 \le\tau[i-1] < c \\ \sigma[i] & = \tau[i-1]+1, \qquad\mbox{ if } c \le\tau[i-1] < 2c \end{align}

for \(c < i \le2c\), and where \(\tau[i]\) is the Faure permutation for base \(c\). The Faure permutations give very small discrepancies (amongst the best known today) for small bases.

Constructor & Destructor Documentation

◆ RadicalInverse()

RadicalInverse ( int  b,
double  x0 
)

Initializes the base of this object to \(b\) and its first value of \(x\) to x0.

Parameters
bBase
x0Initial value of x

Member Function Documentation

◆ getFaureLemieuxPermutation()

static void getFaureLemieuxPermutation ( int  coordinate,
int []  pi 
)
static

Computes the permutations as proposed in [58]  \(\sigma_b\) of the set \(\{0, …, b - 1\}\) and puts it in array pi.

Parameters
coordinatethe coordinate
pian array of size at least b, to be filled with the permutation

◆ getFaurePermutation()

static void getFaurePermutation ( int  b,
int []  pi 
)
static

Computes the Faure permutation [60]  \(\sigma_b\) of the set \(\{0, …, b - 1\}\) and puts it in array pi.

See the description in the introduction above.

Parameters
bthe base
pian array of size at least b, to be filled with the permutation

◆ getPrimes()

static int [] getPrimes ( int  n)
static

Provides an elementary method for obtaining the first \(n\) prime numbers larger than 1.

Creates and returns an array that contains these numbers. This is useful for determining the prime bases for the different coordinates of the Halton sequence and Hammersley nets.

Parameters
nnumber of prime numbers to return
Returns
an array with the first n prime numbers

◆ integerRadicalInverse()

static int integerRadicalInverse ( int  b,
int  i 
)
static

Computes the integer radical inverse of \(i\) in base \(b\), equal to \(b^k \psi_b(i)\) if \(i\) has \(k\) \(b\)-ary digits.

Parameters
bbase used for the operation
ithe value for which the integer radical inverse will be computed
Returns
the integer radical inverse of i in base b

◆ nextRadicalInverse() [1/2]

static double nextRadicalInverse ( double  invb,
double  x 
)
static

A fast method that incrementally computes the radical inverse \(x_{i+1}\) in base \(b\) from \(x_i\) = x = \(\psi_b(i)\), using addition with rigthward carry.

The parameter invb is equal to \(1/b\). Using long incremental streams (i.e., calling this method several times in a row) cause increasing inaccuracy in \(x\). Thus the user should recompute the radical inverse directly by calling radicalInverse every once in a while (i.e. in every few thousand calls).

Parameters
invb\(1/b\) where \(b\) is the base
xthe inverse \(x_i\)
Returns
the radical inverse \(x_{i+1}\)

◆ nextRadicalInverse() [2/2]

double nextRadicalInverse ( )

A fast method that incrementally computes the radical inverse \(x_{i+1}\) in base \(b\) from \(x_i\) = \(\psi_b(i)\), using addition with rigthward carry as described in [235] .

Since using long incremental streams (i.e., calling this method several times in a row) cause increasing inaccuracy in \(x\), the method recomputes the radical inverse directly from \(i\) by calling radicalInverse once in every 1000 calls.

Returns
the radical inverse \(x_{i+1}\)

◆ nextRadicalInverseDigits()

static int nextRadicalInverseDigits ( int  b,
int  k,
int  idigits[] 
)
static

Given the \(k\) digits of the integer radical inverse of \(i\) in bdigits, in base \(b\), this method replaces them by the digits of the integer radical inverse of \(i+1\) and returns their number.

The array must be large enough to hold this new number of digits.

Parameters
bbase
kinitial number of digits in arrays
idigitsdigits of integer radical inverse
Returns
new number of digits in arrays

◆ permutedRadicalInverse()

static double permutedRadicalInverse ( int  b,
int []  pi,
long  i 
)
static

Computes the radical inverse of \(i\) in base \(b\), where the digits are permuted using the permutation \(\pi\).

If \(i=\sum_{r=0}^{k-1} a_r b^r\), the method will compute and return

\[ x = \sum_{r=0}^{k-1} \pi[a_r] b^{-r-1}. \]

Parameters
bbase \(b\) used for the operation
pian array of length at least b containing the permutation used during the computation
ithe value for which the radical inverse will be computed
Returns
the radical inverse of i in base b

◆ radicalInverse()

static double radicalInverse ( int  b,
long  i 
)
static

Computes the radical inverse of \(i\) in base \(b\).

If \(i=\sum_{r=0}^{k-1} a_r b^r\), the method computes and returns

\[ x = \sum_{r=0}^{k-1} a_r b^{-r-1}. \]

Parameters
bbase used for the operation
ithe value for which the radical inverse will be computed
Returns
the radical inverse of i in base b

◆ radicalInverseInteger()

static int radicalInverseInteger ( int  b,
double  x 
)
static

Computes the radical inverse of \(x\) in base \(b\).

If \(x\) has more decimals in base \(b\) than \(\log_b\)(Long.MAX_VALUE), it is truncated to its minimum precision in base \(b\). If \(x=\sum_{r=0}^{k-1} a_r b^{-r-1}\), the method computes and returns

\[ i = \sum_{r=0}^{k-1} a_r b^r. \]

Parameters
bbase used for the operation
xthe value for which the radical inverse will be computed
Returns
the radical inverse of x in base b

◆ radicalInverseLong()

static long radicalInverseLong ( int  b,
double  x 
)
static

Computes the radical inverse of \(x\) in base \(b\).

If \(x\) has more decimals in base \(b\) than \(\log_b\)(Long.MAX_VALUE), it is truncated to its minimum precision in base \(b\). If \(x=\sum_{r=0}^{k-1} a_r b^{-r-1}\), the method computes and returns

\[ i = \sum_{r=0}^{k-1} a_r b^r. \]

Parameters
bbase used for the operation
xthe value for which the radical inverse will be computed
Returns
the radical inverse of x in base b

◆ reverseDigits()

static void reverseDigits ( int  k,
int  bdigits[],
int  idigits[] 
)
static

A fast method that incrementally computes the radical inverse \(x_{i+1}\) in base \(b\) from \(x_i =\) x, using the Struckmeier’s algorithm described in [227] .

It uses a small precomputed table of values \(L_k = 1 - b^{-k}\). The method returns the next radical inverse \(x_{i+1} = x_i + (b + 1 - b^k) / b^k\), where \(L_{k-1} \le x < L_k\).

Remarks
Richard: This method can work only if it is reprogrammed with integers. With floating-point numbers, unavoidable accumulating rounding errors will sooner or later lead to choosing the wrong interval, after which, all subsequent x’s will be completely wrong.
Parameters
bbase
xthe inverse \(x_i\)
Returns
the radical inverse \(x_{i+1}\) Given the \(k\) \(b\)-ary digits of \(i\) in bdigits, returns the \(k\) digits of the integer radical inverse of \(i\) in idigits. This simply reverses the order of the digits.
Parameters
knumber of digits in arrays
bdigitsdigits in original order
idigitsdigits in reverse order

The documentation for this class was generated from the following file: