lib: Make random number generation thread safe

[marisn]

As random file generation keeps internal state in a static variable, it is prone to failures as multiple threads would attempt to update sate simultaneously. This PR adds a mutex lock around internal state update calls thus preventing clashes when called from multiple threads.

Contains code to be merged as PR #6479
Also contains necessary changes in Make/CMake to compile gis lib with pthread support.

Should fix issues observed in r.mapcalc after thread support was merged #5742

[wenzeslaus]

The code changes look straightforward. Do you have an idea about performance? I agree we need the change in any case, but at least knowing whether we need a better solution in the future would be good.

[wenzeslaus]

Good news. This makes the failing rand function tests in #6476 pass at least on Linux.

[marisn]

The code changes look straightforward. Do you have an idea about performance? I agree we need the change in any case, but at least knowing whether we need a better solution in the future would be good.

I haven't tested but locked part is small and should be fast (just update, no IO), thus I would expect impact on speed to be negligible.

[wenzeslaus]

I think the locking needs to expand to reading of the shared state, it is both write and read which need to be isolated.

More generally, I'm afraid we didn't finalize the r.mapcalc's rand behavior with multiple threads (and seeds). (Or did we?) If r.mapcalc is using these functions directly as before parallelization, the random numbers placed in the raster will depend on the order of threads, i.e., what rows will come first, but I did not review the relevant code now. (An ultimate fix there would be different "context" objects for the random numbers, per row (for maximum reproducibility) or per thread (for at least seed+nprocs reproducibility).

[wenzeslaus]

I haven't tested but locked part is small and should be fast (just update, no IO), thus I would expect impact on speed to be negligible.

No IO, but it is called for every cell in r.mapcalc, no? Still better slow with overall faster code than a segfault, and rand() is already much slower than let's say a constant (I did not test recently). If it is slow, we could use a solution with an explicit state (context) which may allow us to avoid locking, while having a separate locking version of the API.

[marisn]

This is the worst case scenario – code calls random number generation and discards produced numbers (thread creation + locking overhead):

Total random calls per run: 100,000,000
Runs per thread count: 5

Threads	Avg. Time (s)	Overhead vs. 1 Thread
1	0.8202	+0.0000s (+0.00%)
2	3.3390	+2.5187s (+307.07%)
4	3.2877	+2.4674s (+300.81%)
6	3.5505	+2.7303s (+332.86%)
8	4.3752	+3.5549s (+433.40%)
10	4.9465	+4.1263s (+503.05%)
12	5.7690	+4.9488s (+603.32%)
16	6.2315	+5.4113s (+659.71%)

---

[marisn]

No IO, but it is called for every cell in r.mapcalc, no? Still better slow with overall faster code than a segfault, and rand() is already much slower than let's say a constant (I did not test recently). If it is slow, we could use a solution with an explicit state (context) which may allow us to avoid locking, while having a separate locking version of the API.

RN generation now is always serial – outcome order is guaranteed. Consumption is another topic. I lean towards just stating that for fully reproducible outcome one should set parameters for single thread execution.
Generally it is out of scope of this PR and should be discussed separately with more people involved.

[wenzeslaus]

A random sequence can contain duplicates, so your test is failing, but it is because you don't have to (can't) have that requirement.

AssertionError: 10000 != 9996 : Duplicate values found in multi-threaded run, indicating a race condition.

I think only the comparison of sorted lists gives you a true idea of what happened.

[wenzeslaus]

...outcome order is guaranteed.

From library perspective, but never in the actual usage.

Consumption is another topic.

So I would say just don't overstate the outcome order guarantees in the doc.

I lean towards just stating that for fully reproducible outcome one should set parameters for single thread execution.

That would be my choice at this point as well.

Generally it is out of scope of this PR and should be discussed separately with more people involved.

I agree.

[wenzeslaus]

Commit 38ed9dc updates the PR to use C11 atomics instead of a pthread mutex to fix #5742 (thread-safe G_lrand48/G_mrand48/G_drand48). The generated sequence for a given seed is unchanged and now pinned by the regression tests from #7633 and #7634, which pass on my fork and should run also here.

Details: why atomics, how it works, and what was verified (AI-generated)

Three problems the mutex approach ran into:

1. The default builds never got the fix. HAVE_PTHREAD is only defined with --with-pthread, which is independent of --with-openmp and off by default. Our own macOS CI and both Docker images build with --with-openmp but without --with-pthread, so exactly the builds with the r.mapcalc: Support parallel computing by OpenMP #5742 race would have compiled the mutex out. The atomics gate (__STDC_VERSION__/__STDC_NO_ATOMICS__) is answered by the compiler, with no configure involvement.

