Skip to content

githubgrasp/genran

Repository files navigation

genran

A non-Gaussian random field generator. Produces 1D, 2D, or 3D scalar fields with a prescribed marginal CDF (Gaussian, Weibull, or grafted Weibull–Gaussian) and a prescribed Gaussian-shaped autocorrelation. The underlying method is the spectral representation of Shinozuka & Deodatis; optional iterative refinement (Shields–Deodatis–Bocchini 2011) recovers the prescribed autocorrelation for non-Gaussian targets.

Originally written in 2006; rewritten and modernised since.

Distributed under the MIT licence.


Build

make            # produces ./genran.exe
make test       # runs the regression harness (must pass)
make validate   # empirical demonstration of the SDB-2011 iteration
make clean

Plain g++; no external dependencies. Tested on macOS (clang/Apple toolchain) and Linux (gcc).


Quick start

./genran.exe test1D.in field1D.dat

reads the simulation parameters from test1D.in, writes the generated field to field1D.dat, and writes a mean sdev/mean summary to stat.dat in the current directory.

The repository ships with test{1,2,3}D.in covering 1D Weibull and 2D/3D Gaussian configurations, plus iteration and grafted-CDF variants (test*_iterate.in, test1D_grafted*.in).


Input file reference

Input files are .ini-style: one key = value per line, : accepted as a synonym for =. Lines that don't begin with a recognised key are ignored, so you can scatter comments freely. Keys are case-insensitive on the prefix.

Required, all configurations

Key Meaning
numberOfDimensions 1, 2, or 3
realNumber1 (and realNumber2, realNumber3 for nDim ≥ 2) Number of grid points along each axis. Each must be a positive power of two (radix-2 FFT requirement).
xOrigin (and yOrigin, zOrigin) Origin of the output coordinate system in user units.
autoLength1 (and autoLength2, autoLength3) Autocorrelation length la along each axis, in user spatial units. The autocorrelation kernel is R(τ) = σ² exp(-(τ/b)²) with b = 2·la/√π.
converter Spatial-grid step inverse: samples per user-length-unit. The total domain along axis i is (realNumber_i - 1)/converter.
typeOfCDF Marginal distribution: 1 Gaussian, 2 Weibull, 3 grafted Weibull–Gaussian.

Distribution parameters

For typeOfCDF = 1 (Gaussian):

Key Meaning
mean Target mean.
cov Target coefficient of variation, σ/μ.

For typeOfCDF = 2 (Weibull):

Key Meaning
scaling Weibull scale parameter.
modulus Weibull shape parameter (the modulus).

The Weibull mean and variance are computed analytically from these.

For typeOfCDF = 3 (grafted Weibull–Gaussian, see Bažant & Pang 2006):

Key Meaning
scaling Weibull scale of the lower tail.
modulus Weibull modulus of the lower tail.
graft Grafting probability (probability mass below the splice).
cov Coefficient of variation entering the rescaling formula.

The Gaussian-core parameters and the rescaling factor are computed at runtime from the four user inputs.

Optional

Key Default Meaning
ranint time-based Negative integer seed. Provide a fixed value for reproducible runs. Non-negative or absent values fall back to -time(NULL).
padding 1 Power of two. Simulates on a padding × N grid in each dimension and crops the user's requested corner. Suppresses the FFT's implicit-periodicity boundary contamination at the cost of padding^nDim memory and FFT work. padding = 2 is a reasonable default if boundary statistics matter.
iterate 0 Set to 1 to run the SDB-2011 iteration that recovers the prescribed autocorrelation for non-Gaussian targets. Without it, non-Gaussian fields have the correct marginal but a distorted autocorrelation (the memoryless mapping squashes correlations). Has no effect for typeOfCDF = 1.
iterMaxIter 100 Hard ceiling on the SDB-2011 loop.
iterTolerance 0.01 The loop terminates once the per-iteration improvement in relative spectral error drops below this many percentage points.
iterBeta 1.4 Update exponent (Eq. 18 of SDB-2011). The paper recommends 1.3–1.5; smaller values are slower and more accurate, larger values risk instability.
iterHermiteOrder 30 Number of Hermite terms in the Mehler expansion of the translation kernel. Increase if the Σ c_k²/k! value printed at start-up does not match the target σ² closely.

Output format

The field is written to the file given as argv[2]. Coordinates are in the user's spatial units (i.e. converted back from the simulator's internal scaled coordinates).

1D:

x_0   value
x_1   value
...

2D: the first line carries N1 N2 dx1 dx2; subsequent lines are

x_k   y_l   value

with l varying fastest, blank lines between successive k.

3D: first line N1 N2 N3 dx1 dx2 dx3; subsequent lines

x_k   y_l   z_m   value

