From 9bea3f3e0380e12f74c74a9d0ee70c63b4645694 Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Sun, 8 Nov 2026 22:28:24 +0100 Subject: [PATCH 1/8] Add a multithreaded Geant4 mode to o2-sim This adds the --isMT mode, in which one sim process runs Geant4 with one transport thread per sim worker. - Each Geant4 worker thread asks the primary server for its sub-events and sends its hits to the merger, as a worker process does. - The VMC application is cloned per thread, with its own stack, random generator, external decayer and data channel. - The threads share one request channel to the primary server. - Their data channels are closed before the process exits, so that the merger receives all queued hits. - A worker leaves its event loop once the server has no more work, after an empty event. - Detector clones get their own hit buffers, and the worker volume map points at the worker's own modules. - The media energy cuts and the worker-only Geant4 commands are applied on the worker threads. - The cached magnetic field values are invalidated at the start of each track. - ShmManager takes a mutex for its allocations, and the single sim process occupies the whole shared region. - 16 detectors are enabled for MT; EMCAL, HMPID, ZDC, FOCAL, FT0 and FV0 needed fixes to their copy constructor or destructor. - --isMT stops with an error for other engines, for Geant4 built without multithreading, for VecGeom navigation, and for detectors not yet enabled for it. https://its.cern.ch/jira/browse/O2-7305 https://github.com/vmc-project/geant4_vmc/pull/92 Co-Authored-By: Claude Opus 5.5 --- Common/Utils/include/CommonUtils/ShmManager.h | 8 +- Common/Utils/src/ShmManager.cxx | 12 +- .../FairRootBase/include/FairMCApplication.h | 13 +- .../FairRootBase/include/FairModule.h | 4 + .../FairRootBase/src/FairMCApplication.cxx | 32 ++- Detectors/Base/CMakeLists.txt | 8 + .../Base/include/DetectorsBase/Detector.h | 7 + .../include/DetectorsBase/MTWorkerHooks.h | 52 +++++ Detectors/Base/src/MTWorkerHooks.cxx | 23 +++ Detectors/Base/src/Stack.cxx | 10 +- Detectors/Base/test/testDetImplClone.cxx | 98 ++++++++++ .../include/CPVSimulation/Detector.h | 1 + .../include/EMCALSimulation/Detector.h | 1 + Detectors/EMCAL/simulation/src/Detector.cxx | 4 +- .../include/FDDSimulation/Detector.h | 1 + .../include/FT0Simulation/Detector.h | 1 + Detectors/FIT/FT0/simulation/src/Detector.cxx | 3 + .../include/FV0Simulation/Detector.h | 1 + Detectors/FIT/FV0/simulation/src/Detector.cxx | 3 - .../include/FOCALSimulation/Detector.h | 1 + Detectors/FOCAL/simulation/src/Detector.cxx | 7 +- .../include/HMPIDSimulation/Detector.h | 1 + Detectors/HMPID/simulation/src/Detector.cxx | 3 +- .../include/ITSSimulation/Detector.h | 2 + .../include/MFTSimulation/Detector.h | 1 + .../include/MCHSimulation/Detector.h | 1 + .../include/MIDSimulation/Detector.h | 1 + .../include/PHOSSimulation/Detector.h | 1 + .../include/TOFSimulation/Detector.h | 1 + .../include/TPCSimulation/Detector.h | 1 + .../include/TRDSimulation/Detector.h | 1 + .../include/ZDCSimulation/Detector.h | 4 + Detectors/ZDC/simulation/src/Detector.cxx | 31 ++- Detectors/gconfig/CMakeLists.txt | 2 +- Detectors/gconfig/commonConfig.C | 9 +- Detectors/gconfig/g4Config.C | 50 +++++ Detectors/gconfig/g4config.in | 1 - .../gconfig/include/SimSetup/G4FieldCache.h | 23 +++ Detectors/gconfig/src/G4FieldCache.cxx | 34 ++++ Steer/include/Steer/O2MCApplication.h | 37 +++- Steer/include/Steer/O2MCApplicationBase.h | 11 ++ Steer/src/O2MCApplication.cxx | 183 +++++++++++++++++- run/O2SimDevice.cxx | 29 +-- run/O2SimDeviceRunner.cxx | 145 ++++++++++++++ run/o2sim_parallel.cxx | 10 +- 45 files changed, 810 insertions(+), 62 deletions(-) create mode 100644 Detectors/Base/include/DetectorsBase/MTWorkerHooks.h create mode 100644 Detectors/Base/src/MTWorkerHooks.cxx create mode 100644 Detectors/Base/test/testDetImplClone.cxx create mode 100644 Detectors/gconfig/include/SimSetup/G4FieldCache.h create mode 100644 Detectors/gconfig/src/G4FieldCache.cxx diff --git a/Common/Utils/include/CommonUtils/ShmManager.h b/Common/Utils/include/CommonUtils/ShmManager.h index e7174b909b82c..bb0b8d1350575 100644 --- a/Common/Utils/include/CommonUtils/ShmManager.h +++ b/Common/Utils/include/CommonUtils/ShmManager.h @@ -22,6 +22,7 @@ #include #include #include +#include #if !defined(__CINT__) && !defined(__MAKECINT__) && !defined(__ROOTCLING__) && !defined(__CLING__) #include @@ -70,7 +71,8 @@ class ShmManager // create the local segment // this will occupy a subregion of an already created global shared mem segment - void occupySegment(); + // (all of it for the single process of multithreaded Geant4, whose threads share it) + void occupySegment(bool wholeRegion = false); // simply attaches to the global segment bool attachToGlobalSegment(); @@ -88,7 +90,7 @@ class ShmManager // returns if pointer is part of the shm region under control of this manager bool isPointerOk(void* ptr) const { - return mBufferPtr && getPointerOffset(ptr) < SHMPOOLSIZE; + return mBufferPtr && getPointerOffset(ptr) < mPoolSize; } // returns if shared mem setup is correctly setup/operational @@ -114,6 +116,8 @@ class ShmManager ShmMetaInfo* mSegInfoPtr = nullptr; // pointing to the meta information object bool mIsMaster = false; // true if the manager who allocated the region bool mIsOperational = false; + size_t mPoolSize = SHMPOOLSIZE; // size of the region occupied by this process + std::mutex mAllocMutex; // the boost buffer has no locking of its own; threads of multithreaded Geant4 share it // helper function void* tryAttach(bool& success); size_t getPointerOffset(void* ptr) const { return (size_t)((char*)ptr - (char*)mBufferPtr); } diff --git a/Common/Utils/src/ShmManager.cxx b/Common/Utils/src/ShmManager.cxx index 3ed863dcc96c3..32c11097bb1ef 100644 --- a/Common/Utils/src/ShmManager.cxx +++ b/Common/Utils/src/ShmManager.cxx @@ -178,7 +178,7 @@ bool ShmManager::attachToGlobalSegment() return b; } -void ShmManager::occupySegment() +void ShmManager::occupySegment(bool wholeRegion) { LOG(info) << "OCCUPYING A SEGMENT IN A SHARED REGION"; // get information about the global shared mem in which to occupy a region @@ -197,15 +197,17 @@ void ShmManager::occupySegment() if (b) { // read meta information in the segment to determine the id of this segment auto info = static_cast(addr); - const int segmentcounter = info->counter.fetch_add(1); + const int npools = wholeRegion ? (info->allocedbytes - sizeof(ShmMetaInfo)) / SHMPOOLSIZE : 1; + const int segmentcounter = info->counter.fetch_add(npools); LOG(info) << "SEGMENTCOUNT " << segmentcounter; const auto offset_in_bytes = sizeof(ShmMetaInfo) + segmentcounter * SHMPOOLSIZE; mBufferPtr = (void*)(((char*)addr) + offset_in_bytes); + mPoolSize = npools * SHMPOOLSIZE; - assert((unsigned long long)((char*)mBufferPtr - (char*)addr) + SHMPOOLSIZE <= info->allocedbytes); + assert((unsigned long long)((char*)mBufferPtr - (char*)addr) + mPoolSize <= info->allocedbytes); - boostmanagedbuffer = new boost::interprocess::wmanaged_external_buffer(create_only, mBufferPtr, SHMPOOLSIZE); + boostmanagedbuffer = new boost::interprocess::wmanaged_external_buffer(create_only, mBufferPtr, mPoolSize); boostallocator = new boost::interprocess::allocator( boostmanagedbuffer->get_segment_manager()); @@ -225,6 +227,7 @@ ShmManager::~ShmManager() void* ShmManager::getmemblock(size_t size) { void* addr = nullptr; + std::lock_guard lock(mAllocMutex); try { addr = (void*)boostallocator->allocate(size).get(); } catch (const std::exception& e) { @@ -235,6 +238,7 @@ void* ShmManager::getmemblock(size_t size) void ShmManager::freememblock(void* ptr, size_t s) { + std::lock_guard lock(mAllocMutex); boostallocator->deallocate((char*)ptr, s); } diff --git a/Detectors/3rdparty/FairRootBase/include/FairMCApplication.h b/Detectors/3rdparty/FairRootBase/include/FairMCApplication.h index 6ee3dd0488f15..e2a201c554845 100644 --- a/Detectors/3rdparty/FairRootBase/include/FairMCApplication.h +++ b/Detectors/3rdparty/FairRootBase/include/FairMCApplication.h @@ -57,6 +57,7 @@ enum class FairMCApplicationState * @since 12.01.04 */ +class FairModule; class FairMCApplication : public TVirtualMCApplication { public: @@ -304,9 +305,11 @@ class FairMCApplication : public TVirtualMCApplication ClassDef(FairMCApplication, 4); - private: - /** Protected copy constructor */ + protected: + /** Copy constructor, used to create the application of a worker thread */ FairMCApplication(const FairMCApplication&); + + private: /** Protected assignment operator */ FairMCApplication& operator=(const FairMCApplication&); @@ -314,6 +317,12 @@ class FairMCApplication : public TVirtualMCApplication Bool_t fGeometryIsInitialized; static FairMCApplication* fgMasterInstance; + + protected: + /** True for the application on the master thread, false for a worker clone */ + Bool_t IsMasterApplication() const { return !fgMasterInstance || fgMasterInstance == this; } + /** The worker clone of a module of the master application (the module itself on the master) */ + FairModule* OwnModule(FairModule* module) const; }; // inline functions diff --git a/Detectors/3rdparty/FairRootBase/include/FairModule.h b/Detectors/3rdparty/FairRootBase/include/FairModule.h index fdf8eae1267bd..ae6b330117357 100644 --- a/Detectors/3rdparty/FairRootBase/include/FairModule.h +++ b/Detectors/3rdparty/FairRootBase/include/FairModule.h @@ -159,6 +159,10 @@ class FairModule : public TNamed Bool_t fGeoSaved; //! flag for initialisation TVirtualMC* fMC; //! cahed pointer to MC (available only after initialization) + public: + /** Sets the cached MC, for a module cloned for a worker thread before its MC existed */ + void SetMC(TVirtualMC* mc) { fMC = mc; } + ClassDef(FairModule, 4); }; diff --git a/Detectors/3rdparty/FairRootBase/src/FairMCApplication.cxx b/Detectors/3rdparty/FairRootBase/src/FairMCApplication.cxx index b84f5c9b278c5..f56cdc33f13e9 100644 --- a/Detectors/3rdparty/FairRootBase/src/FairMCApplication.cxx +++ b/Detectors/3rdparty/FairRootBase/src/FairMCApplication.cxx @@ -613,6 +613,20 @@ TVirtualMCApplication* FairMCApplication::CloneForWorker() const return workerApplication; } +//_____________________________________________________________________________ +FairModule* FairMCApplication::OwnModule(FairModule* module) const +{ + // the worker's clone of a module, found by name + if (!module || IsMasterApplication() || fModules->IndexOf(module) >= 0) { + return module; + } + auto clone = static_cast(fModules->FindObject(module->GetName())); + if (!clone) { + LOG(fatal) << "FairMCApplication: no worker clone for module " << module->GetName(); + } + return clone; +} + //_____________________________________________________________________________ void FairMCApplication::InitOnWorker() { @@ -1014,7 +1028,8 @@ void FairMCApplication::InitGeometry() while ((obj = fModIter->Next())) { detector = dynamic_cast(obj); module = dynamic_cast(obj); - if (module) { + if (module && IsMasterApplication()) { + // media cuts are set once on the master; multithreaded Geant4 workers share them module->SetSpecialPhysicsCuts(); } if (detector) { @@ -1101,7 +1116,7 @@ void FairMCApplication::InitGeometry() } FairVolume* fNewV = new FairVolume(fv->GetName(), id); fNewV->setModId(fv->getModId()); - fNewV->SetModule(fv->GetModule()); + fNewV->SetModule(OwnModule(fv->GetModule())); fNewV->setCopyNo(fN->GetNumber()); fNewV->setMCid(id); fVolMap.insert(pair(id, fNewV)); @@ -1109,13 +1124,22 @@ void FairMCApplication::InitGeometry() } else { FairVolume* fNewV = new FairVolume(fv->GetName(), id); fNewV->setModId(fv->getModId()); - fNewV->SetModule(fv->GetModule()); + fNewV->SetModule(OwnModule(fv->GetModule())); fNewV->setCopyNo(1); fNewV->setMCid(id); fVolMap.insert(pair(id, fNewV)); } - } else { + } else if (IsMasterApplication()) { fVolMap.insert(pair(id, fv)); + } else { + // the volume object is shared with the master; a worker needs its own, pointing at its own module + FairVolume* fNewV = new FairVolume(fv->GetName(), id); + fNewV->setModId(fv->getModId()); + fNewV->SetModule(OwnModule(fv->GetModule())); + fNewV->setCopyNo(fv->getCopyNo()); + fNewV->setMCid(id); + fNewV->setGeoNode(fv->getGeoNode()); + fVolMap.insert(pair(id, fNewV)); } } // end off loop Fill sensitive volumes diff --git a/Detectors/Base/CMakeLists.txt b/Detectors/Base/CMakeLists.txt index 2025c6855c4db..e225a2d975ef3 100644 --- a/Detectors/Base/CMakeLists.txt +++ b/Detectors/Base/CMakeLists.txt @@ -35,6 +35,7 @@ o2_add_library(DetectorsBase src/Stack.cxx src/VMCSeederService.cxx src/SimRandom.cxx + src/MTWorkerHooks.cxx src/GlobalParams.cxx src/O2Tessellated.cxx src/TGeoGeometryUtils.cxx @@ -97,6 +98,13 @@ o2_add_test( PUBLIC_LINK_LIBRARIES O2::DetectorsBase LABELS detectorsbase) +o2_add_test( + DetImplClone + SOURCES test/testDetImplClone.cxx + COMPONENT_NAME DetectorsBase + PUBLIC_LINK_LIBRARIES O2::DetectorsBase + LABELS detectorsbase) + o2_add_test( O2Tessellated SOURCES test/testO2Tessellated.cxx diff --git a/Detectors/Base/include/DetectorsBase/Detector.h b/Detectors/Base/include/DetectorsBase/Detector.h index b110112ba2c60..9ed345e6827c5 100644 --- a/Detectors/Base/include/DetectorsBase/Detector.h +++ b/Detectors/Base/include/DetectorsBase/Detector.h @@ -219,6 +219,9 @@ class Detector : public FairDetector // ... for things that should be setup in each simulation worker separately) virtual void initializeLate() = 0; + /// whether this detector may run in multithreaded Geant4, with one clone per worker thread + virtual bool isMTSafe() const { return false; } + /// helper wrapper function to register a geometry volume given by name with FairRoot /// @returns The MonteCarlo ID for the volume int registerSensitiveVolumeAndGetVolID(std::string const& name); @@ -334,6 +337,10 @@ class DetImpl : public o2::base::Detector // offer same constructors as base using Detector::Detector; + /// copy for a worker thread of multithreaded Geant4; the copy creates its own hit buffers in initializeLate + DetImpl(const DetImpl& rhs) : Detector(rhs) {} + DetImpl& operator=(const DetImpl&) = delete; + // default implementation for getHitBranchNames std::string getHitBranchNames(int probe) const override { diff --git a/Detectors/Base/include/DetectorsBase/MTWorkerHooks.h b/Detectors/Base/include/DetectorsBase/MTWorkerHooks.h new file mode 100644 index 0000000000000..36bf1975640e0 --- /dev/null +++ b/Detectors/Base/include/DetectorsBase/MTWorkerHooks.h @@ -0,0 +1,52 @@ +// Copyright 2019-2026 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef ALICEO2_BASE_MTWORKERHOOKS_H +#define ALICEO2_BASE_MTWORKERHOOKS_H + +#include +#include +#include + +class TVirtualMCDecayer; + +namespace o2::data +{ +class PrimaryChunk; +} + +namespace o2::base +{ + +/// What the simulation setup provides to the transport threads of multithreaded Geant4. +/// Installed once on the master thread before transport starts; each function is called on a worker thread. +struct MTWorkerHooks { + /// Returns the next chunk of primaries for the calling thread, or nullptr when there is no more work + std::function()> fetchChunk; + /// Returns the calling thread's channel to the hit merger + std::function dataChannel; + /// Creates the external decayer of the calling thread (empty if there is none) + std::function makeDecayer; + /// Ends the event loop of the calling thread after the current event + std::function stopWorker; + /// Called on a worker thread at the start of each run, after its user actions exist + std::function beginRunOnWorker; + /// Called once on the master thread after initialisation and before the first transport + std::function beforeTransport; + /// Called on a worker thread at the start of each track + std::function beginTrack; + + static MTWorkerHooks& instance(); +}; + +} // namespace o2::base + +#endif diff --git a/Detectors/Base/src/MTWorkerHooks.cxx b/Detectors/Base/src/MTWorkerHooks.cxx new file mode 100644 index 0000000000000..85df356d67b1c --- /dev/null +++ b/Detectors/Base/src/MTWorkerHooks.cxx @@ -0,0 +1,23 @@ +// Copyright 2019-2026 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#include "DetectorsBase/MTWorkerHooks.h" + +namespace o2::base +{ + +MTWorkerHooks& MTWorkerHooks::instance() +{ + static MTWorkerHooks hooks; + return hooks; +} + +} // namespace o2::base diff --git a/Detectors/Base/src/Stack.cxx b/Detectors/Base/src/Stack.cxx index e7caca2982318..5f6b8c6cf791c 100644 --- a/Detectors/Base/src/Stack.cxx +++ b/Detectors/Base/src/Stack.cxx @@ -120,6 +120,7 @@ Stack::Stack(const Stack& rhs) : FairGenericStack(rhs), mStack(), mParticles(), + mTrackIDtoParticlesEntry(rhs.mTrackIDtoParticlesEntry.size(), -1), mTracks(nullptr), mIndexMap(), mIndexOfCurrentTrack(-1), @@ -131,10 +132,15 @@ Stack::Stack(const Stack& rhs) mNumberOfPrimariesPopped(0), mStoreMothers(rhs.mStoreMothers), mStoreSecondaries(rhs.mStoreSecondaries), + mPruneKinematics(rhs.mPruneKinematics), mMinHits(rhs.mMinHits), mEnergyCut(rhs.mEnergyCut), - mTrackRefs(new std::vector), - mIsG4Like(rhs.mIsG4Like) + mCleanupThreshold(rhs.mCleanupThreshold), + mIsG4Like(rhs.mIsG4Like), + mIsExternalMode(rhs.mIsExternalMode), + mDoTrackSeeding(rhs.mDoTrackSeeding), + mTransportPrimary(rhs.mTransportPrimary), + mTrackRefs(new std::vector) { LOG(debug) << "copy constructor called"; mTracks = new std::vector(); diff --git a/Detectors/Base/test/testDetImplClone.cxx b/Detectors/Base/test/testDetImplClone.cxx new file mode 100644 index 0000000000000..ab619f5ae99fe --- /dev/null +++ b/Detectors/Base/test/testDetImplClone.cxx @@ -0,0 +1,98 @@ +// Copyright 2019-2026 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file testDetImplClone.cxx +/// \author Sandro Wenzel (CERN) +/// \brief Checks that a worker clone of a detector gets its own hit buffers + +#define BOOST_TEST_MODULE Test DetImplClone +#define BOOST_TEST_MAIN +#define BOOST_TEST_DYN_LINK +#include + +#include "DetectorsBase/Detector.h" +#include "SimulationDataFormat/BaseHits.h" +#include + +namespace o2::test +{ +using TestHit = o2::BasicXYZEHit; + +// A minimal detector with one hit vector in shared-memory style buffers +class CloneTestDetector : public o2::base::DetImpl +{ + public: + CloneTestDetector() : o2::base::DetImpl("ITS", true) {} + CloneTestDetector(const CloneTestDetector& rhs) = default; + + Bool_t ProcessHits(FairVolume*) override { return false; } + void Register() override {} + void Reset() override {} + void InitializeO2Detector() override {} + + std::vector* getHits(int i) const { return i == 0 ? mHits : nullptr; } + bool initialized() const { return mInitialized; } + void* cachedBuffer(int b) const { return mCachedPtr[b].empty() ? nullptr : mCachedPtr[b][0]; } + const void* busyFlag(int b) const { return mShmBusy[b]; } + + std::vector* mHits = nullptr; +}; +} // namespace o2::test + +namespace o2::base +{ +template <> +struct UseShm { + static constexpr bool value = true; +}; +} // namespace o2::base + +// the test detector has no ROOT dictionary, so the ClassDef members of its base are defined here +template <> +TClass* o2::base::DetImpl::Class() +{ + return nullptr; +} +template <> +void o2::base::DetImpl::Streamer(TBuffer&) +{ +} + +BOOST_AUTO_TEST_CASE(clone_gets_own_hit_buffers) +{ + o2::test::CloneTestDetector master; + master.initializeLate(); + master.BeginEvent(); + BOOST_REQUIRE(master.getHits(0) != nullptr); + + std::unique_ptr clone(static_cast(master.CloneModule())); + BOOST_CHECK(!clone->initialized()); + for (int b = 0; b < 3; ++b) { + BOOST_CHECK(clone->cachedBuffer(b) == nullptr); + BOOST_CHECK(clone->busyFlag(b) == nullptr); + } + + clone->initializeLate(); + clone->BeginEvent(); + BOOST_REQUIRE(clone->getHits(0) != nullptr); + BOOST_CHECK(clone->getHits(0) != master.getHits(0)); + for (int b = 0; b < 3; ++b) { + BOOST_CHECK(clone->cachedBuffer(b) != master.cachedBuffer(b)); + BOOST_CHECK(clone->busyFlag(b) != master.busyFlag(b)); + } + + // hits of the clone do not reach the master, and the master survives the clone + clone->getHits(0)->emplace_back(); + BOOST_CHECK(master.getHits(0)->empty()); + clone.reset(); + master.getHits(0)->emplace_back(); + BOOST_CHECK_EQUAL(master.getHits(0)->size(), 1); +} diff --git a/Detectors/CPV/simulation/include/CPVSimulation/Detector.h b/Detectors/CPV/simulation/include/CPVSimulation/Detector.h index 3bd41f89c7ea2..3595d74df2d6d 100644 --- a/Detectors/CPV/simulation/include/CPVSimulation/Detector.h +++ b/Detectors/CPV/simulation/include/CPVSimulation/Detector.h @@ -70,6 +70,7 @@ class Detector : public o2::base::DetImpl /// Initializing detector /// void InitializeO2Detector() final; + bool isMTSafe() const override { return true; } /// /// Processing hit creation in the CPV crystalls diff --git a/Detectors/EMCAL/simulation/include/EMCALSimulation/Detector.h b/Detectors/EMCAL/simulation/include/EMCALSimulation/Detector.h index 2edee5aabe153..72efbbb846d05 100644 --- a/Detectors/EMCAL/simulation/include/EMCALSimulation/Detector.h +++ b/Detectors/EMCAL/simulation/include/EMCALSimulation/Detector.h @@ -71,6 +71,7 @@ class Detector : public o2::base::DetImpl /// \brief Initializing detector void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } /// \brief Processing hit creation in the EMCAL scintillator volume /// \param v Current sensitive volume diff --git a/Detectors/EMCAL/simulation/src/Detector.cxx b/Detectors/EMCAL/simulation/src/Detector.cxx index 4e6c7ced5c3e8..f2f29eb2cc0be 100644 --- a/Detectors/EMCAL/simulation/src/Detector.cxx +++ b/Detectors/EMCAL/simulation/src/Detector.cxx @@ -84,9 +84,9 @@ Detector::Detector(const Detector& rhs) mBirkC0(rhs.mBirkC0), mBirkC1(rhs.mBirkC1), mBirkC2(rhs.mBirkC2), - mSensitive(), + mSensitive(rhs.mSensitive), mSMVolumeID(), - mSMVolNames(), + mSMVolNames(rhs.mSMVolNames), mVolumeIDScintillator(-1), mHits(o2::utils::createSimVector()), mGeometry(rhs.mGeometry), diff --git a/Detectors/FIT/FDD/simulation/include/FDDSimulation/Detector.h b/Detectors/FIT/FDD/simulation/include/FDDSimulation/Detector.h index f25a7a18cfb4c..9647aace614b9 100644 --- a/Detectors/FIT/FDD/simulation/include/FDDSimulation/Detector.h +++ b/Detectors/FIT/FDD/simulation/include/FDDSimulation/Detector.h @@ -44,6 +44,7 @@ class Detector : public o2::base::DetImpl ~Detector() override; void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } Bool_t ProcessHits(FairVolume* v = nullptr) override; diff --git a/Detectors/FIT/FT0/simulation/include/FT0Simulation/Detector.h b/Detectors/FIT/FT0/simulation/include/FT0Simulation/Detector.h index 75f3c7a1f188f..3ec22d16c0296 100644 --- a/Detectors/FIT/FT0/simulation/include/FT0Simulation/Detector.h +++ b/Detectors/FIT/FT0/simulation/include/FT0Simulation/Detector.h @@ -71,6 +71,7 @@ class Detector : public o2::base::DetImpl /// Initialization of the detector is done here void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } /// This method is called for each step during simulation (see FairMCApplication::Stepping()) Bool_t ProcessHits(FairVolume* v) override; diff --git a/Detectors/FIT/FT0/simulation/src/Detector.cxx b/Detectors/FIT/FT0/simulation/src/Detector.cxx index de7e72e5f1709..813f38c41d712 100644 --- a/Detectors/FIT/FT0/simulation/src/Detector.cxx +++ b/Detectors/FIT/FT0/simulation/src/Detector.cxx @@ -38,6 +38,7 @@ #include "FT0Simulation/Detector.h" #include "DetectorsBase/Stack.h" #include "DetectorsBase/SimRandom.h" +#include using namespace o2::ft0; using o2::ft0::Geometry; @@ -55,6 +56,8 @@ Detector::Detector(Bool_t Active) Detector::Detector(const Detector& rhs) : o2::base::DetImpl(rhs), mIdSens1(rhs.mIdSens1), mPMTeff(rhs.mPMTeff), mHits(o2::utils::createSimVector()), mTrackIdTop(-1), mTrackIdMCPtop(-1) { + // the channel lookup is built with the optical geometry, which only the master constructs + std::copy(std::begin(rhs.mSim2LUT), std::end(rhs.mSim2LUT), mSim2LUT); } Detector::~Detector() diff --git a/Detectors/FIT/FV0/simulation/include/FV0Simulation/Detector.h b/Detectors/FIT/FV0/simulation/include/FV0Simulation/Detector.h index fdff96b6eafa8..e6301733f181e 100644 --- a/Detectors/FIT/FV0/simulation/include/FV0Simulation/Detector.h +++ b/Detectors/FIT/FV0/simulation/include/FV0Simulation/Detector.h @@ -50,6 +50,7 @@ class Detector : public o2::base::DetImpl /// Initializes the detector (adds sensitive volume) void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } /// This method is called for each step during simulation (see FairMCApplication::Stepping()) Bool_t ProcessHits(FairVolume* v = nullptr) override; diff --git a/Detectors/FIT/FV0/simulation/src/Detector.cxx b/Detectors/FIT/FV0/simulation/src/Detector.cxx index a91c599f30020..92a3158623611 100644 --- a/Detectors/FIT/FV0/simulation/src/Detector.cxx +++ b/Detectors/FIT/FV0/simulation/src/Detector.cxx @@ -45,9 +45,6 @@ Detector::~Detector() if (mHits) { o2::utils::freeSimVector(mHits); // delete mHits; } - if (mGeometry) { - delete mGeometry; - } } Detector::Detector(Bool_t isActive) diff --git a/Detectors/FOCAL/simulation/include/FOCALSimulation/Detector.h b/Detectors/FOCAL/simulation/include/FOCALSimulation/Detector.h index 5f6bed3a037b0..b890f75465cd6 100644 --- a/Detectors/FOCAL/simulation/include/FOCALSimulation/Detector.h +++ b/Detectors/FOCAL/simulation/include/FOCALSimulation/Detector.h @@ -81,6 +81,7 @@ class Detector : public o2::base::DetImpl /// \brief Initializing detector void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } /// \brief Processing hit creation in the FOCAL sensitive volume /// \param v Current sensitive volume diff --git a/Detectors/FOCAL/simulation/src/Detector.cxx b/Detectors/FOCAL/simulation/src/Detector.cxx index 9c89baa8f271b..21edc4564b89f 100644 --- a/Detectors/FOCAL/simulation/src/Detector.cxx +++ b/Detectors/FOCAL/simulation/src/Detector.cxx @@ -49,7 +49,12 @@ Detector::Detector(bool active, std::string geofilename) } Detector::Detector(const Detector& rhs) - : o2::base::DetImpl(rhs) + : o2::base::DetImpl(rhs), + mHits(o2::utils::createSimVector()), + mCurrentSuperparent(nullptr), + mCurrentTrack(-1), + mCurrentPrimaryID(-1), + mCurrentParentID(-1) { mGeometry = rhs.mGeometry; mMedSensHCal = rhs.mMedSensHCal; diff --git a/Detectors/HMPID/simulation/include/HMPIDSimulation/Detector.h b/Detectors/HMPID/simulation/include/HMPIDSimulation/Detector.h index 04280748de6e3..6aaf372e4f09c 100644 --- a/Detectors/HMPID/simulation/include/HMPIDSimulation/Detector.h +++ b/Detectors/HMPID/simulation/include/HMPIDSimulation/Detector.h @@ -38,6 +38,7 @@ class Detector : public o2::base::DetImpl } void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } bool ProcessHits(FairVolume* v) override; o2::hmpid::HitType* AddHit(float x, float y, float z, float time, float energy, Int_t trackId, Int_t detId); void GenFee(float qtot); diff --git a/Detectors/HMPID/simulation/src/Detector.cxx b/Detectors/HMPID/simulation/src/Detector.cxx index 8780a93afb56c..539ab28568561 100644 --- a/Detectors/HMPID/simulation/src/Detector.cxx +++ b/Detectors/HMPID/simulation/src/Detector.cxx @@ -42,7 +42,8 @@ namespace hmpid Detector::Detector(Bool_t active) : o2::base::DetImpl("HMP", active), mHits(new std::vector) {} -Detector::Detector(const Detector& other) : mSensitiveVolumes(other.mSensitiveVolumes), +Detector::Detector(const Detector& other) : o2::base::DetImpl(other), + mSensitiveVolumes(other.mSensitiveVolumes), mHits(new std::vector) {} void Detector::InitializeO2Detector() diff --git a/Detectors/ITSMFT/ITS/simulation/include/ITSSimulation/Detector.h b/Detectors/ITSMFT/ITS/simulation/include/ITSSimulation/Detector.h index 57301ac4babd0..48460b6e7c389 100644 --- a/Detectors/ITSMFT/ITS/simulation/include/ITSSimulation/Detector.h +++ b/Detectors/ITSMFT/ITS/simulation/include/ITSSimulation/Detector.h @@ -229,6 +229,8 @@ class Detector : public o2::base::DetImpl /// Returns the number of layers Int_t getNumberOfLayers() const { return mNumberLayers; } + bool isMTSafe() const override { return true; } + GeometryTGeo* mGeometryTGeo; //! access to geometry details protected: diff --git a/Detectors/ITSMFT/MFT/simulation/include/MFTSimulation/Detector.h b/Detectors/ITSMFT/MFT/simulation/include/MFTSimulation/Detector.h index 3c44705f747c0..fd5bf8b5b0a88 100644 --- a/Detectors/ITSMFT/MFT/simulation/include/MFTSimulation/Detector.h +++ b/Detectors/ITSMFT/MFT/simulation/include/MFTSimulation/Detector.h @@ -62,6 +62,7 @@ class Detector : public o2::base::DetImpl /// Initialization of the detector is done here void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } /// This method is called for each step during simulation (see FairMCApplication::Stepping()) Bool_t ProcessHits(FairVolume* v = nullptr) override; diff --git a/Detectors/MUON/MCH/Simulation/include/MCHSimulation/Detector.h b/Detectors/MUON/MCH/Simulation/include/MCHSimulation/Detector.h index 2c378e0a8e805..d1a04f2cf18a5 100644 --- a/Detectors/MUON/MCH/Simulation/include/MCHSimulation/Detector.h +++ b/Detectors/MUON/MCH/Simulation/include/MCHSimulation/Detector.h @@ -34,6 +34,7 @@ class Detector : public o2::base::DetImpl Bool_t ProcessHits(FairVolume* v = nullptr) override; void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } void Register() override; diff --git a/Detectors/MUON/MID/Simulation/include/MIDSimulation/Detector.h b/Detectors/MUON/MID/Simulation/include/MIDSimulation/Detector.h index f44093fce0bfa..39d76d3c4e91e 100644 --- a/Detectors/MUON/MID/Simulation/include/MIDSimulation/Detector.h +++ b/Detectors/MUON/MID/Simulation/include/MIDSimulation/Detector.h @@ -29,6 +29,7 @@ class Detector : public o2::base::DetImpl Detector(bool active = true); void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } bool ProcessHits(FairVolume* vol = nullptr) override; diff --git a/Detectors/PHOS/simulation/include/PHOSSimulation/Detector.h b/Detectors/PHOS/simulation/include/PHOSSimulation/Detector.h index cba3a1a7e9c6a..2d177c59afb32 100644 --- a/Detectors/PHOS/simulation/include/PHOSSimulation/Detector.h +++ b/Detectors/PHOS/simulation/include/PHOSSimulation/Detector.h @@ -77,6 +77,7 @@ class Detector : public o2::base::DetImpl /// Initializing detector /// void InitializeO2Detector() final; + bool isMTSafe() const override { return true; } /// /// Processing hit creation in the PHOS crystalls diff --git a/Detectors/TOF/simulation/include/TOFSimulation/Detector.h b/Detectors/TOF/simulation/include/TOFSimulation/Detector.h index 1dfd730efe716..04652a41c0217 100644 --- a/Detectors/TOF/simulation/include/TOFSimulation/Detector.h +++ b/Detectors/TOF/simulation/include/TOFSimulation/Detector.h @@ -82,6 +82,7 @@ class Detector : public o2::base::DetImpl ~Detector() override; void InitializeO2Detector() final; + bool isMTSafe() const override { return true; } Bool_t ProcessHits(FairVolume* v = nullptr) final; diff --git a/Detectors/TPC/simulation/include/TPCSimulation/Detector.h b/Detectors/TPC/simulation/include/TPCSimulation/Detector.h index 1f5c2840ff9f5..d937f2054383c 100644 --- a/Detectors/TPC/simulation/include/TPCSimulation/Detector.h +++ b/Detectors/TPC/simulation/include/TPCSimulation/Detector.h @@ -82,6 +82,7 @@ class Detector : public o2::base::DetImpl /** Initialization of the detector is done here */ void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } /** this method is called for each step during simulation * (see FairMCApplication::Stepping()) diff --git a/Detectors/TRD/simulation/include/TRDSimulation/Detector.h b/Detectors/TRD/simulation/include/TRDSimulation/Detector.h index bde79532f387e..9a10af8462e1f 100644 --- a/Detectors/TRD/simulation/include/TRDSimulation/Detector.h +++ b/Detectors/TRD/simulation/include/TRDSimulation/Detector.h @@ -35,6 +35,7 @@ class Detector : public o2::base::DetImpl Detector(Bool_t active = true); ~Detector() override; void InitializeO2Detector() override; + bool isMTSafe() const override { return true; } bool ProcessHits(FairVolume* v = nullptr) override; void Register() override; std::vector* getHits(int iColl) const diff --git a/Detectors/ZDC/simulation/include/ZDCSimulation/Detector.h b/Detectors/ZDC/simulation/include/ZDCSimulation/Detector.h index 11fb95001b981..451aad05cc64d 100644 --- a/Detectors/ZDC/simulation/include/ZDCSimulation/Detector.h +++ b/Detectors/ZDC/simulation/include/ZDCSimulation/Detector.h @@ -82,6 +82,7 @@ class Detector : public o2::base::DetImpl #endif void InitializeO2Detector() final; + bool isMTSafe() const override { return true; } Bool_t ProcessHits(FairVolume* v = nullptr) final; @@ -140,6 +141,9 @@ class Detector : public o2::base::DetImpl void resetHitIndices(); + /// Builds the fast simulation models, if they are enabled + void setupFastSim(); + // common function for hit creation (can be called from multiple interfaces) bool createOrAddHit(int detector, int sector, diff --git a/Detectors/ZDC/simulation/src/Detector.cxx b/Detectors/ZDC/simulation/src/Detector.cxx index 198044b34118b..39ef4155389b6 100644 --- a/Detectors/ZDC/simulation/src/Detector.cxx +++ b/Detectors/ZDC/simulation/src/Detector.cxx @@ -69,7 +69,31 @@ Detector::Detector(Bool_t active) mMediumPMCid = -1; // minus for unitialized mMediumPMQid = -2; // different to PMC in any case resetHitIndices(); + setupFastSim(); +} +//_____________________________________________________________________________ +Detector::Detector(const Detector& rhs) + : o2::base::DetImpl(rhs), + mXImpact(-999, -999, -999), + // the medium ids are set while the geometry is built, which only the master does + mMediumPMCid(rhs.mMediumPMCid), + mMediumPMQid(rhs.mMediumPMQid), + mHits(new std::vector), + mNeutronResponseImage(rhs.mNeutronResponseImage), + mProtonResponseImage(rhs.mProtonResponseImage) +{ + mTrackEta = 999; + mTotLightPMC = 0; + mTotLightPMQ = 0; + resetHitIndices(); + // each copy has its own fast simulation models, since running a model changes its state + setupFastSim(); +} + +//_____________________________________________________________________________ +void Detector::setupFastSim() +{ #ifdef ZDC_FASTSIM_ONNX // If FastSim module was disabled, log appropriate message // otherwise check if all necessary parameters were passed, if so try build objects @@ -127,13 +151,6 @@ Detector::Detector(Bool_t active) #endif } -//_____________________________________________________________________________ -Detector::Detector(const Detector& rhs) - : o2::base::DetImpl(rhs), - mHits(new std::vector) -{ -} - //_____________________________________________________________________________ #ifdef ZDC_FASTSIM_ONNX Detector::~Detector() diff --git a/Detectors/gconfig/CMakeLists.txt b/Detectors/gconfig/CMakeLists.txt index d35bbed3deb36..afdbab663b744 100644 --- a/Detectors/gconfig/CMakeLists.txt +++ b/Detectors/gconfig/CMakeLists.txt @@ -27,7 +27,7 @@ if(TGeo2VecGeom_FOUND) endif() endif() -set(G4SETUP_SOURCES src/G4Config.cxx src/G4RunConfiguration.cxx src/G4LocalFieldConstruction.cxx +set(G4SETUP_SOURCES src/G4Config.cxx src/G4RunConfiguration.cxx src/G4LocalFieldConstruction.cxx src/G4FieldCache.cxx src/VecGeomNavigation.cxx src/O2MonopolePhysics.cxx) if(G4SETUP_WITH_VECGEOM) list(APPEND G4SETUP_SOURCES src/VecGeomG4Map.cxx src/VecGeomChecks.cxx src/VecGeomG4NavigatorBase.cxx diff --git a/Detectors/gconfig/commonConfig.C b/Detectors/gconfig/commonConfig.C index e218da39b8461..99e495b182dd7 100644 --- a/Detectors/gconfig/commonConfig.C +++ b/Detectors/gconfig/commonConfig.C @@ -45,9 +45,14 @@ void stackSetup(T* vmc, R* run) */ } -void decayerSetup(TVirtualMC* vmc) +TVirtualMCDecayer* makeDecayer() { auto decayer = new o2::eventgen::DecayerPythia8; decayer->Init(); - vmc->SetExternalDecayer(decayer); + return decayer; +} + +void decayerSetup(TVirtualMC* vmc) +{ + vmc->SetExternalDecayer(makeDecayer()); } diff --git a/Detectors/gconfig/g4Config.C b/Detectors/gconfig/g4Config.C index f57a364f908d8..6173066309d98 100644 --- a/Detectors/gconfig/g4Config.C +++ b/Detectors/gconfig/g4Config.C @@ -68,12 +68,22 @@ R__LOAD_LIBRARY(libgeant4vmc) #include #include "SimSetup/G4RunConfiguration.h" #include "SimSetup/VecGeomNavigation.h" +#include "SimSetup/G4FieldCache.h" #include "DetectorsBase/VMCSeederService.h" #include "CLHEP/Random/Random.h" #include "TRandom.h" +#include "DetectorsBase/MTWorkerHooks.h" +#include "G4RunManager.hh" +#include "G4Threading.hh" +#include "G4UImanager.hh" +#include "G4RegionStore.hh" +#include "G4ProductionCuts.hh" #endif #include "commonConfig.C" +// skips particles unknown to Geant4 (e.g. some ions) with a warning instead of a crash +constexpr const char* kSkipUnknownParticles = "/mcPrimaryGenerator/skipUnknownParticles true"; + // Configuration macro for Geant4 VirtualMC void Config() { @@ -108,6 +118,11 @@ void Config() //Geant4 VMC 3.x Bool_t mtMode = FairRunSim::Instance()->IsMT(); +#ifndef G4MULTITHREADED + if (mtMode) { + LOG(fatal) << "Multithreaded mode (--isMT) needs Geant4 built with GEANT4_BUILD_MULTITHREADED=ON"; + } +#endif Bool_t specialStacking = true; // leads to default stack behaviour in which new primaries are only started if the previous // one and all of its secondaries have been transported // any other choice is dangerously inconsistent with the FinishPrimary() interface of VMCApp @@ -152,6 +167,11 @@ void Config() } /// avoid the use of G4BACKTRACE (it seems to inferfere with process logic in o2-sim) setenv("G4BACKTRACE", "none", 1); + if (mtMode) { + // one transport thread per worker, read by the TGeant4 constructor + const char* nthreads = getenv("ALICE_NSIMWORKERS"); + setenv("G4FORCENUMBEROFTHREADS", nthreads ? nthreads : "1", 1); + } /// Create the G4 VMC TGeant4* geant4 = new TGeant4("TGeant4", "The Geant4 Monte Carlo", runConfiguration); std::cout << "Geant4 has been created." << std::endl; @@ -165,6 +185,33 @@ void Config() CLHEP::HepRandom::setTheSeeds(seeds); }); + if (mtMode) { + if (g4Params.navmode == o2::conf::EG4Nav::kVecGeom) { + LOG(fatal) << "G4.navmode=kVecGeom is not yet supported in multithreaded mode (--isMT)"; + } + // each thread takes one sub-event at a time + geant4->ProcessGeantCommand("/run/eventModulo 1"); + // a thread without more work ends its own event loop after the current (empty) event + o2::base::MTWorkerHooks::instance().stopWorker = []() { G4RunManager::GetRunManager()->AbortRun(true); }; + // each worker thread creates its own external decayer + o2::base::MTWorkerHooks::instance().makeDecayer = &makeDecayer; + // the field seen by a track must not depend on the previous track of the thread + o2::base::MTWorkerHooks::instance().beginTrack = &o2::simsetup::invalidateFieldCache; + // the primary generator exists only on the worker threads, so its commands are applied there + o2::base::MTWorkerHooks::instance().beginRunOnWorker = []() { G4UImanager::GetUIpointer()->ApplyCommand(kSkipUnknownParticles); }; + // makes the first run apply the media energy cuts + o2::base::MTWorkerHooks::instance().beforeTransport = []() { + for (auto region : *G4RegionStore::GetInstance()) { + if (auto cuts = region->GetProductionCuts()) { + for (int i = 0; i < NumberOfG4CutIndex; ++i) { + cuts->SetProductionCut(cuts->GetProductionCut(i), i); + } + } + } + G4UImanager::GetUIpointer()->ApplyCommand("/run/physicsModified"); + }; + } + if (g4Params.navmode == o2::conf::EG4Nav::kVecGeom) { o2::simsetup::installVecGeomNavigator(); } @@ -187,6 +234,9 @@ void Config() /// (verbose level, global range cut, ..) geant4->ProcessGeantMacro(configm1.Data()); } + if (!mtMode) { + geant4->ProcessGeantCommand(kSkipUnknownParticles); + } // Enter in Geant4 Interactive mode // geant4->StartGeantUI(); diff --git a/Detectors/gconfig/g4config.in b/Detectors/gconfig/g4config.in index cd61afd493560..cfc10ee20918b 100644 --- a/Detectors/gconfig/g4config.in +++ b/Detectors/gconfig/g4config.in @@ -33,7 +33,6 @@ /mcPhysics/emModel/setEmModel PAI /mcPhysics/emModel/setRegions TRD_Gas-mix /mcPhysics/emModel/setParticles e- e+ mu- mu+ pi- pi+ kaon- kaon+ proton anti_proton deuteron triton He3 alpha -/mcPrimaryGenerator/skipUnknownParticles true # don't crash when seeing unknown ion etc. (issue warning) # # Precise Msc for EMCAL diff --git a/Detectors/gconfig/include/SimSetup/G4FieldCache.h b/Detectors/gconfig/include/SimSetup/G4FieldCache.h new file mode 100644 index 0000000000000..694ce07f9b5f8 --- /dev/null +++ b/Detectors/gconfig/include/SimSetup/G4FieldCache.h @@ -0,0 +1,23 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef O2_SIMSETUP_G4FIELDCACHE_H_ +#define O2_SIMSETUP_G4FIELDCACHE_H_ + +namespace o2::simsetup +{ + +/// Makes the cached magnetic fields of the calling Geant4 thread forget their last value +void invalidateFieldCache(); + +} // namespace o2::simsetup + +#endif diff --git a/Detectors/gconfig/src/G4FieldCache.cxx b/Detectors/gconfig/src/G4FieldCache.cxx new file mode 100644 index 0000000000000..925b6b461f580 --- /dev/null +++ b/Detectors/gconfig/src/G4FieldCache.cxx @@ -0,0 +1,34 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#include "SimSetup/G4FieldCache.h" + +#include +#include +#include + +namespace o2::simsetup +{ + +void invalidateFieldCache() +{ + // the field caches reuse their value only near the last point, so one query far away invalidates them + constexpr G4double farAway[4] = {1.e12, 1.e12, 1.e12, 0.}; + G4double value[6]; + // the field manager store holds the field managers of the calling thread + for (auto manager : *G4FieldManagerStore::GetInstance()) { + if (auto field = manager->GetDetectorField()) { + field->GetFieldValue(farAway, value); + } + } +} + +} // namespace o2::simsetup diff --git a/Steer/include/Steer/O2MCApplication.h b/Steer/include/Steer/O2MCApplication.h index e43a61ec419b5..575d0082d0f66 100644 --- a/Steer/include/Steer/O2MCApplication.h +++ b/Steer/include/Steer/O2MCApplication.h @@ -24,6 +24,8 @@ #include #include +#include +#include namespace o2 { @@ -39,6 +41,16 @@ class O2MCApplication : public O2MCApplicationBase using O2MCApplicationBase::O2MCApplicationBase; ~O2MCApplication() override = default; + // multithreaded Geant4: each worker thread gets its own copy of the application + O2MCApplication(const O2MCApplication& rhs); + TVirtualMCApplication* CloneForWorker() const override; + void InitOnWorker() override; + void BeginRunOnWorker() override; + void FinishRunOnWorker() override; + + /// a worker thread asks the primary server for its sub-event here, as a worker process does before ProcessEvent + void BeginEvent() override; + // triggers data sending/io void SendData(); @@ -47,6 +59,14 @@ class O2MCApplication : public O2MCApplicationBase /** Define actions at the end of event */ void FinishEvent() override { + if (mIsEndOfWork) { + // the empty event that ends a worker thread's loop: nothing to send + for (auto det : listActiveDetectors) { + det->EndOfEvent(); + } + fStack->Reset(); + return; + } // update the stack fStack->FillTrackArray(); fStack->UpdateTrackIndex(fActiveDetectors); @@ -89,16 +109,31 @@ class O2MCApplication : public O2MCApplicationBase { mPrimaries = p; } + void setPrimaries(std::vector&& p) + { + mPrimaries = std::move(p); + } + + /// whether a sub-event passes the O2SIM_RESTRICT_EVENTPART ("event:part") selection, a debugging aid + static bool isSelectedSubEvent(int eventID, int part); void setSimDataChannel(fair::mq::Channel* channel) { mSimDataChannel = channel; } void setSubEventInfo(o2::data::SubEventInfo* i); std::vector mPrimaries; //! - fair::mq::Channel* mSimDataChannel; //! generic channel on which to send sim data + fair::mq::Channel* mSimDataChannel = nullptr; //! generic channel on which to send sim data o2::data::SubEventInfo* mSubEventInfo = nullptr; //! what are we currently processing? std::vector mActiveO2Detectors; //! active (data taking) o2 detectors + private: + void fetchWorkerChunk(); + bool mIsMTWorker = false; //! this is the application of a multithreaded Geant4 worker + bool mIsEndOfWork = false; //! the current event is the empty one ending the worker's loop + bool mWorkerInitialized = false; //! this worker's detectors and data channel are set up + o2::data::SubEventInfo mWorkerSubEventInfo; //! sub-event being transported by this worker + std::unique_ptr mWorkerRandom; //! this worker's transport random generator + ClassDefOverride(O2MCApplication, 1); //Interface to MonteCarlo application }; diff --git a/Steer/include/Steer/O2MCApplicationBase.h b/Steer/include/Steer/O2MCApplicationBase.h index cdbf99b3b4b01..42c39c1275777 100644 --- a/Steer/include/Steer/O2MCApplicationBase.h +++ b/Steer/include/Steer/O2MCApplicationBase.h @@ -46,6 +46,17 @@ class O2MCApplicationBase : public FairMCApplication initTrackRefHook(); } + /// copy for a worker thread of multithreaded Geant4; the macro hooks are copied, not compiled again + O2MCApplicationBase(const O2MCApplicationBase& rhs) + : FairMCApplication(rhs), + mCutParams(rhs.mCutParams), + mModIdToName(rhs.mModIdToName), + mTrackRefFcn(rhs.mTrackRefFcn), + mKeepStepFcn(rhs.mKeepStepFcn), + mHasStepFilterMacro(rhs.mHasStepFilterMacro) + { + } + ~O2MCApplicationBase() override = default; void Stepping() override; diff --git a/Steer/src/O2MCApplication.cxx b/Steer/src/O2MCApplication.cxx index fc45a870603b4..54f1760f6af96 100644 --- a/Steer/src/O2MCApplication.cxx +++ b/Steer/src/O2MCApplication.cxx @@ -10,6 +10,9 @@ // or submit itself to any jurisdiction. #include +#include +#include +#include #include #include @@ -36,6 +39,8 @@ #include #include "SimConfig/GlobalProcessCutSimParam.h" #include +#include +#include #include // full type: FairField derives from TVirtualMagField #include #include "DetectorsBase/GeometryManagerParam.h" @@ -51,6 +56,9 @@ #include #include #include "DetectorsBase/SimRandom.h" +#include "DetectorsBase/MTWorkerHooks.h" +#include +#include namespace o2 { @@ -185,6 +193,9 @@ bool O2MCApplicationBase::seedsInPreTrack() const void O2MCApplicationBase::PreTrack() { + if (auto& beginTrack = o2::base::MTWorkerHooks::instance().beginTrack) { + beginTrack(); + } if (mCutParams.trackSeed && seedsInPreTrack()) { // Per-track seeding for engines that do not go through // o2::data::Stack::PopNextTrack(). Geant4 is one: it takes primaries via @@ -275,18 +286,24 @@ void O2MCApplicationBase::InitGeometry() // set up the engine seeding here, on the main thread, so that transport threads never JIT-compile it o2::base::VMCSeederService::instance(); // load special cuts which might be given from the outside first. + // the media cuts are set up once, on the master; the worker threads of multithreaded Geant4 share them auto& matMgr = o2::base::MaterialManager::Instance(); - matMgr.loadCutsAndProcessesFromJSON(o2::base::MaterialManager::ESpecial::kTRUE); - matMgr.SetLowEnergyNeutronTransport(mCutParams.lowneut); + if (IsMasterApplication()) { + matMgr.loadCutsAndProcessesFromJSON(o2::base::MaterialManager::ESpecial::kTRUE); + matMgr.SetLowEnergyNeutronTransport(mCutParams.lowneut); + } // During the following, FairModule::SetSpecialPhysicsCuts will be called for each module FairMCApplication::InitGeometry(); - matMgr.writeCutsAndProcessesToJSON(); // now the sensitive volumes are set up in fVolMap and we can query them for (auto e : fVolMap) { // since fVolMap contains multiple entries (if multiple copies), this may // write to the same entry multiple times mSensitiveVolumes[e.first] = e.second->GetName(); } + if (!IsMasterApplication()) { + return; + } + matMgr.writeCutsAndProcessesToJSON(); std::ofstream sensvolfile("MCStepLoggerSenVol.dat"); for (auto e : mSensitiveVolumes) { sensvolfile << e.first << ":" << e.second << "\n"; @@ -1822,7 +1839,20 @@ void O2MCApplicationBase::initTrackRefHook() void O2MCApplication::initLate() { - o2::utils::ShmManager::Instance().occupySegment(); + const bool isMT = fRun && fRun->IsMT(); + if (isMT) { + // detectors are enabled for multithreaded Geant4 one by one; other modules make no hits + for (auto det : listActiveDetectors) { + auto o2det = dynamic_cast(det); + if (o2det && !o2det->isMTSafe()) { + LOG(fatal) << "Detector " << det->GetName() << " is not yet supported in multithreaded mode (--isMT)"; + } + } + } + o2::utils::ShmManager::Instance().occupySegment(isMT); + if (isMT) { + return; // the hit buffers belong to the detector clones of the worker threads + } for (auto det : listActiveDetectors) { if (dynamic_cast(det)) { ((o2::base::Detector*)det)->initializeLate(); @@ -1882,5 +1912,150 @@ void O2MCApplication::SendData() LOG(info) << "sending message with " << simdataparts.Size() << " parts"; mSimDataChannel->Send(simdataparts); } + +namespace +{ +// serialises the set-up of the worker threads of multithreaded Geant4 +std::mutex gWorkerSetupMutex; + +// the sink of a worker thread in multithreaded mode: its output goes to the hit merger with every event +class WorkerSink : public FairSink +{ + public: + Bool_t InitSink() override { return kTRUE; } + void Close() override {} + void Reset() override {} + Sink_Type GetSinkType() override { return kONLINESINK; } + void SetOutTree(TTree*) override {} + void Fill() override {} + Int_t Write(const char*, Int_t, Int_t) override { return 0; } + void RegisterImpl(const char*, const char*, void*) override {} + void RegisterAny(const char*, const std::type_info&, const std::type_info&, void*) override {} + void WriteFolder() override {} + bool CreatePersistentBranchesAny() override { return true; } + void WriteObject(TObject*, const char*, Int_t) override {} + void WriteGeometry() override {} + FairSink* CloneSink() override { return new WorkerSink(); } +}; +} // namespace + +O2MCApplication::O2MCApplication(const O2MCApplication& rhs) : O2MCApplicationBase(rhs) +{ +} + +TVirtualMCApplication* O2MCApplication::CloneForWorker() const +{ + std::lock_guard lock(gWorkerSetupMutex); + // the run object of the worker thread (FairRunSim::Instance() is thread-local), without generator + auto workerRun = new FairRunSim(kFALSE); + workerRun->SetName(fRun->GetName()); + workerRun->SetIsMT(kTRUE); + workerRun->SetSink(new WorkerSink()); + return new O2MCApplication(*this); +} + +void O2MCApplication::InitOnWorker() +{ + std::lock_guard lock(gWorkerSetupMutex); + mIsMTWorker = true; + fRootManager = FairRootManager::Instance(); + fRootManager->InitSink(); + fRun = FairRunSim::Instance(); + fMC = TVirtualMC::GetMC(); + fMC->SetStack(fStack); + fMC->SetMagField(fxField); + // the modules were cloned before this thread's MC existed + for (auto obj : *fModules) { + if (auto module = dynamic_cast(obj)) { + module->SetMC(fMC); + } + } + // the stack's output branches in this thread's FairRootManager + fStack->Register(); + // this thread's own transport random generator + mWorkerRandom = std::make_unique(); + o2::base::setSimRandom(mWorkerRandom.get()); + // this thread's own external decayer; set here, before Geant4 builds the physics of this thread + auto& hooks = o2::base::MTWorkerHooks::instance(); + if (hooks.makeDecayer) { + fMC->SetExternalDecayer(hooks.makeDecayer()); + } +} + +void O2MCApplication::BeginRunOnWorker() +{ + if (auto& hook = o2::base::MTWorkerHooks::instance().beginRunOnWorker) { + hook(); + } +} + +void O2MCApplication::FinishRunOnWorker() +{ + // nothing to finish on a worker: its output went to the hit merger with every event +} + +void O2MCApplication::BeginEvent() +{ + if (mIsMTWorker) { + fetchWorkerChunk(); + } + O2MCApplicationBase::BeginEvent(); +} + +bool O2MCApplication::isSelectedSubEvent(int eventID, int part) +{ + static const auto selection = []() -> std::optional> { + auto env = getenv("O2SIM_RESTRICT_EVENTPART"); + if (!env) { + return std::nullopt; + } + const std::string value(env); + const auto pos = value.find(':'); + if (pos == std::string::npos) { + return std::make_pair(0, 0); + } + return std::make_pair(std::atoi(value.substr(0, pos).c_str()), std::atoi(value.substr(pos + 1).c_str())); + }(); + return !selection || (selection->first == eventID && selection->second == part); +} + +void O2MCApplication::fetchWorkerChunk() +{ + auto& hooks = o2::base::MTWorkerHooks::instance(); + if (!mWorkerInitialized) { // the first sub-event of this thread + mSimDataChannel = hooks.dataChannel(); + for (auto det : listActiveDetectors) { + if (auto o2det = dynamic_cast(det)) { + o2det->initializeLate(); + } + } + mWorkerInitialized = true; + } + auto chunk = hooks.fetchChunk(); + mIsEndOfWork = !chunk; + if (mIsEndOfWork) { + // transport an empty event and leave the event loop after it + mPrimaries.clear(); + mWorkerSubEventInfo = o2::data::SubEventInfo{}; + setSubEventInfo(&mWorkerSubEventInfo); + hooks.stopWorker(); + return; + } + mWorkerSubEventInfo = chunk->mSubEventInfo; + if (isSelectedSubEvent(mWorkerSubEventInfo.eventID, mWorkerSubEventInfo.part)) { + setPrimaries(std::move(chunk->mParticles)); + } else { + setPrimaries(std::vector{}); + LOG(info) << "This chunk will be skipped"; + } + setSubEventInfo(&mWorkerSubEventInfo); + LOG(info) << "Processing " << mPrimaries.size() << " primary particles for event " << mWorkerSubEventInfo.eventID << "/" + << mWorkerSubEventInfo.maxEvents << " part " << mWorkerSubEventInfo.part << "/" << mWorkerSubEventInfo.nparts; + if (!getenv("O2_DISABLE_REPRODUCIBLE_SIM")) { + o2::base::simRandom()->SetSeed(mWorkerSubEventInfo.seed); + o2::base::VMCSeederService::instance().setSeed(); + } +} + } // namespace steer } // namespace o2 diff --git a/run/O2SimDevice.cxx b/run/O2SimDevice.cxx index f2d3f22d7e036..22c3b3b4b4b61 100644 --- a/run/O2SimDevice.cxx +++ b/run/O2SimDevice.cxx @@ -12,6 +12,7 @@ /// @author Sandro Wenzel #include "O2SimDevice.h" +#include "Steer/O2MCApplication.h" #include "../macro/o2sim.C" #include "TVirtualMC.h" #include @@ -22,6 +23,7 @@ #include #include #include +#include #include "DetectorsBase/SimRandom.h" void doLogInfo(int workerID, std::string const& message) @@ -135,28 +137,6 @@ bool O2SimDevice::Kernel(int workerID, fair::mq::Channel& requestchannel, fair:: reproducibleSim = false; } - // Mainly for debugging reasons, we allow to transport - // a specific event + eventpart. This allows to reproduce and debug bugs faster, once - // we know in which precise chunk they occur. The expected format for the environment variable - // is "eventnum:partid". - auto eventselection = getenv("O2SIM_RESTRICT_EVENTPART"); - int focus_on_event = -1; - int focus_on_part = -1; - if (eventselection) { - auto splitString = [](const std::string& str) { - std::pair parts; - size_t pos = str.find(':'); - if (pos != std::string::npos) { - parts.first = str.substr(0, pos); - parts.second = str.substr(pos + 1); - } - return parts; - }; - auto p = splitString(eventselection); - focus_on_event = std::atoi(p.first.c_str()); - focus_on_part = std::atoi(p.second.c_str()); - } - mVMCApp->setSimDataChannel(&dataoutchannel); // we log info with workerID prepended @@ -181,7 +161,8 @@ bool O2SimDevice::Kernel(int workerID, fair::mq::Channel& requestchannel, fair:: << "for event " << info.eventID << "/" << info.maxEvents << " " << "part " << info.part << "/" << info.nparts; - if (eventselection == nullptr || (focus_on_event == info.eventID && focus_on_part == info.part)) { + // for debugging, O2SIM_RESTRICT_EVENTPART="event:part" transports only that sub-event + if (o2::steer::O2MCApplication::isSelectedSubEvent(info.eventID, info.part)) { mVMCApp->setPrimaries(chunk->mParticles); } else { // nothing to transport here @@ -219,7 +200,7 @@ bool O2SimDevice::Kernel(int workerID, fair::mq::Channel& requestchannel, fair:: O2SimDevice::FetchedChunk O2SimDevice::fetchChunk(int workerID, fair::mq::Channel& requestchannel) { - static int counter = 0; + static std::atomic counter{0}; // shared by the transport threads in multithreaded mode fair::mq::MessagePtr request(requestchannel.NewSimpleMessage(PrimaryChunkRequest{workerID, -1, counter++})); // <-- don't need content; channel means -> give primaries fair::mq::Parts reply; diff --git a/run/O2SimDeviceRunner.cxx b/run/O2SimDeviceRunner.cxx index 11da3792d0ac5..45d35d0b080be 100644 --- a/run/O2SimDeviceRunner.cxx +++ b/run/O2SimDeviceRunner.cxx @@ -22,6 +22,8 @@ #include #include #include +#include +#include #include #include #include @@ -35,6 +37,10 @@ #include #endif #include "SimPublishChannelHelper.h" +#include "DetectorsBase/MTWorkerHooks.h" +#include +#include +#include #include "rapidjson/document.h" #include "rapidjson/stringbuffer.h" @@ -265,6 +271,133 @@ bool waitForControlInput(int workerID) return true; } +// Multithreaded Geant4: one process whose Geant4 worker threads each ask the primary server for sub-events +// and send their hits to the merger, as the forked worker processes do otherwise. +std::atomic gServerDone{false}; +std::atomic gNextWorkerID{0}; +std::atomic gChunkRequests{0}; // requests for work in the current run + +// The channels of the transport threads, kept here so that they can be closed before the process exits +struct MTChannels { + std::shared_ptr factory; + std::mutex setupMutex; // serialises the channel creation in the shared factory + std::mutex requestMutex; // one request channel for all threads (with one per thread, only one got work) + std::unique_ptr request; + std::vector> data; // one per thread +}; +MTChannels gMTChannels; + +fair::mq::Channel* makeWorkerChannel(std::unique_ptr& owner, char const* name, char const* type, std::string const& address) +{ + std::lock_guard lock(gMTChannels.setupMutex); + owner = std::make_unique(name, type, gMTChannels.factory); + if (!owner->Connect(address) || !owner->Validate()) { + LOG(fatal) << "Could not connect channel " << name << " of a transport thread to " << address; + } + return owner.get(); +} + +// Closes the channels of the transport threads; with linger -1 this waits until the merger has taken all queued hits +void closeMTChannels() +{ + auto& hooks = o2::base::MTWorkerHooks::instance(); + hooks.fetchChunk = nullptr; + hooks.dataChannel = nullptr; + gMTChannels.request.reset(); + gMTChannels.data.clear(); + gMTChannels.factory.reset(); +} + +// The Geant4 worker threads already exist when this is called; they use these hooks from the first real ProcessRun. +void installMTWorkerHooks(std::string const& serveraddress, std::string const& mergeraddress) +{ + gMTChannels.factory = fair::mq::TransportFactory::CreateTransportFactory("zeromq"); + auto& hooks = o2::base::MTWorkerHooks::instance(); + hooks.fetchChunk = [serveraddress]() -> std::unique_ptr { + thread_local int workerID = gNextWorkerID++; + ++gChunkRequests; + std::lock_guard lock(gMTChannels.requestMutex); + if (!gMTChannels.request) { + makeWorkerChannel(gMTChannels.request, "primary-get", "req", serveraddress); + } + using Status = o2::devices::O2SimDevice::ChunkStatus; + while (true) { + auto fetched = o2::devices::O2SimDevice::fetchChunk(workerID, *gMTChannels.request); + if (fetched.status == Status::Work) { + return std::move(fetched.chunk); + } + if (fetched.status == Status::Stop) { + gServerDone = true; + return nullptr; + } + // Retry and EndOfWork: ask again, as a worker process does; the server answers as idle once done + } + }; + hooks.dataChannel = [mergeraddress]() -> fair::mq::Channel* { + thread_local fair::mq::Channel* channel = nullptr; + if (!channel) { + std::unique_ptr owner; + channel = makeWorkerChannel(owner, "simdata", "push", mergeraddress); + channel->GetSocket().SetLinger(-1); // keep queued hits until the merger has taken them + std::lock_guard lock(gMTChannels.setupMutex); + gMTChannels.data.push_back(std::move(owner)); + } + return channel; + }; +} + +[[noreturn]] void runMT(std::string const& serveraddress, std::string const& mergeraddress, fair::mq::Channel& pubchannel) +{ + installMTWorkerHooks(serveraddress, mergeraddress); + auto sim = getDevice(); + sim->lateInit(); + auto pushchannel = o2::simpubsub::createPUBChannel(o2::simpubsub::getPublishAddress("o2sim-workerinternal", gDriverProcess), "push"); + std::thread([&pubchannel]() { + auto collectorchannel = o2::simpubsub::createPUBChannel(o2::simpubsub::getPublishAddress("o2sim-workerinternal", gDriverProcess), "pull"); + std::unique_ptr msg(collectorchannel.NewMessage()); + while (true) { + if (collectorchannel.Receive(msg) > 0) { + std::string text(reinterpret_cast(msg->GetData()), msg->GetSize()); + o2::simpubsub::publishMessage(pubchannel, text); + } + } + }).detach(); + o2::simpubsub::publishMessage(pushchannel, o2::simpubsub::simStatusString("WORKER0", "STATUS", "SETUP COMPLETED")); + + // a run ends once every thread was told there is no more work; it is repeated if the events ran out first + constexpr int kEventBlock = 100000; + auto vmc = TVirtualMC::GetMC(); + if (auto& beforeTransport = o2::base::MTWorkerHooks::instance().beforeTransport) { + beforeTransport(); + } + auto& conf = o2::conf::SimConfig::Instance(); + bool more = true; + while (more) { + gServerDone = false; + while (!gServerDone) { + gChunkRequests = 0; + vmc->ProcessRun(kEventBlock); + if (gChunkRequests == 0) { + LOG(fatal) << "No transport thread asked for work in this run; Geant4 did not start it"; + } + } + doLogInfo(0, "simulation is done"); + if (conf.asService()) { + LOG(info) << "IN SERVICE MODE WAITING"; + o2::simpubsub::publishMessage(pushchannel, o2::simpubsub::simStatusString("WORKER0", "STATUS", "AWAITING INPUT")); + more = waitForControlInput(0); + usleep(100); + } else { + o2::simpubsub::publishMessage(pushchannel, o2::simpubsub::simStatusString("WORKER0", "STATUS", "TERMINATING")); + LOG(info) << "FINISHING"; + more = false; + } + } + o2::SimSetup::shutdown(); + closeMTChannels(); + _exit(0); +} + int main(int argc, char* argv[]) { struct sigaction act; @@ -374,6 +507,12 @@ int main(int argc, char* argv[]) // many processes. Here we are not using fair::mq::Devices and just setup // some channels manually and do our own runloop. + // in multithreaded mode, Geant4 starts its worker threads while the simulation is initialised + const bool mtMode = getenv("ALICE_O2SIM_MT") != nullptr; + if (mtMode) { + ROOT::EnableThreadSafety(); + } + // we init the simulation first std::unique_ptr simrun; // TODO: take the addresses from somewhere else @@ -382,6 +521,12 @@ int main(int argc, char* argv[]) return 1; } + if (mtMode) { + gMasterProcess = getpid(); + gDriverProcess = getppid(); + runMT(serveraddress, mergeraddress, pubchannel); + } + // should be factored out? unsigned int nworkers = std::max(1u, std::thread::hardware_concurrency() / 2); auto f = getenv("ALICE_NSIMWORKERS"); diff --git a/run/o2sim_parallel.cxx b/run/o2sim_parallel.cxx index 08e3e24935d0c..38112df1fddee 100644 --- a/run/o2sim_parallel.cxx +++ b/run/o2sim_parallel.cxx @@ -533,6 +533,14 @@ int main(int argc, char* argv[]) int nworkers = conf.getNSimWorkers(); setenv("ALICE_NSIMWORKERS", std::to_string(nworkers).c_str(), 1); LOG(info) << "Running with " << nworkers << " sim workers "; + if (conf.getIsMT()) { + if (conf.getMCEngine() != "TGeant4") { + LOG(fatal) << "Multithreaded mode (--isMT) is only available with TGeant4"; + } + // one worker process with a Geant4 transport thread per sim worker + setenv("ALICE_O2SIM_MT", "ON", 1); + LOG(info) << "Multithreaded Geant4 with " << nworkers << " transport threads"; + } o2::utils::ShmManager::Instance().createGlobalSegment(nworkers); @@ -817,7 +825,7 @@ int main(int argc, char* argv[]) } // sum the Geant4 scoring meshes written by the individual workers - if (!errored && o2::conf::mergeG4ScoringDumps(".", conf.getNSimWorkers()) < 0) { + if (!errored && o2::conf::mergeG4ScoringDumps(".", conf.getIsMT() ? 1 : conf.getNSimWorkers()) < 0) { errored = true; } From a2996d129495252987c81f41b3e97ddeb4fe759e Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Sun, 8 Nov 2026 20:40:26 +0100 Subject: [PATCH 2/8] Reset the FT0 photon track ids at the end of each event This fixes hits that depended on the previous event of the same worker. - mTrackIdTop and mTrackIdMCPtop kept their values across events. - A photon in the next event with the same track id then skipped the photo-electron registration. - Both are now reset in Reset(). Co-Authored-By: Claude Opus 5.5 --- Detectors/FIT/FT0/simulation/src/Detector.cxx | 2 ++ 1 file changed, 2 insertions(+) diff --git a/Detectors/FIT/FT0/simulation/src/Detector.cxx b/Detectors/FIT/FT0/simulation/src/Detector.cxx index 813f38c41d712..3a498f8f1c839 100644 --- a/Detectors/FIT/FT0/simulation/src/Detector.cxx +++ b/Detectors/FIT/FT0/simulation/src/Detector.cxx @@ -1218,6 +1218,8 @@ void Detector::Reset() mHits->clear(); } mTrackRefLastTrackID = -1; + mTrackIdTop = -1; + mTrackIdMCPtop = -1; } void Detector::CreateMaterials() From 75cb49963090eb61ca24c03a22470e0b13f84faa Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Sun, 8 Nov 2026 21:41:01 +0100 Subject: [PATCH 3/8] Look up the PDG record once per thread and PDG code in the Stack This removes most TDatabasePDG lookups from Stack::PushTrack. - Each lookup takes ROOT's global lock, which serialises the threads of multithreaded Geant4. - A TParticle is now copied from a per-thread prototype, which keeps the PDG record and the mass. - Unknown PDG codes still use the full TParticle constructor. Co-Authored-By: Claude Opus 5.5 --- Detectors/Base/src/Stack.cxx | 30 +++++++++++++++++++++++++++++- 1 file changed, 29 insertions(+), 1 deletion(-) diff --git a/Detectors/Base/src/Stack.cxx b/Detectors/Base/src/Stack.cxx index 5f6b8c6cf791c..78b7f0baa4264 100644 --- a/Detectors/Base/src/Stack.cxx +++ b/Detectors/Base/src/Stack.cxx @@ -36,6 +36,7 @@ #include #include // for NULL #include +#include #include "DetectorsBase/SimRandom.h" using std::cout; @@ -43,6 +44,33 @@ using std::endl; using std::pair; using namespace o2::data; +namespace +{ +// a TParticle copied from a per-thread prototype per PDG code, avoiding the ROOT-locked PDG lookup +TParticle makeParticle(int pdg, int status, int mother1, int mother2, int daughter1, int daughter2, + double px, double py, double pz, double e, double vx, double vy, double vz, double t) +{ + thread_local std::unordered_map prototypes; + auto it = prototypes.find(pdg); + if (it == prototypes.end()) { + TParticle p(pdg, status, mother1, mother2, daughter1, daughter2, px, py, pz, e, vx, vy, vz, t); + if (p.GetPDG(1)) { // unknown codes may still be added to TDatabasePDG later + prototypes.emplace(pdg, p); + } + return p; + } + TParticle p(it->second); // the copy keeps the PDG record and the mass + p.SetStatusCode(status); + p.SetFirstMother(mother1); + p.SetLastMother(mother2); + p.SetFirstDaughter(daughter1); + p.SetLastDaughter(daughter2); + p.SetMomentum(px, py, pz, e); + p.SetProductionVertex(vx, vy, vz, t); + return p; +} +} // namespace + // small helper function to append to vector at arbitrary position template void insertInVector(std::vector& v, I index, T e) @@ -219,7 +247,7 @@ void Stack::PushTrack(Int_t toBeDone, Int_t parentId, Int_t pdgCode, Double_t px // Int_t daughter1Id = -1; // Int_t daughter2Id = -1; Int_t iStatus = (proc == kPPrimary) ? is : trackId; - TParticle p(pdgCode, iStatus, parentId, secondparentId, daughter1Id, daughter2Id, px, py, pz, e, vx, vy, vz, time); + TParticle p = makeParticle(pdgCode, iStatus, parentId, secondparentId, daughter1Id, daughter2Id, px, py, pz, e, vx, vy, vz, time); p.SetPolarisation(polx, poly, polz); p.SetWeight(weight); p.SetUniqueID(proc); // using the unique ID to transfer process ID From 73f8f9daf0c31cc85a7f8b7750666cf49dcdd686 Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Mon, 9 Nov 2026 10:45:17 +0100 Subject: [PATCH 4/8] Add a test that multithreaded o2-sim gives the same output for 1 and N threads This adds the script o2-sim-mt-equality.sh and a ctest that runs it. - The script runs the same seeded workload with --isMT -j 1 and -j N and compares the kinematics and hits per event and branch. - It reports the first difference and exits with 1 when the outputs differ. - A sequential run with the same seed is the reference for the number of hits, so that a detector clone without hits fails. - The ctest replaces the disabled o2sim_G4_mt test and is registered only for Geant4 built with multithreading. Usage: o2-sim-mt-equality.sh -n 5 -j 4 --field 0 HALL MAG DIPO COMP PIPE ABSO SHIL TPC https://its.cern.ch/jira/browse/O2-7305 Co-Authored-By: Claude Opus 5.5 --- run/CMakeLists.txt | 36 +++++----- run/o2-sim-mt-equality.sh | 141 ++++++++++++++++++++++++++++++++++++++ 2 files changed, 157 insertions(+), 20 deletions(-) create mode 100755 run/o2-sim-mt-equality.sh diff --git a/run/CMakeLists.txt b/run/CMakeLists.txt index 5fda996400ed2..9ad9705360856 100644 --- a/run/CMakeLists.txt +++ b/run/CMakeLists.txt @@ -230,26 +230,21 @@ set_tests_properties(o2sim_G4 G4) set_property(TEST o2sim_G4 APPEND PROPERTY ENVIRONMENT ${G4ENV}) -# # note that the MT is currently only supported in the non FairMQ version -# # switched off -# o2_add_test_command(NAME o2sim_G4_mt -# WORKING_DIRECTORY ${SIMTESTDIR} -# TIMEOUT 400 -# COMMAND $ -# COMMAND_LINE_ARGS -n -# 1 -# -e -# TGeant4 -# --isMT -# on -# -o -# o2simG4MT -# ENVIRONMENT "${SIMENV}" -# LABELS "g4;sim;long") -#set_tests_properties(o2sim_G4_mt -# PROPERTIES PASS_REGULAR_EXPRESSION -# "Macro finished succesfully") -#set_property(TEST o2sim_G4_mt APPEND PROPERTY ENVIRONMENT ${G4ENV}) +# multithreaded Geant4 must give the same output for 1 and 2 transport threads +if(Geant4_multithreaded_FOUND) + o2_add_test_command(NAME o2sim_G4_mt_equality + WORKING_DIRECTORY ${SIMTESTDIR} + TIMEOUT 600 + COMMAND ${CMAKE_SOURCE_DIR}/run/o2-sim-mt-equality.sh + COMMAND_LINE_ARGS --workdir + ${SIMTESTDIR}/mt_equality + ENVIRONMENT "${SIMENV}" + LABELS "g4;sim;long") + set_tests_properties(o2sim_G4_mt_equality + PROPERTIES PASS_REGULAR_EXPRESSION + "EQUALITY TEST PASSED") + set_property(TEST o2sim_G4_mt_equality APPEND PROPERTY ENVIRONMENT ${G4ENV}) +endif() o2_add_test(CheckStackG4 @@ -360,6 +355,7 @@ set_tests_properties(o2sim_G4_checklogs PROPERTIES FIXTURES_REQUIRED G4) endif() +install(PROGRAMS o2-sim-mt-equality.sh DESTINATION ${CMAKE_INSTALL_BINDIR}) install(FILES o2-sim-client.py PERMISSIONS GROUP_READ GROUP_EXECUTE OWNER_EXECUTE OWNER_WRITE OWNER_READ WORLD_EXECUTE WORLD_READ DESTINATION ${CMAKE_INSTALL_BINDIR}) install(DIRECTORY SimExamples/ DESTINATION examples diff --git a/run/o2-sim-mt-equality.sh b/run/o2-sim-mt-equality.sh new file mode 100755 index 0000000000000..834b2156e6d87 --- /dev/null +++ b/run/o2-sim-mt-equality.sh @@ -0,0 +1,141 @@ +#!/usr/bin/env bash +# Checks that multithreaded Geant4 in o2-sim gives the same output for 1 and N transport threads. +# +# Usage: o2-sim-mt-equality.sh [-n events] [-j threads] [-g generator] [--field kG] [--seed seed] +# [--seeding track|subevent] [--tolerance fraction] [--workdir dir] [module ...] +# +# Runs the same seeded workload with --isMT -j 1 and --isMT -j N and compares the kinematics and hit files +# event by event and branch by branch. A sequential run with the same seed is the reference for the number of +# hits of each detector, which must not be zero and must agree within the tolerance (default 0.2). Exits with +# 1 at the first difference. Without modules, the passive geometry with ITS is used. + +nev=2 +nthreads=2 +gen=pythia8pp +field=-5 +seed=4711 +seeding=track +tolerance=0.2 +workdir= +modules=() +while [ $# -gt 0 ]; do + case $1 in + -n) nev=$2; shift 2 ;; + -j) nthreads=$2; shift 2 ;; + -g) gen=$2; shift 2 ;; + --field) field=$2; shift 2 ;; + --seed) seed=$2; shift 2 ;; + --seeding) seeding=$2; shift 2 ;; + --tolerance) tolerance=$2; shift 2 ;; + --workdir) workdir=$2; shift 2 ;; + -h|--help) sed -n '2,10p' $0; exit 0 ;; + *) modules+=("$1"); shift ;; + esac +done +[ ${#modules[@]} -eq 0 ] && modules=(HALL MAG DIPO COMP PIPE ABSO SHIL ITS) +case $seeding in + track) seedconf="SimCutParams.trackSeed=true" ;; + subevent) seedconf="SimCutParams.trackSeed=false" ;; + *) echo "unknown seeding '$seeding', use track or subevent"; exit 1 ;; +esac +o2sim=o2-sim +[ -x "$O2_ROOT/bin/o2-sim" ] && o2sim=$O2_ROOT/bin/o2-sim +[ -z "$workdir" ] && workdir=$(mktemp -d ${TMPDIR:-/tmp}/o2sim-mt-equality-XXXXXX) +mkdir -p $workdir && cd $workdir || exit 1 + +# content hash and number of entries per event and branch of the o2sim tree, independent of the basket layout +cat > hashtree.C <<'EOF' +void hashtree(const char* fname) +{ + TFile f(fname); + auto t = (TTree*)f.Get("o2sim"); + if (!t) { + return; + } + std::vector> objs; + objs.reserve(t->GetListOfBranches()->GetEntries()); // the branch addresses point into objs + for (auto b : *t->GetListOfBranches()) { + auto br = (TBranch*)b; + TClass* cl = nullptr; + EDataType dt; + br->GetExpectedType(cl, dt); + if (!cl) { + continue; + } + objs.emplace_back(cl, cl->New(), br->GetName()); + br->SetAddress(&std::get<1>(objs.back())); + } + for (Long64_t i = 0; i < t->GetEntries(); ++i) { + t->GetEntry(i); + for (auto& [cl, obj, name] : objs) { + TBufferFile buf(TBuffer::kWrite); + buf.WriteObjectAny(obj, cl); + TMD5 md5; + md5.Update((UChar_t*)buf.Buffer(), buf.Length()); + md5.Final(); + printf("HASH event %lld branch %s %s\n", i, name.c_str(), md5.AsString()); + if (auto proxy = cl->GetCollectionProxy()) { + TVirtualCollectionProxy::TPushPop helper(proxy, obj); + printf("COUNT event %lld branch %s %u\n", i, name.c_str(), proxy->Size()); + } + } + } +} +EOF + +run() { # + local dir=$1 + shift + mkdir -p $dir && cd $dir + $o2sim "$@" -n $nev -g $gen -e TGeant4 --seed $seed --field $field -m ${modules[@]} --configKeyValues "$seedconf" > sim.log 2>&1 + local rc=$? + if [ $rc -ne 0 ] || ! grep -q "SIMULATION RETURNED SUCCESFULLY" sim.log; then + echo "o2-sim in $dir failed (exit code $rc), see $workdir/$dir" + exit 1 + fi + for f in o2sim_Kine.root o2sim_Hits*.root; do + [ -f $f ] || continue + root -l -b -q "../hashtree.C(\"$f\")" 2>/dev/null > $f.out + grep "^HASH" $f.out | sed "s/^HASH/$f/" > $f.hashes + grep "^COUNT" $f.out | awk '{n += $NF} END {print n + 0}' > $f.count + done + cd .. +} + +echo "o2-sim MT equality: ${modules[*]}, $nev events of $gen, field $field kG, $seeding seeding, 1 against $nthreads threads, in $workdir" +run seq -j 1 +run mt1 --isMT 1 -j 1 +run mt$nthreads --isMT 1 -j $nthreads +status=0 +for h in mt1/*.hashes; do + f=$(basename $h) + file=${f%.hashes} + if [ ! -s $h ]; then + echo "$file: no entries" + status=1 + elif cmp -s $h mt$nthreads/$f; then + echo "$file: identical ($(wc -l < $h) event-branch hashes)" + else + echo "$file: DIFFERS, first at $(diff $h mt$nthreads/$f | grep -m1 '^<' | cut -d' ' -f3-6)" + status=1 + fi +done +[ $(ls mt$nthreads/*.hashes | wc -l) -eq $(ls mt1/*.hashes | wc -l) ] || { echo "the runs wrote different output files"; status=1; } +# the sequential run as reference for the number of hits, which catches a detector clone that makes no hits +for c in seq/o2sim_Hits*.root.count; do + [ -f "$c" ] || continue + file=$(basename ${c%.count}) + nseq=$(cat $c) + nmt=$(cat mt1/$file.count 2>/dev/null || echo 0) + if [ "$nmt" -eq 0 ]; then + echo "$file: no hits in MT (sequential: $nseq)" + status=1 + elif ! awk -v a=$nseq -v b=$nmt -v t=$tolerance 'BEGIN {exit !(a > 0 && (a - b < 0 ? b - a : a - b) <= t * a)}'; then + echo "$file: $nmt hits in MT against $nseq in sequential mode, beyond the tolerance $tolerance" + status=1 + else + echo "$file: $nmt hits in MT, $nseq in sequential mode" + fi +done +[ $status -eq 0 ] && echo "EQUALITY TEST PASSED" || echo "EQUALITY TEST FAILED" +exit $status From ee9abf30fcb3699edac5311619b89bd7e0b2cf98 Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Wed, 11 Nov 2026 16:58:46 +0100 Subject: [PATCH 5/8] Rebuild the geometry voxels after alignment This fixes a crash of o2-sim in multithreaded Geant4 mode when the alignment is applied. - TGeoPhysicalNode::Align only flags the voxels and assembly bounding boxes of the mother volumes for a rebuild. - The rebuild then happens lazily during navigation, concurrently in several transport threads, and corrupts the heap. - They are now rebuilt once at the end of MisalignGeometry, before transport starts. https://its.cern.ch/jira/browse/O2-7305 Co-Authored-By: Claude Opus 5.5 --- Steer/src/O2MCApplication.cxx | 19 +++++++++++++++++++ 1 file changed, 19 insertions(+) diff --git a/Steer/src/O2MCApplication.cxx b/Steer/src/O2MCApplication.cxx index 54f1760f6af96..b3bbe0e59223e 100644 --- a/Steer/src/O2MCApplication.cxx +++ b/Steer/src/O2MCApplication.cxx @@ -46,6 +46,7 @@ #include "DetectorsBase/GeometryManagerParam.h" #include #include +#include #include #include #include @@ -357,6 +358,24 @@ bool O2MCApplicationBase::MisalignGeometry() // performs possible optimizations (shape replacements on the runtime geometry) fixTGeoRuntimeShapes(); + // alignment leaves voxels and assembly boxes to be rebuilt lazily during navigation, which + // is a data race between Geant4 transport threads; rebuild them here once + TIter nextvol(gGeoManager->GetListOfVolumes()); + while (auto vol = static_cast(nextvol())) { + if (vol->IsAssembly()) { + vol->GetShape()->ComputeBBox(); + } + } + nextvol.Reset(); + int nvoxelized = 0; + while (auto vol = static_cast(nextvol())) { + if (auto voxels = vol->GetVoxels(); voxels && voxels->NeedRebuild()) { + voxels->Voxelize(); + ++nvoxelized; + } + } + LOG(info) << "Rebuilt the voxels of " << nvoxelized << " volumes after alignment"; + // return original return value of misalignment procedure return true; } From 395018ef646a3ed2a13f4a2ad2c09c26c6efff55 Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Wed, 11 Nov 2026 17:02:24 +0100 Subject: [PATCH 6/8] Add a test of o2-sim in multithreaded Geant4 mode This adds the test o2sim_G4_MT, the o2sim_G4 setup with --isMT 1. It runs only with a multithreaded Geant4. https://its.cern.ch/jira/browse/O2-7305 Co-Authored-By: Claude Opus 5.5 --- run/CMakeLists.txt | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/run/CMakeLists.txt b/run/CMakeLists.txt index 9ad9705360856..d0e7a3439597a 100644 --- a/run/CMakeLists.txt +++ b/run/CMakeLists.txt @@ -244,6 +244,36 @@ if(Geant4_multithreaded_FOUND) PROPERTIES PASS_REGULAR_EXPRESSION "EQUALITY TEST PASSED") set_property(TEST o2sim_G4_mt_equality APPEND PROPERTY ENVIRONMENT ${G4ENV}) + + # the o2sim_G4 setup in multithreaded Geant4 mode + o2_add_test_command(NAME o2sim_G4_MT + WORKING_DIRECTORY ${SIMTESTDIR} + TIMEOUT 400 + COMMAND $ + COMMAND_LINE_ARGS --isMT + 1 + -n + 2 + -j + 2 + -e + TGeant4 + -o + o2simG4MT + --chunkSize + 2 + --skipModules + MFT ZDC + --seed + 15946057944514955802 + --configKeyValues + "align-geom.mDetectors=none" + ENVIRONMENT "${SIMENV}" + LABELS "g4;sim;long") + set_tests_properties(o2sim_G4_MT + PROPERTIES PASS_REGULAR_EXPRESSION + "SIMULATION RETURNED SUCCESFULLY") + set_property(TEST o2sim_G4_MT APPEND PROPERTY ENVIRONMENT ${G4ENV}) endif() From 030e95221ecbe6376710db2e37efc4c7d5afd7fa Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Wed, 11 Nov 2026 17:04:44 +0100 Subject: [PATCH 7/8] Add the copyright header to commonConfig.C and remove trailing spaces --- Detectors/gconfig/commonConfig.C | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/Detectors/gconfig/commonConfig.C b/Detectors/gconfig/commonConfig.C index 99e495b182dd7..c554cbbdc800b 100644 --- a/Detectors/gconfig/commonConfig.C +++ b/Detectors/gconfig/commonConfig.C @@ -1,3 +1,14 @@ +// Copyright 2020-2022 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + // common piece of code to setup stack and register // with VMC instances @@ -25,7 +36,7 @@ void stackSetup(T* vmc, R* run) // register the stack as an observer on FinishPrimary events (managed by Cave) bool foundCave = false; auto modules = run->GetListOfModules(); - + if( strcmp(vmc->GetName(), "TGeant4")==0 ) { // there is no way to get the module by name, so we have // to iterate through the complete list From 7cf50977a2a89b709ec9f7bdf80505b862fe535d Mon Sep 17 00:00:00 2001 From: Sandro Wenzel Date: Wed, 11 Nov 2026 20:11:46 +0100 Subject: [PATCH 8/8] Do not run the multithreaded Geant4 tests on macOS --- run/CMakeLists.txt | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/run/CMakeLists.txt b/run/CMakeLists.txt index d0e7a3439597a..53223a9925ba9 100644 --- a/run/CMakeLists.txt +++ b/run/CMakeLists.txt @@ -230,8 +230,8 @@ set_tests_properties(o2sim_G4 G4) set_property(TEST o2sim_G4 APPEND PROPERTY ENVIRONMENT ${G4ENV}) -# multithreaded Geant4 must give the same output for 1 and 2 transport threads -if(Geant4_multithreaded_FOUND) +# multithreaded Geant4 must give the same output for 1 and 2 transport threads (not run on macOS) +if(Geant4_multithreaded_FOUND AND NOT APPLE) o2_add_test_command(NAME o2sim_G4_mt_equality WORKING_DIRECTORY ${SIMTESTDIR} TIMEOUT 600