2. Build system changes are no longer needed. The Makefile and CMake changes are dropped entirely (also resolving the CMake review comment). As a side effect this avoids defining HAVE_PTHREAD on MSVC, where find_package(Threads) succeeds with Win32 threads but <pthread.h> does not exist.

3. Performance. The mutex serializes every call and adds contention: in a microbenchmark, threaded runs were slower than single-threaded at every thread count. The atomic version keeps the exact same serialized sequence at about half the cost when several threads generate numbers at the same time, with near zero single-threaded overhead. (A single shared generator still cannot scale with threads; per-thread streams remain future work.)

How it works: the generator is the standard drand48 linear congruential generator (LCG): the next state is computed from the current one as X' = (A*X + B) mod 2^48. The previous code kept those 48 bits split across three 16-bit variables and advanced them with hand-written multi-word arithmetic, a leftover from times before 64-bit integers were portable; that split is what made concurrent updates able to interleave and corrupt the state. Now the whole state is a single 48-bit value in one atomic_uint_least64_t, and advancing it is one atomic compare-and-swap (CAS): a thread reads the state, computes the successor with the formula above, and atomically replaces the state only if it still holds the value that was read; otherwise it retries with the fresh value. In other words, the update optimistically assumes the state was not touched by another thread in the meantime and just recomputes when it was — unlike a lock, which pessimistically excludes all other threads ahead of time and pays that cost on every call whether a conflict occurs or not. Each successful swap is exactly one generator step, and because the caller computes its return value from the successor it installed rather than by re-reading shared variables, concurrent callers consume distinct, consecutive steps of the same sequence. Without C11 atomics the code falls back to the previous plain implementation (not thread-safe, i.e. unchanged behavior on such platforms; RFC 7 allows optional C11 features with a fallback).

The same C11-atomics-with-fallback approach is used in #7350; if both land, the capability test can be factored into a shared macro.

Reproducibility: the sequence of values is serialized and seed-reproducible, but which thread (and so which raster cell) receives which value depends on scheduling. Fully reproducible results still require single-threaded execution; this is now stated in the file documentation. The concern about G_fatal_error() under a held lock in the mutex-based implementation no longer applies (there is no lock).

Test cleanup: the duplicate-values assertion is removed (a correct random sequence can contain duplicates, as the observed 10000 != 9996 failure showed; the sorted-lists comparison is the valid check), the AI authorship note moved to the PR, and the banner comments are gone.

The 64-bit reformulation and its verification were AI-assisted (Claude). Verification included a standalone build of this file compared against the previous implementation for 100k steps and the reference sequences for five seeds, an 8-thread multiset check (1M draws), and the regression tests from #7633 and #7634 run against a GRASS build using this file.

[wenzeslaus]

@nilason This is now using your approach from #7350 while not touching any build files (unlike the pthreads locking version). Please review.

I changed the code to address my previous review.

[wenzeslaus]

This is passing tests from #7633 and #7634, i.e. generating same numbers as before. I approve, but I'm also responsible for the new code, so someone else should give the formal approval, too.

[nilason]

Look good to me, the atomic ops are reasonable.

[marisn]

Lets merge this, but remember that this is the first step, as we still have other issues to tackle:

• reproducibility in threaded environments;
• MSVC;
• fallback if compiled with older compiler (do we block older compilers?)

[nilason]

Lets merge this, but remember that this is the first step, as we still have other issues to tackle:

• reproducibility in threaded environments;

#7884 is an enhancement built upon this PR. I'm not sure if that answers your question.

• MSVC;

MSVC do have support for atomic operations, but with a bit different API compared to C11. I plan to suggest to add a GRASS wrapper API for atomic operations similar to https://github.com/python/cpython/pull/108701/changes. That will hopefully make implementation of atomic operations with broader support easier in GRASS.

• fallback if compiled with older compiler (do we block older compilers?)

cpython do have a wrapper for GCC's internal functions which match C11. This will enable support for GCC 4.8+ according to https://github.com/python/cpython/blob/main/Include/cpython/pyatomic_gcc.h.
As this PR currently stands, atomic operations (and thus thread safe random number generation) is enabled only if the compiler supports it.

[wenzeslaus]

With available C11 atomics, main is now without race although potentially slower (price for correctness). Follow-up in #7884 for a more complete fix.