From 8a94c90851c3af1a072f61a7a1dfc763eb70b7b9 Mon Sep 17 00:00:00 2001 From: Gert Van Hoey Date: Sat, 8 Aug 2026 09:26:27 +0000 Subject: [PATCH 1/4] Introduce a thread-barrier in GateTimeSorter to reduce memory consumption in MT simulations --- .../digitizer/GateCoincidenceSorterActor.cpp | 9 + .../digitizer/GateCoincidenceSorterActor.h | 2 + .../digitizer/GateDigitizerDeadTimeActor.cpp | 8 + .../digitizer/GateDigitizerDeadTimeActor.h | 4 +- .../digitizer/GateDigitizerPileupActor.cpp | 8 + .../digitizer/GateDigitizerPileupActor.h | 4 +- .../opengate_lib/digitizer/GateTimeSorter.cpp | 223 +++++++++++++++++- .../opengate_lib/digitizer/GateTimeSorter.h | 29 +++ opengate/actors/digitizers.py | 36 +++ 9 files changed, 317 insertions(+), 6 deletions(-) diff --git a/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.cpp b/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.cpp index 5bfa79c637..6b882537d5 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.cpp +++ b/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.cpp @@ -103,6 +103,13 @@ void GateCoincidenceSorterActor::InitializeUserInfo(py::dict &user_info) { if (py::len(user_info) > 0 && user_info.contains("sorting_time")) { fSortingTime = DictGetDouble(user_info, "sorting_time"); // nanoseconds } + if (py::len(user_info) > 0 && user_info.contains("sorting_buffer_size")) { + fSortingBufferSyncThreshold = DictGetInt(user_info, "sorting_buffer_size"); + } + if (py::len(user_info) > 0 && user_info.contains("thread_sync_enabled")) { + fThreadSyncEnabled = DictGetBool(user_info, "thread_sync_enabled"); + } + fGroupVolumeDepth = -1; fInputDigiCollectionName = DictGetStr(user_info, "input_digi_collection"); } @@ -183,6 +190,8 @@ void GateCoincidenceSorterActor::BeginOfRunActionMasterThread(int run_id) { fTimeSorter->Init(fInputDigiCollection); fTimeSorter->SetSortingWindow(fSortingTime); fTimeSorter->SetMaxSize(fClearEveryNEvents); + fTimeSorter->SetThreadSyncEnabled(fThreadSyncEnabled); + fTimeSorter->SetBufferThreadSyncThreshold(fSortingBufferSyncThreshold); // Create temporary storage, which stores time-sorted digis until a // sufficiently large GlobalTime range is covered to allow for coincidence diff --git a/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.h b/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.h index 921b88734a..3b0223f459 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.h +++ b/core/opengate_core/opengate_lib/digitizer/GateCoincidenceSorterActor.h @@ -63,6 +63,8 @@ class GateCoincidenceSorterActor : public GateVDigitizerWithOutputActor { TransaxialPlane fTransaxialPlane{TransaxialPlane::XY}; int fGroupVolumeDepth; double fSortingTime; + int fSortingBufferSyncThreshold; + bool fThreadSyncEnabled; struct TemporaryStorage { TemporaryStorage(GateDigiCollection *input, GateDigiCollection *output, diff --git a/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.cpp b/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.cpp index 963f4aacdf..401d43f403 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.cpp +++ b/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.cpp @@ -39,6 +39,12 @@ void GateDigitizerDeadTimeActor::InitializeUserInfo(py::dict &user_info) { if (py::len(user_info) > 0 && user_info.contains("sorting_time")) { fSortingTime = DictGetDouble(user_info, "sorting_time"); // nanoseconds } + if (py::len(user_info) > 0 && user_info.contains("sorting_buffer_size")) { + fSortingBufferSyncThreshold = DictGetInt(user_info, "sorting_buffer_size"); + } + if (py::len(user_info) > 0 && user_info.contains("thread_sync_enabled")) { + fThreadSyncEnabled = DictGetBool(user_info, "thread_sync_enabled"); + } fGroupVolumeDepth = -1; fInputDigiCollectionName = DictGetStr(user_info, "input_digi_collection"); } @@ -49,6 +55,8 @@ void GateDigitizerDeadTimeActor::BeginOfRunActionMasterThread(int run_id) { fTimeSorter->Init(fInputDigiCollection); fTimeSorter->SetSortingWindow(fSortingTime); fTimeSorter->SetMaxSize(fClearEveryNEvents); + fTimeSorter->SetThreadSyncEnabled(fThreadSyncEnabled); + fTimeSorter->SetBufferThreadSyncThreshold(fSortingBufferSyncThreshold); auto &outputIter = fTimeSorter->OutputIterator(); outputIter.TrackAttribute("GlobalTime", &fTimeSorterOutputTime); diff --git a/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.h b/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.h index b1ce464ec5..e47370c534 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.h +++ b/core/opengate_core/opengate_lib/digitizer/GateDigitizerDeadTimeActor.h @@ -53,8 +53,10 @@ class GateDigitizerDeadTimeActor : public GateVDigitizerWithOutputActor { // User parameters double fDeadTime; DeadTimePolicy fPolicy; - double fSortingTime; int fGroupVolumeDepth; + double fSortingTime; + int fSortingBufferSyncThreshold; + bool fThreadSyncEnabled; std::unique_ptr fTimeSorter; std::map fVolumeEndOfDeadTimeInterval; diff --git a/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.cpp b/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.cpp index 03d3d1ae1b..56084c4ea1 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.cpp +++ b/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.cpp @@ -66,6 +66,12 @@ void GateDigitizerPileupActor::InitializeUserInfo(py::dict &user_info) { if (py::len(user_info) > 0 && user_info.contains("sorting_time")) { fSortingTime = DictGetDouble(user_info, "sorting_time"); // nanoseconds } + if (py::len(user_info) > 0 && user_info.contains("sorting_buffer_size")) { + fSortingBufferSyncThreshold = DictGetInt(user_info, "sorting_buffer_size"); + } + if (py::len(user_info) > 0 && user_info.contains("thread_sync_enabled")) { + fThreadSyncEnabled = DictGetBool(user_info, "thread_sync_enabled"); + } fGroupVolumeDepth = -1; fInputDigiCollectionName = DictGetStr(user_info, "input_digi_collection"); } @@ -76,6 +82,8 @@ void GateDigitizerPileupActor::BeginOfRunActionMasterThread(int run_id) { fTimeSorter->Init(fInputDigiCollection); fTimeSorter->SetSortingWindow(fSortingTime); fTimeSorter->SetMaxSize(fClearEveryNEvents); + fTimeSorter->SetThreadSyncEnabled(fThreadSyncEnabled); + fTimeSorter->SetBufferThreadSyncThreshold(fSortingBufferSyncThreshold); auto &outputIter = fTimeSorter->OutputIterator(); outputIter.TrackAttribute("GlobalTime", &fTimeSorterOutputTime); diff --git a/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.h b/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.h index 0aba836383..118502df74 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.h +++ b/core/opengate_core/opengate_lib/digitizer/GateDigitizerPileupActor.h @@ -60,8 +60,10 @@ class GateDigitizerPileupActor : public GateVDigitizerWithOutputActor { TimeWindowPolicy fTimeWindowPolicy; AttributePolicy fAttributePolicy; PositionAttributePolicy fPositionAttributePolicy; - double fSortingTime; int fGroupVolumeDepth; + double fSortingTime; + int fSortingBufferSyncThreshold; + bool fThreadSyncEnabled; // Struct for storing digis in one particular volume which belong to the same // time window. diff --git a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp index 9c27600e24..bfed7bb165 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp +++ b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp @@ -14,6 +14,10 @@ #include #include +// Static variable that allows a time sorter to determine whether it is the most +// upstream time sorter instance in the sequence of actors in the simulation. +std::atomic GateTimeSorter::sMostUpstreamInstance{nullptr}; + GateTimeSorter::GateTimeSorter(const std::string &name) : fName(name) { fNumWorkingThreads = std::max(1, G4Threading::GetNumberOfRunningWorkerThreads()); @@ -22,6 +26,16 @@ GateTimeSorter::GateTimeSorter(const std::string &name) : fName(name) { std::make_unique(fNumWorkingThreads); } +GateTimeSorter::~GateTimeSorter() { + // Reset the static so that subsequent simulation runs can elect a new + // upstream instance. Without this, sMostUpstreamInstance would keep pointing + // to a destroyed object and the compare-and-swap (CAS) in IsFirstUpstream() + // would never succeed again (expected == nullptr would always fail). + GateTimeSorter *expected = this; + sMostUpstreamInstance.compare_exchange_strong( + expected, nullptr, std::memory_order_acq_rel, std::memory_order_relaxed); +} + void GateTimeSorter::Init(GateDigiCollection *input) { fInputCollection = input; @@ -129,19 +143,37 @@ void GateTimeSorter::SetMaxSize(size_t maxSize) { fMaxSize = maxSize; } +void GateTimeSorter::SetBufferThreadSyncThreshold(size_t size) { + if (fProcessingStarted) { + Fatal("SetBufferThreadSyncThreshold() cannot be called after Ingest() has " + "been called."); + } + fThreadSync.activationThreshold = size; +} + +void GateTimeSorter::SetThreadSyncEnabled(bool enabled) { + if (fProcessingStarted) { + Fatal("SetThreadSyncEnabled() cannot be called after Ingest() has " + "been called."); + } + fThreadSync.enabled = enabled; +} + void GateTimeSorter::OnEndOfEventAction(std::function work) { // This method is intended to be called by an actor in its EndOfEventAction() // method. The work function provided by the actor may then be called for // downstream processing of the time-sorted digis. - // There are two stages: + // There are three stages: // 1. The time sorter ingests the digis that are provided by the actor by // copying them into an ingestion buffer. This ingestion must be synchronized // with a mutex, because all threads must copy their digis sequentially into // the same buffer. - // 2. The actual time-sorting, followed by the work provided by the actor, + // 2. From time to time, the threads in a multi-threaded simulation require + // synchronization to reduce memory consumption in the time sorter. + // 3. The actual time-sorting, followed by the work provided by the actor, // happens next when a certain number of ingestions has happened since the // previous time sorting. This condition avoid incurring the overhead of stage - // 2 after every ingestion of (typically very few) digis. + // 3 after every ingestion of (typically very few) digis. // Phase 1 @@ -152,6 +184,25 @@ void GateTimeSorter::OnEndOfEventAction(std::function work) { // Phase 2 + // GlobalTime divergence between threads leads to increased memory consumption + // in the time sorter. A thread synchronization barrier reduces GlobalTime + // divergence by occasionally forcing the faster progressing threads to wait + // for the slower ones to catch up. The thread synchronization should be + // active in only one time sorter instance, to avoid deadlock. It should be + // active in the most upstream time sorter instance, since only there the + // GlobalTime observed by a thread is related to the primary events simulated + // by that same thread (digis may switch from one thread's input collection to + // another thread's output collection while being processed in the time + // sorter). The first time sorter to ingest at least one digi will be marked + // as the most upstream time sorter. + + if (ThreadSyncRequired()) { + SetupBarrierIfNeeded(); + WaitAtBarrierIfNeeded(); + } + + // Phase 3 + // If the number of ingestions has not yet reached the threshold, then // increment the counter and return. constexpr int numIngestionsPerProcessCall = 10; @@ -174,7 +225,12 @@ void GateTimeSorter::OnEndOfEventAction(std::function work) { std::memory_order_acquire, std::memory_order_relaxed)) { Process(); // executes time-sorting logic - work(); // executes the work provided by the actor + if (ThreadSyncRequired()) { + fThreadSync.sortedIndicesSize.store(fSortedIndicesA->size(), + std::memory_order_release); + fThreadSync.barrierSetupAllowed.store(true, std::memory_order_release); + } + work(); // executes the work provided by the actor fProcessingOngoing.store(false, std::memory_order_release); } } @@ -190,6 +246,12 @@ void GateTimeSorter::OnEndOfRunAction( // anyThreadWork. // lastThreadWork allows the calling actor to execute logic that is intended // to run after the GateTimeSorter has finalized all digi sorting. + + if (ThreadSyncRequired()) { + // Disable thread synchronization at the end of a run. + fThreadSync.barrierBypassed.store(true, std::memory_order_release); + fThreadSync.barrierConditionVariable.notify_all(); + } if (fNumActiveWorkingThreads.fetch_sub(1, std::memory_order_acq_rel) <= 1) { Process(); Flush(); @@ -278,6 +340,151 @@ bool GateTimeSorter::Ingest() { return true; } +bool GateTimeSorter::IsFirstUpstream() { + // The most upstream time sorter is the one that has succeeded in replacing + // the nullptr value with its own this pointer. + if (sMostUpstreamInstance.load(std::memory_order_relaxed) == this) { + return true; + } + // Try to replace nullptr by the this pointer value, and return true if it + // succeeded. + GateTimeSorter *expected = nullptr; + return sMostUpstreamInstance.compare_exchange_strong( + expected, this, std::memory_order_acq_rel, std::memory_order_relaxed); +} + +bool GateTimeSorter::ThreadSyncRequired() { + return fNumWorkingThreads > 1 && fThreadSync.enabled && IsFirstUpstream(); +} + +void GateTimeSorter::SetupBarrierIfNeeded() { + // * barrierSetupAllowed: set to true as soon as Process() has been executed + // after all threads have reached the barrier, indicating that + // sortedIndicesSize has been reduced. + // * barrierSetupClaimed: ensures that only one thread can set up the barrier. + // * barrierSetupComplete: indicates that the barrier is active. + + auto &ts = fThreadSync; + // As soon as the number of buffered digis reaches the activation threshold, + // one thread activates the barrier. + if (ts.barrierSetupAllowed.load(std::memory_order_relaxed) && + !ts.barrierSetupClaimed.load(std::memory_order_relaxed) && + ts.sortedIndicesSize.load(std::memory_order_relaxed) >= + fThreadSync.activationThreshold) { + bool expected = false; + if (ts.barrierSetupClaimed.compare_exchange_strong( + expected, true, std::memory_order_acq_rel, + std::memory_order_relaxed)) { + // Determine the highest GlabalTime value across all threads. + auto maxIt = std::max_element( + fMaxGlobalTimePerThread.get(), + fMaxGlobalTimePerThread.get() + fNumWorkingThreads, + [](const PaddedAtomicDouble &a, const PaddedAtomicDouble &b) { + return a.value.load() < b.value.load(); + }); + // Store target before the release on fBarrierSetupComplete so that + // threads that acquire fBarrierSetupComplete == true are guaranteed + // to see it. + const double maxTime = maxIt->value.load(); + ts.barrierGlobalTimeTarget.store(maxTime, std::memory_order_relaxed); + std::cout << "barrierGlobalTimeTarget " + << ts.barrierGlobalTimeTarget.load() << " sortedIndicesSize " + << ts.sortedIndicesSize.load(std::memory_order_relaxed) << "\n"; + ts.barrierSetupComplete.store(true, std::memory_order_release); + } + } +} + +void GateTimeSorter::WaitAtBarrierIfNeeded() { + auto &ts = fThreadSync; + // If the barrier has been set up and this thread's GlobalTime has reached + // the current target value, then wait until every other thread has also + // reached it. + if (ts.barrierSetupComplete.load(std::memory_order_acquire)) { + const int tid = std::max(0, G4Threading::G4GetThreadId()); + const double threadTime = + fMaxGlobalTimePerThread[tid].value.load(std::memory_order_relaxed); + const double targetTime = + ts.barrierGlobalTimeTarget.load(std::memory_order_acquire); + + if (threadTime >= targetTime) { + std::unique_lock cvLock(ts.barrierConditionVariableMutex); + + // Re-check under lock: the barrier may have been released while we + // were between the outer fBarrierSetupComplete check and here. + if (ts.barrierSetupComplete.load(std::memory_order_relaxed)) { + const int generation = + ts.barrierGeneration.load(std::memory_order_relaxed); + const int numArrived = + ts.numThreadsAtBarrier.fetch_add(1, std::memory_order_relaxed) + 1; + + // TODO remove this logging + if (numArrived == 1) { + auto [minIt, maxIt] = std::minmax_element( + fMaxGlobalTimePerThread.get(), + fMaxGlobalTimePerThread.get() + fNumWorkingThreads, + [](const PaddedAtomicDouble &a, const PaddedAtomicDouble &b) { + return a.value.load() < b.value.load(); + }); + std::cout << "[" << fName << "] Barrier hit with divergence " + << (maxIt->value.load() - minIt->value.load()) << " at " + << maxIt->value.load() << "\n"; + } + // TODO remove until here + + if (numArrived >= fNumWorkingThreads) { + // Last thread to arrive: reset state and release all waiters. + ts.numThreadsAtBarrier.store(0, std::memory_order_relaxed); + ts.barrierSetupComplete.store(false, std::memory_order_relaxed); + ts.barrierSetupClaimed.store(false, std::memory_order_relaxed); + ts.barrierSetupAllowed.store(false, std::memory_order_relaxed); + + // Reset sorting window to allow number of buffered digis to decrease + // in the next call to Process(). Do this before unlocking so woken + // threads see the updated value during Process(). + auto [minIt, maxIt] = std::minmax_element( + fMaxGlobalTimePerThread.get(), + fMaxGlobalTimePerThread.get() + fNumWorkingThreads, + [](const PaddedAtomicDouble &a, const PaddedAtomicDouble &b) { + return a.value.load() < b.value.load(); + }); + fSortingWindow.store(fMinimumSortingWindow + maxIt->value.load() - + minIt->value.load()); + + ts.barrierGeneration.fetch_add(1, std::memory_order_relaxed); + + // TODO remove this logging + { + auto [minIt, maxIt] = std::minmax_element( + fMaxGlobalTimePerThread.get(), + fMaxGlobalTimePerThread.get() + fNumWorkingThreads, + [](const PaddedAtomicDouble &a, const PaddedAtomicDouble &b) { + return a.value.load() < b.value.load(); + }); + std::cout << "[" << fName << "] Barrier lifted with divergence " + << (maxIt->value.load() - minIt->value.load()) << "\n"; + } + std::cout << "fBarrierGlobalTimeTarget " + << ts.barrierGlobalTimeTarget.load() << "\n"; + // TODO remove until here + + cvLock.unlock(); + // Last thread has arrived so all can resume their work. + ts.barrierConditionVariable.notify_all(); + } else { + // Park the thread until the barrier is released or the run ends. + // cvLock is already held; wait() releases it while parked. + ts.barrierConditionVariable.wait(cvLock, [&] { + return ts.barrierGeneration.load(std::memory_order_relaxed) != + generation || + ts.barrierBypassed.load(std::memory_order_relaxed); + }); + } + } + } + } +} + void GateTimeSorter::Process() { // Processes all digis from the ingestion buffer by copying and sorting them // according to GlobalTime. Next, copies the oldest sorted digis to the @@ -335,6 +542,7 @@ void GateTimeSorter::Process() { if (fMostRecentTimeDeparted.has_value() && (digiTime < *fMostRecentTimeDeparted)) { ++fNumDroppedDigi; + std::cout << "[" << fName << "] drop " << fNumDroppedDigi << "\n"; fMaxDropDelta = std::max(fMaxDropDelta, *fMostRecentTimeDeparted - digiTime); if (!fSortingWindowWarningIssued) { @@ -426,6 +634,8 @@ void GateTimeSorter::Prune() { // 2. Sorted collection A is cleared. // 3. The two collections and sorted index queues are swapped. + const auto numBefore = fSortedCollectionA->GetSize(); + // Step 1 GateDigiAttributesFiller transferFiller( fSortedCollectionA, fSortedCollectionB, @@ -449,4 +659,9 @@ void GateTimeSorter::Prune() { // Step 3 std::swap(fSortedCollectionA, fSortedCollectionB); std::swap(fSortedIndicesA, fSortedIndicesB); + + std::cout << "Prune [" << fName << "] " << numBefore << " " + << fSortedCollectionA->GetSize() << "\n"; + std::cout << "[" << fName << "] thread " + << std::max(0, G4Threading::G4GetThreadId()) << "\n"; } diff --git a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h index 4efc2014a3..756fc54b78 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h +++ b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h @@ -7,6 +7,7 @@ #include #include #include +#include #include #include @@ -16,11 +17,14 @@ class GateDigiAttributesFiller; class GateTimeSorter { public: GateTimeSorter(const std::string &name); + ~GateTimeSorter(); void Init(GateDigiCollection *input); void SetSortingWindow(double duration); void SetMaxSize(size_t size); + void SetBufferThreadSyncThreshold(size_t size); + void SetThreadSyncEnabled(bool enabled); void OnEndOfEventAction(std::function work); void OnEndOfRunAction(std::function anyThreadWork, @@ -32,6 +36,10 @@ class GateTimeSorter { private: bool Ingest(); + bool IsFirstUpstream(); + bool ThreadSyncRequired(); + void SetupBarrierIfNeeded(); + void WaitAtBarrierIfNeeded(); void Process(); void Flush(); @@ -102,6 +110,27 @@ class GateTimeSorter { double fMaxDropDelta{}; std::optional fMostRecentTimeArrived; std::optional fMostRecentTimeDeparted; + + static std::atomic sMostUpstreamInstance; + + struct ThreadSync { + bool enabled{true}; + size_t activationThreshold{50'000}; + + std::atomic barrierSetupAllowed{true}; + std::atomic barrierSetupClaimed{false}; + std::atomic barrierSetupComplete{false}; + std::atomic barrierBypassed{false}; + std::atomic barrierGlobalTimeTarget{0.0}; + std::atomic numThreadsAtBarrier{0}; + std::atomic barrierGeneration{0}; + + std::atomic sortedIndicesSize{0}; + std::mutex barrierConditionVariableMutex; + std::condition_variable barrierConditionVariable; + }; + + ThreadSync fThreadSync; }; #endif // GateTimeSorter_h diff --git a/opengate/actors/digitizers.py b/opengate/actors/digitizers.py index 67c5455a06..6b21a6627b 100644 --- a/opengate/actors/digitizers.py +++ b/opengate/actors/digitizers.py @@ -620,6 +620,18 @@ class DigitizerDeadTimeActor(DigitizerWithRootOutput, g4.GateDigitizerDeadTimeAc "doc": "Time interval during which digis are buffered for time-sorting", }, ), + "thread_sync_enabled": ( + True, + { + "doc": "Thread synchronization is applied to limit sorting buffer size", + }, + ), + "sorting_buffer_size": ( + 5e4, + { + "doc": "Sorting buffer size (digis) at which thread synchronization is applied", + }, + ), "skip_attributes": ( [], { @@ -742,6 +754,18 @@ class DigitizerPileupActor(DigitizerWithRootOutput, g4.GateDigitizerPileupActor) "doc": "Time interval during which digis are buffered for time-sorting", }, ), + "thread_sync_enabled": ( + True, + { + "doc": "Thread synchronization is applied to limit sorting buffer size", + }, + ), + "sorting_buffer_size": ( + 5e4, + { + "doc": "Sorting buffer size (digis) at which thread synchronization is applied", + }, + ), "skip_attributes": ( [], { @@ -1467,6 +1491,18 @@ class CoincidenceSorterActor(DigitizerWithRootOutput, g4.GateCoincidenceSorterAc "doc": "Time interval during which digis are buffered for time-sorting", }, ), + "thread_sync_enabled": ( + True, + { + "doc": "Thread synchronization is applied to limit sorting buffer size", + }, + ), + "sorting_buffer_size": ( + 5e4, + { + "doc": "Sorting buffer size (digis) at which thread synchronization is applied", + }, + ), "skip_attributes": ( [], { From 8386365a14ffd89a28bde436fb375084836f448f Mon Sep 17 00:00:00 2001 From: Gert Van Hoey Date: Sat, 8 Aug 2026 15:45:32 +0000 Subject: [PATCH 2/4] Improve upstream actor detection --- .../opengate_lib/digitizer/GateTimeSorter.cpp | 17 +++++++++++++---- .../opengate_lib/digitizer/GateTimeSorter.h | 1 + 2 files changed, 14 insertions(+), 4 deletions(-) diff --git a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp index bfed7bb165..211c679301 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp +++ b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp @@ -343,14 +343,23 @@ bool GateTimeSorter::Ingest() { bool GateTimeSorter::IsFirstUpstream() { // The most upstream time sorter is the one that has succeeded in replacing // the nullptr value with its own this pointer. - if (sMostUpstreamInstance.load(std::memory_order_relaxed) == this) { + if (fIsFirstUpstream.load(std::memory_order_relaxed)) { return true; } - // Try to replace nullptr by the this pointer value, and return true if it + // Once any instance is elected, there is no point retrying the CAS. + if (sMostUpstreamInstance.load(std::memory_order_relaxed) != nullptr) { + return false; + } + // Try to replace nullptr by the this pointer value and return true if it // succeeded. GateTimeSorter *expected = nullptr; - return sMostUpstreamInstance.compare_exchange_strong( - expected, this, std::memory_order_acq_rel, std::memory_order_relaxed); + if (sMostUpstreamInstance.compare_exchange_strong( + expected, this, std::memory_order_acq_rel, + std::memory_order_relaxed)) { + fIsFirstUpstream.store(true, std::memory_order_release); + return true; + } + return false; } bool GateTimeSorter::ThreadSyncRequired() { diff --git a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h index 756fc54b78..d426ad8947 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h +++ b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.h @@ -112,6 +112,7 @@ class GateTimeSorter { std::optional fMostRecentTimeDeparted; static std::atomic sMostUpstreamInstance; + std::atomic fIsFirstUpstream{false}; struct ThreadSync { bool enabled{true}; From 6680ccbec1ff595dcf5a405e3159a0039390c667 Mon Sep 17 00:00:00 2001 From: Gert Van Hoey Date: Sat, 8 Aug 2026 15:04:14 +0000 Subject: [PATCH 3/4] Remove debug logging from GateTimeSorter --- .../opengate_lib/digitizer/GateTimeSorter.cpp | 41 ------------------- 1 file changed, 41 deletions(-) diff --git a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp index 211c679301..b6e2336cda 100644 --- a/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp +++ b/core/opengate_core/opengate_lib/digitizer/GateTimeSorter.cpp @@ -396,9 +396,6 @@ void GateTimeSorter::SetupBarrierIfNeeded() { // to see it. const double maxTime = maxIt->value.load(); ts.barrierGlobalTimeTarget.store(maxTime, std::memory_order_relaxed); - std::cout << "barrierGlobalTimeTarget " - << ts.barrierGlobalTimeTarget.load() << " sortedIndicesSize " - << ts.sortedIndicesSize.load(std::memory_order_relaxed) << "\n"; ts.barrierSetupComplete.store(true, std::memory_order_release); } } @@ -427,20 +424,6 @@ void GateTimeSorter::WaitAtBarrierIfNeeded() { const int numArrived = ts.numThreadsAtBarrier.fetch_add(1, std::memory_order_relaxed) + 1; - // TODO remove this logging - if (numArrived == 1) { - auto [minIt, maxIt] = std::minmax_element( - fMaxGlobalTimePerThread.get(), - fMaxGlobalTimePerThread.get() + fNumWorkingThreads, - [](const PaddedAtomicDouble &a, const PaddedAtomicDouble &b) { - return a.value.load() < b.value.load(); - }); - std::cout << "[" << fName << "] Barrier hit with divergence " - << (maxIt->value.load() - minIt->value.load()) << " at " - << maxIt->value.load() << "\n"; - } - // TODO remove until here - if (numArrived >= fNumWorkingThreads) { // Last thread to arrive: reset state and release all waiters. ts.numThreadsAtBarrier.store(0, std::memory_order_relaxed); @@ -461,22 +444,6 @@ void GateTimeSorter::WaitAtBarrierIfNeeded() { minIt->value.load()); ts.barrierGeneration.fetch_add(1, std::memory_order_relaxed); - - // TODO remove this logging - { - auto [minIt, maxIt] = std::minmax_element( - fMaxGlobalTimePerThread.get(), - fMaxGlobalTimePerThread.get() + fNumWorkingThreads, - [](const PaddedAtomicDouble &a, const PaddedAtomicDouble &b) { - return a.value.load() < b.value.load(); - }); - std::cout << "[" << fName << "] Barrier lifted with divergence " - << (maxIt->value.load() - minIt->value.load()) << "\n"; - } - std::cout << "fBarrierGlobalTimeTarget " - << ts.barrierGlobalTimeTarget.load() << "\n"; - // TODO remove until here - cvLock.unlock(); // Last thread has arrived so all can resume their work. ts.barrierConditionVariable.notify_all(); @@ -551,7 +518,6 @@ void GateTimeSorter::Process() { if (fMostRecentTimeDeparted.has_value() && (digiTime < *fMostRecentTimeDeparted)) { ++fNumDroppedDigi; - std::cout << "[" << fName << "] drop " << fNumDroppedDigi << "\n"; fMaxDropDelta = std::max(fMaxDropDelta, *fMostRecentTimeDeparted - digiTime); if (!fSortingWindowWarningIssued) { @@ -643,8 +609,6 @@ void GateTimeSorter::Prune() { // 2. Sorted collection A is cleared. // 3. The two collections and sorted index queues are swapped. - const auto numBefore = fSortedCollectionA->GetSize(); - // Step 1 GateDigiAttributesFiller transferFiller( fSortedCollectionA, fSortedCollectionB, @@ -668,9 +632,4 @@ void GateTimeSorter::Prune() { // Step 3 std::swap(fSortedCollectionA, fSortedCollectionB); std::swap(fSortedIndicesA, fSortedIndicesB); - - std::cout << "Prune [" << fName << "] " << numBefore << " " - << fSortedCollectionA->GetSize() << "\n"; - std::cout << "[" << fName << "] thread " - << std::max(0, G4Threading::G4GetThreadId()) << "\n"; } From 16de9b3a88b0459eefebabf88f721c98e257d995 Mon Sep 17 00:00:00 2001 From: Gert Van Hoey Date: Sun, 20 Sep 2026 12:07:06 +0000 Subject: [PATCH 4/4] Update developer documentation --- .../developer_guide_coincidence_sorter_mt.rst | 179 ++++++++++++++++-- 1 file changed, 161 insertions(+), 18 deletions(-) diff --git a/docs/source/developer_guide/developer_guide_coincidence_sorter_mt.rst b/docs/source/developer_guide/developer_guide_coincidence_sorter_mt.rst index dc2e4114a9..6143f4d172 100644 --- a/docs/source/developer_guide/developer_guide_coincidence_sorter_mt.rst +++ b/docs/source/developer_guide/developer_guide_coincidence_sorter_mt.rst @@ -32,17 +32,24 @@ Component Overview ``GateTimeSorter`` (``GateTimeSorter.h/.cpp``) A self-contained merge-sorter. It is not a GATE actor but a helper object - owned by ``GateCoincidenceSorterActor``. It manages several internal digi - collections (see `Data-Flow Diagram`_) and exposes three operations that - callers use via the higher-level wrappers ``OnEndOfEventAction()`` and - ``OnEndOfRunAction()``: + owned by a digitizer actor. ``GateDigitizerDeadTimeActor``, + ``GateDigitizerPileupActor``, and ``GateCoincidenceSorterActor`` each have + their own ``GateTimeSorter`` instance, and several of these actors can be + chained one after another in the same digitizer pipeline. It manages + several internal digi collections (see `Data-Flow Diagram`_) and exposes + operations that callers use via the higher-level wrappers + ``OnEndOfEventAction()`` and ``OnEndOfRunAction()``: * ``Ingest()`` — copies the calling thread's digis into a shared ingestion buffer. Always mutex-protected; deliberately minimal. + * A thread-synchronization barrier (optional, see + `Phase 2 — Thread Synchronization Barrier`_) that periodically forces + threads to catch up with each other, bounding the memory growth caused + by ``GlobalTime`` divergence between threads. Active only in the + *most upstream* ``GateTimeSorter`` instance of the pipeline. * ``Process()`` — sorts the buffered digis and drains the oldest ones to - the output collection. Only one thread executes this at a time, and - only the thread currently tracking the highest ``GlobalTime`` is - selected (see `Phase 2 — Sorting and Output`_). + the output collection. Only one thread executes this at a time, + selected via a compare-and-swap race (see `Phase 3 — Sorting and Output`_). * ``Flush()`` — drains all remaining sorted digis at end-of-run, without applying the sorting window. @@ -113,6 +120,11 @@ memory-reclamation paths. | swapped | +--------------------------------------------------------------------+ +Not shown above: between ``Ingest()`` and the drain step, the (optional) +thread-synchronization barrier periodically shrinks ``fSortingWindow`` from +the outside, which is what keeps ``fSortedCollectionA``/``fSortedIndicesA`` +bounded in the long run (see `Phase 2 — Thread Synchronization Barrier`_). + GateTimeSorter: Design and Threading Model ------------------------------------------ @@ -142,7 +154,71 @@ window extension. Threads therefore block each other for only a very short time. -Phase 2 — Sorting and Output +Phase 2 — Thread Synchronization Barrier +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +``GlobalTime`` divergence between threads is the main driver of memory growth +in ``GateTimeSorter``: ``fSortingWindow`` extends to cover the full spread of +``fMaxGlobalTimePerThread`` (see `The Sorting Window`_), and a wide window +keeps digis buffered in ``fSortedCollectionA``/``fSortedIndicesA`` for longer. +A thread-synchronization barrier bounds this growth by occasionally forcing +threads that have progressed further to wait for the slower ones to catch up, +which shrinks the observed divergence and, in turn, ``fSortingWindow``. + +This mechanism must run in only *one* place in the digitizer pipeline — if +two chained actors both maintained a barrier, a thread parked at the +downstream barrier could prevent another thread from ever reaching the +upstream one, causing a deadlock. It is therefore active only when all of +the following hold, checked by ``ThreadSyncRequired()``: + +* the simulation is multi-threaded (``fNumWorkingThreads > 1``); +* it is enabled by the user (``thread_sync_enabled``, default ``true``); +* this instance is the *most upstream* ``GateTimeSorter`` in the pipeline — + determined by ``IsFirstUpstream()``, which has every constructed + ``GateTimeSorter`` race to replace the static + ``sMostUpstreamInstance`` pointer (initially ``nullptr``) with its own + ``this`` via CAS. The first instance to ingest a digi wins the race and + every subsequent call short-circuits via a cached ``fIsFirstUpstream`` + flag. The destructor resets ``sMostUpstreamInstance`` back to ``nullptr`` + if it still points to the instance being destroyed, so a new simulation run + can elect a fresh upstream instance. + +Each call to ``OnEndOfEventAction()`` (after ``Ingest()``) calls +``SetupBarrierIfNeeded()`` followed by ``WaitAtBarrierIfNeeded()``: + +``SetupBarrierIfNeeded()`` + Once ``fSortedIndicesA`` has grown to ``sorting_buffer_size`` digis + (``fThreadSync.activationThreshold``, default 50 000) since the barrier + was last released, one thread wins a CAS on ``barrierSetupClaimed``, + computes the highest ``GlobalTime`` currently reached by any thread + (``fMaxGlobalTimePerThread``), stores it as ``barrierGlobalTimeTarget``, + and publishes the barrier with a release store to ``barrierSetupComplete``. + +``WaitAtBarrierIfNeeded()`` + Any thread whose own ``fMaxGlobalTimePerThread[tid]`` has already reached + ``barrierGlobalTimeTarget`` blocks on ``barrierConditionVariable`` until + every working thread has arrived. The last thread to arrive resets the + barrier state (``numThreadsAtBarrier``, ``barrierSetupComplete``, + ``barrierSetupClaimed``, ``barrierSetupAllowed``), recomputes + ``fSortingWindow`` from the now much smaller spread of + ``fMaxGlobalTimePerThread`` — this is what actually lets + ``Process()`` drain more digis on the next call — advances + ``barrierGeneration``, and notifies all waiters. Waiting threads recheck + their captured ``barrierGeneration`` value against the current one to + guard against spurious and lost wakeups. + +``barrierSetupAllowed`` is only set back to ``true`` once ``Process()`` has +run after a barrier release (see below), which prevents a new barrier from +being set up again at (approximately) the same target before the reduced +``fSortingWindow`` has had a chance to take effect. + +At ``OnEndOfRunAction()``, before decrementing ``fNumActiveWorkingThreads``, +the calling thread sets ``barrierBypassed`` and notifies all waiters, so that +any thread still parked at the barrier when the run ends is released instead +of deadlocking against threads that have already left the event loop. + + +Phase 3 — Sorting and Output ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ After returning from ``Ingest()``, ``OnEndOfEventAction()`` increments the @@ -170,6 +246,13 @@ If the CAS succeeds, the thread calls ``Process()`` and ``work()``: .. code-block:: cpp Process(); // executes time-sorting logic + if (ThreadSyncRequired()) { + // Publish the (now smaller) buffer size and re-allow a new barrier + // to be set up on a subsequent call. + fThreadSync.sortedIndicesSize.store(fSortedIndicesA->size(), + std::memory_order_release); + fThreadSync.barrierSetupAllowed.store(true, std::memory_order_release); + } work(); // actor lambda: ProcessTimeSortedSingles + DetectCoincidences fProcessingOngoing.store(false, std::memory_order_release); @@ -223,8 +306,12 @@ thread observes a new maximum ``GlobalTime`` difference across threads: + max(fMaxGlobalTimePerThread) - min(fMaxGlobalTimePerThread)) -The window grows monotonically and never shrinks, so a temporary burst of -thread divergence leaves a permanent safety margin. +Outside of the barrier, the window grows monotonically and never shrinks, so +a temporary burst of thread divergence leaves a permanent safety margin. The +only place where ``fSortingWindow`` is reduced is the release of the +thread-synchronization barrier (see `Phase 2 — Thread Synchronization Barrier`_), +which recomputes it from the divergence observed *after* threads have caught +up with each other. End of Run — Flush @@ -248,7 +335,11 @@ reads and writes ``fIngestionBufferB``, ``fSortedCollectionA``, and ``fMostRecentTimeDeparted`` **without holding any mutex**, and reads ``fSortingWindow`` via an atomic load (``fSortingWindow`` is ``std::atomic`` because ``Ingest()`` -writes to it concurrently). This section explains why that is safe. +and the thread-synchronization barrier both write to it concurrently). This +section explains why that is safe. The thread-synchronization barrier itself +does use blocking primitives (a mutex and a condition variable); its +thread-safety is discussed separately in +`Thread-Safety of the Synchronization Barrier`_. Exclusive Ownership of ``fIngestionBufferB`` After the Swap @@ -328,6 +419,45 @@ Concurrent access is safe because: store, immediately visible to any subsequent atomic load in ``Process()``. +Thread-Safety of the Synchronization Barrier +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +Unlike the rest of ``GateTimeSorter``, the barrier (see +`Phase 2 — Thread Synchronization Barrier`_) genuinely blocks threads using a +mutex and a condition variable. This is acceptable because it only runs in +the single, most-upstream ``GateTimeSorter`` instance and only after a large +number of digis has accumulated, so it is far off the per-event hot path. + +* **Upstream election.** ``sMostUpstreamInstance`` is a + ``std::atomic``, initialised to ``nullptr``. Every + instance's first call to ``IsFirstUpstream()`` races to CAS it from + ``nullptr`` to ``this``; the winner caches the result in + ``fIsFirstUpstream`` so later calls avoid the atomic entirely. The + destructor CASes the pointer back from ``this`` to ``nullptr``, but only if + it still holds this instance's address — protecting against a later + instance having already claimed the slot for a subsequent run. +* **Single barrier setup per cycle.** ``barrierSetupClaimed`` is a + compare-and-swap gate: only the thread that flips it from ``false`` to + ``true`` computes ``barrierGlobalTimeTarget`` and publishes the barrier. + The release store to ``barrierSetupComplete`` happens after the target is + written, and every thread reads ``barrierSetupComplete`` with acquire + ordering before reading the target, so the target is always seen correctly + once the barrier is visible. +* **No lost or spurious wakeups.** Waiting threads capture the current + ``barrierGeneration`` before parking and re-check it inside the condition + variable's predicate. Because the last thread to arrive increments + ``barrierGeneration`` under the same mutex used for waiting, a wakeup is + only acted upon once the generation has actually advanced (or the run has + ended, see below), which rules out both premature wakeups and wakeups lost + between the predicate check and the call to ``wait()``. +* **No deadlock at end of run.** ``OnEndOfRunAction()`` sets + ``barrierBypassed`` and notifies all waiters before any thread decrements + ``fNumActiveWorkingThreads``. Waiting threads also wake up when + ``barrierBypassed`` becomes true, so a thread that reaches the barrier for + the last time it will ever call ``OnEndOfEventAction()`` cannot be left + parked forever. + + GateCoincidenceSorterActor: Consuming the Sorted Stream ------------------------------------------------------- @@ -335,7 +465,13 @@ GateCoincidenceSorterActor: Consuming the Sorted Stream (sorting window = ``fSortingTime``) and two ``TemporaryStorage`` objects (``fCurrentStorage`` and ``fFutureStorage``). Each ``TemporaryStorage`` holds a ``GateDigiCollection`` of time-sorted singles together with the -``GlobalTime`` of the earliest and latest digi it currently contains. +``GlobalTime`` of the earliest and latest digi it currently contains. It also +forwards the user-facing ``thread_sync_enabled`` and ``sorting_buffer_size`` +parameters to the time sorter via ``SetThreadSyncEnabled()`` and +``SetBufferThreadSyncThreshold()`` (see +`Phase 2 — Thread Synchronization Barrier`_). ``GateDigitizerDeadTimeActor`` +and ``GateDigitizerPileupActor`` expose the same two parameters and wire them +up identically for their own ``GateTimeSorter`` instance. At each ``EndOfEventAction``, the actor passes a lambda to ``fTimeSorter->OnEndOfEventAction()``. When the time sorter decides that the @@ -395,14 +531,21 @@ Thread Lifecycle Summary - All workers - ``Ingest()`` — acquires ``fIngestionMutex`` briefly (copy, atomic store to per-thread max, occasional ``fSortingWindow`` extension). Returns - early if no digis; otherwise increment ``fNumIngestions`` (relaxed - atomic). If threshold reached: attempt CAS on ``fProcessingOngoing`` - (atomic, non-blocking). If CAS succeeds: ``Process()`` — acquires - ``fIngestionMutex`` for buffer swap only, then fully lock-free — - followed by actor lambda (lock-free). + early if no digis. If this is the most-upstream, multi-threaded, + thread-sync-enabled time sorter: ``SetupBarrierIfNeeded()`` / + ``WaitAtBarrierIfNeeded()`` — one thread may configure the barrier + (CAS), then any thread that has reached the target ``GlobalTime`` may + block on a condition variable until all threads have arrived. + Otherwise increment ``fNumIngestions`` (relaxed atomic). If threshold + reached: attempt CAS on ``fProcessingOngoing`` (atomic, non-blocking). + If CAS succeeds: ``Process()`` — acquires ``fIngestionMutex`` for + buffer swap only, then fully lock-free — followed by actor lambda + (lock-free). * - ``EndOfRunAction`` - All workers - - Decrement ``fNumActiveWorkingThreads`` (atomic). Last thread: + - If thread sync is active: set ``barrierBypassed`` and notify all + waiters, so no thread remains parked at the barrier. + Decrement ``fNumActiveWorkingThreads`` (atomic). Last thread: ``Process()`` then ``Flush()`` (both single-threaded at this point) then ``lastThreadWork()``. All threads: ``anyThreadWork()``.