with m fastest, then l, then k.

A second file stat.dat is always written in the current directory and contains a single line mean sdev/mean.


Algorithms

Spectral representation

Sample functions are generated by

f(x) = √2 · Σ_n A(κ_n) · cos(κ_n · x + φ_n)

with random phases φ_n uniformly distributed on [0, 2π] and amplitudes A(κ_n) = √(2 S(κ_n) Δκ). The wavenumber sums are evaluated via FFT. Multi-dimensional fields use 2^(D−1) sign combinations of the wavenumber vector, following Shinozuka & Deodatis (1996).

Translation to a non-Gaussian marginal

Each Gaussian sample g(x) is mapped to the target distribution by

X(x) = F_NG⁻¹(Φ(g(x)/σ_G)),

where F_NG⁻¹ is the inverse of the prescribed marginal CDF and σ_G is the standard deviation of the underlying Gaussian field. This memoryless transform preserves the marginal exactly, but distorts the autocorrelation.

SDB-2011 iteration (iterate = 1)

To recover the prescribed autocorrelation, the underlying Gaussian PSDF S_G(κ) is iterated using

S_G⁽ⁱ⁺¹⁾(κ) = (S_T(κ) / S_NG⁽ⁱ⁾(κ))^β · S_G⁽ⁱ⁾(κ),

where S_NG⁽ⁱ⁾ is the spectrum of the non-Gaussian field implied by S_G⁽ⁱ⁾ through translation theory. S_NG⁽ⁱ⁾ is computed via Mehler's formula

R_NG(τ) = Σ_{k≥1} c_k² · ξ(τ)^k / k!,

avoiding sample-function generation in the inner loop. The Hermite coefficients c_k are computed once per run by quadrature against the standard normal density.

Spectral discretisation

The cutoff wavenumber is hard-coded to κ_max = 1.8π (in scaled coordinates), with Δκ = κ_max / (N/2). The DC and Nyquist bins are zeroed by convention.

References

  • Shinozuka & Jan, J. Sound and Vibration 25, 111–128 (1972).
  • Shinozuka & Deodatis, Appl. Mech. Rev. 44, 191–204 (1991).
  • Shinozuka & Deodatis, Appl. Mech. Rev. 49, 29–53 (1996).
  • Deodatis & Micaletti, J. Eng. Mech. 127, 1284–1295 (2001).
  • Shields, Deodatis & Bocchini, Probab. Eng. Mech. 26, 511–519 (2011).
  • Kim & Shields, Comput. Struct. 161, 31–42 (2015).
  • Bažant & Pang, J. Mech. Phys. Solids 55, 91–131 (2007).

Tests

make test runs the regression harness in tests/run.sh. Each scenario (see the cases enumerated at the bottom of the script) is generated with a fixed ranint seed; the program's stat.dat and a SHA256 of the field file are diffed against committed golden values in tests/golden/. All scenarios must pass.

The harness assumes std::mt19937_64 is byte-identical for a given seed across machines (it is by C++ standard). FFT round-off can in principle differ across toolchains; if a fresh machine reports hash mismatches without algorithm changes, regenerate the goldens once and recommit.

make validate runs an empirical sanity check: it generates a long 1D Weibull field with iterate=0 and iterate=1 (modulus = 1.5), computes the empirical autocorrelation of each, and reports the maximum deviation from the prescribed Gaussian target. The iterated version is expected to beat the un-iterated one by at least 2×; in practice the gap is 3–5× for moderately skewed distributions.

tools/autocorr is the small standalone helper used by the validation script. It can also be run directly to inspect the autocorrelation of any 1D field file.


Repository layout

genran.C                    command-line driver, input parsing, output writing
randomfieldgenerator.C      orchestration, SDB-2011 iteration, mapping
gaussianfieldgenerator.C    spectral-representation Gaussian field synthesis
fourier.C                   radix-2 Cooley-Tukey FFT
help.C                      moments, RNG, gamma, normal CDF / inverse CDF,
                            grafted-CDF rescaling
format.C                    .ini-style input parsing
random.h                    public function declarations
test{1,2,3}D.in             stock 1D Weibull, 2D Gaussian, 3D Gaussian inputs
test*_iterate.in            iterate-enabled variants
test1D_grafted*.in          grafted Weibull-Gaussian inputs
tests/run.sh                regression harness
tests/golden/               committed golden statistics and field hashes
tests/validate_iteration.sh empirical demonstration of the iteration's effect
tools/autocorr.cpp          autocorrelation tool used by the validation
Literature/                 reference PDFs (untracked)

Licence

MIT.

About

Random field generator

Resources

Stars

Watchers

Forks

Releases

No releases published

Packages

 
 
 

Contributors