pf_random_permutation Interface

public interface pf_random_permutation

Fills perm with the whole permutation of 1 .. size(perm) under seed.

The bulk form of pf_random_perm_at, and an identity rather than a second algorithm: perm(k) equals pf_random_perm_at(seed, size(perm), k) for every k, so the two may be mixed freely and a prefix of one is a prefix of the other. That is the same relationship pf_random_at and pf_random_fill_draws already have, and it is what keeps the two forms from drifting apart under a later optimisation.

perm is a rank-1 integer(int32) or integer(int64) array, intent(out); seed is integer(int64). A zero-sized perm is a defined no-op.

threads is optional and changes only how fast the array is filled, never what is in it -- and that is provable here rather than merely intended, which is unusual for a threading argument. perm(k) is a pure function of (seed, size(perm), k) and reads nothing another element writes, so the 1-thread and N-thread results are bit-identical and a program may switch thread counts while debugging without its permutation changing under it. Absent means automatic: as many threads as OpenMP offers, capped by parquet_set_random_threads, and 1 inside an OpenMP parallel region, since a nested region is the caller's business. An explicit threads= is honoured there too.

A work floor applies to both, and it is not a tuning detail. Below parquet_set_random_parallel_min_elements (default 1000) elements per thread the call runs serially, because threading a small permutation is measurably slower than not threading it -- so threads=8 on a 100-element array is honoured by being declined.

It is about 2.4x cheaper per element than calling the elemental form size(perm) times (machine B, gfortran 14.2.1, --profile release: 9.9 ns against 24.1 ns at m = 10**8), for two structural reasons rather than any tuning. Enumerating the domain in order hands the loop the (l, r) split that the scalar form has to recover with an integer division; and, the larger of the two, the width rule and the key schedule are derived once per call instead of once per element, which a pure elemental function has no way to avoid. The values are unchanged by either -- only the route to them is.

Module Procedures