FAURE_NEXT

Computes a shuffled Faure sequence.

Required Arguments

STATE — An IMSL_FAURE pointer containing the structure created by the call to FAURE_INIT. The structure contains information about the sequence. The structure should be freed using FAURE_FREE after it is no longer needed. (Input/Output)

NEXT_PT — Vector of length NDIM containing the next point in the shuffled Faure sequence, where NDIM is the dimension of the hyper-rectangle specified in FAURE_INIT. (Output)

Optional Arguments

IMSL_RETURN_SKIP — Returns the current point in the sequence. The sequence can be restarted by calling FAURE_INIT using this value for NSKIP, and using the same value for NDIM. (Input)

FORTRAN 90 Interface

Generic: CALL FAURE_NEXT (STATE, NEXT_PT [, …])

Specific: The specific interface names are S_FAURE_NEXT and D_FAURE_NEXT.

Description

The functions FAURE_INIT and FAURE_NEXT are used to generate shuffled Faure sequence of low discrepancy n-dimensional points. Low discrepancy series fill an n-dimensional cube more uniformly than pseudo-random sequences, and are used in multivariate quadrature, simulation, and global optimization. Because of this uniformity, use of low discrepancy series is generally more efficient than pseudo-random series for multivariate Monte Carlo methods. See the IMSL routine QMC (in Chapter 4, “Integration and Differentiation”) for a discussion of quasi-Monte Carlo quadrature based on low discrepancy series.

Discrepancy measures the deviation from uniformity of a point set.

The discrepancy of the point set , is defined

 

 

where the supremum is over all subsets of [0, 1]d of the form

 

λ is the Lebesque measure, and is the number of the xj contained in E.

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

 

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 NBASE defaults to the smallest prime greater than or equal to the dimension.

The generalized Faure sequence x1, x2, …, is computed as follows:

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

 

where ai(n) are integers, .

The j-th coordinate of xn is

 

The generator matrix for the series, , is defined to be

 

and is an element of the Pascal matrix,

 

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 FAURE_INIT is used to create a structure that holds the state of the sequence. Each call to FAURE_NEXT returns the next point in the sequence and updates the IMSL_FAURE structure. The final call to FAURE_FREE frees data items, stored in the structure, that were allocated by FAURE_INIT.

 

use faure_int

implicit none

type (s_imsl_faure), pointer :: state

real(kind(1e0)) :: x(3)

integer,parameter :: ndim=3

integer :: k

! CREATE THE STRUCTURE THAT HOLDS

! THE STATE OF THE SEQUENCE.

call faure_init(ndim, state)

! GET THE NEXT POINT IN THE SEQUENCE

do k=1,5

call faure_next(state, x)

write(*,'(3F15.3)') x(1), x(2) , x(3)

enddo

! FREE DATA ITEMS STORED IN

! state STRUCTURE

call faure_free(state)

end

Output

 

0.334 0.493 0.064

0.667 0.826 0.397

0.778 0.270 0.175

0.111 0.604 0.509

0.445 0.937 0.842