faureNextPoint¶

Computes a shuffled Faure sequence.

Synopsis¶

faureSequenceInit (ndim)

faureNextPoint (state)

faureSequenceFree (state)

Required Arguments for faureSequenceInit¶

int ndim (Input): The dimension of the hyper-rectangle.

Return Value for faureSequenceInit¶

Returns a structure that contains information about the sequence. The structure should be freed using faureSequenceFree after it is no longer needed.

Required Arguments for faureNextPoint¶

structure state (Input/Output): Structure created by a call to faureSequenceInit.

Return Value for faureNextPoint¶

Returns the next point in the shuffled Faure sequence. To release this space, use faureSequenceFree.

Required Arguments for faureSequenceFree¶

structure state (Input/Output): Structure created by a call to faureSequenceInit.

Optional Arguments¶

base, int (Input)

The base of the Faure sequence.

Default: The smallest prime greater than or equal to ndim.

skip, int (Input)

The number of points to be skipped at the beginning of the Faure sequence.

Default: \(\lfloor\mathtt{base}^{m/2-1} \rfloor\), where \(m=\lfloor\log \mathrm{B}/\log\mathtt{base} \rfloor\) and B is the largest representable integer.

returnSkip (Output)

The current point in the sequence. The sequence can be restarted by initializing a new sequence using this value for skip, and using the same dimension for ndim.

Description¶

Discrepancy measures the deviation from uniformity of a point set.

The discrepancy of the point set \(x_1,\ldots,x_n\in[0,1]^d,d\geq 1\), is

\[D_n^{(d)} = \sup_E \left| \frac{A(E; n)}{n} - \lambda(E)\right|,\]

where the supremum is over all subsets of \(\left[0,1\right]^d\) of the form

\[E = \left[0,t_1\right) \times \ldots \times \left[0,t_d\right), 0 \leq t_j \leq 1, 1 \leq j \leq d\]

λ is the Lebesque measure, and \(A (E; n)\) is the number of the \(x_j\) contained in E.

The sequence \(x_1\), \(x_2\), …of points \([0,1]^d\) is a low-discrepancy sequence if there exists a constant c(d), depending only on d, such that

\[D_n^{(d)} \leq c(d) \frac{(\log n)^d}{n}\]

for all n>1.

Generalized Faure sequences can be defined for any prime base b≥d. The lowest bound for the discrepancy is obtained for the smallest prime b≥d, so the optional argument base defaults to the smallest prime greater than or equal to the dimension.

The generalized Faure sequence \(x_1\), \(x_2\), …, is computed as follows:

Write the positive integer n in its b-ary expansion,

\[n = \sum_{i=0}^{\infty} a_i(n)b^i\]

where \(a_i(n)\) are integers, \(0\leq a_i (n)<b\).

The j-th coordinate of \(x_n\) is

\[x_n^{(j)} = \sum_{k=0}^{\infty}\sum_{d=0}^{\infty} c_{kd}^{(j)}a_d(n)b^{-k-1}, 1 \leq j \leq d\]

The generator matrix for the series, \(c_{kd}^{(j)}\), is defined to be

\[c_{kd}^{(j)} = j^{d-k}c_{kd}\]

and \(c_{kd}\) is an element of the Pascal matrix,

\[\begin{split}c_{kd} = \begin{cases} \tfrac{d!}{c!(d-c)!} & k \leq d \\ \hfil 0 & k > d \end{cases}\end{split}\]

It is faster to compute a shuffled Faure sequence than to compute the Faure sequence itself. It can be shown that this shuffling preserves the low-discrepancy property.

The shuffling used is the b-ary Gray code. The function G(n) maps the positive integer n into the integer given by its b-ary expansion.

The sequence computed by this function is x(G(n)), where x is the generalized Faure sequence.

Example¶

In this example, five points in the Faure sequence are computed. The points are in the three-dimensional unit cube.

Note that faureSequenceInit (see Synopsis) is used to create a structure that holds the state of the sequence. Each call to faureNextPoint returns the next point in the sequence and updates the structure. The final call to faureSequenceFree (see Synopsis) frees data items, stored in the structure, that were allocated by faureSequenceInit.

from __future__ import print_function
from numpy import *
from pyimsl.math.faureSequenceInit import faureSequenceInit
from pyimsl.math.faureNextPoint import faureNextPoint
from pyimsl.math.faureSequenceFree import faureSequenceFree

ndim = 3
state = faureSequenceInit(ndim)

for k in range(0, 5):
    x = faureNextPoint(state)
    print("%10.3f %10.3f  %10.3f" % (x[0], x[1], x[2]))

faureSequenceFree(state)

Output¶

334      0.493       0.064
667      0.826       0.397
778      0.270       0.175
111      0.604       0.509
445      0.937       0.842