Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // Class G4GeometryManager implementation << 27 // 26 // 28 // 26.07.95, P.Kent - Initial version, includi << 27 // $Id: G4GeometryManager.cc 67975 2013-03-13 10:19:44Z gcosmo $ 29 // 12.06.24, J.Apostolakis - Added parallel op << 28 // >> 29 // class G4GeometryManager >> 30 // >> 31 // Implementation >> 32 // >> 33 // Author: >> 34 // 26.07.95 P.Kent Initial version, including optimisation Build 30 // ------------------------------------------- 35 // -------------------------------------------------------------------- 31 36 32 #include <iomanip> 37 #include <iomanip> 33 << 34 #include "G4ios.hh" << 35 #include "G4Timer.hh" 38 #include "G4Timer.hh" 36 #include "G4GeometryManager.hh" 39 #include "G4GeometryManager.hh" 37 #include "G4SystemOfUnits.hh" 40 #include "G4SystemOfUnits.hh" 38 #include "G4Threading.hh" << 41 >> 42 #ifdef G4GEOMETRY_VOXELDEBUG >> 43 #include "G4ios.hh" >> 44 #endif 39 45 40 // Needed for building optimisations 46 // Needed for building optimisations 41 // 47 // 42 #include "G4LogicalVolumeStore.hh" 48 #include "G4LogicalVolumeStore.hh" 43 #include "G4VPhysicalVolume.hh" 49 #include "G4VPhysicalVolume.hh" 44 #include "G4SmartVoxelHeader.hh" 50 #include "G4SmartVoxelHeader.hh" 45 #include "voxeldefs.hh" 51 #include "voxeldefs.hh" 46 52 47 // Needed for setting the extent for tolerance 53 // Needed for setting the extent for tolerance value 48 // 54 // 49 #include "G4GeometryTolerance.hh" 55 #include "G4GeometryTolerance.hh" 50 #include "G4SolidStore.hh" 56 #include "G4SolidStore.hh" 51 #include "G4VSolid.hh" 57 #include "G4VSolid.hh" 52 58 53 // Needed for parallel optimisation << 54 #include "G4AutoLock.hh" << 55 << 56 namespace // Data structures / mutexes for pa << 57 { << 58 // Mutex to obtain a volume to optimise << 59 G4Mutex obtainVolumeMutex = G4MUTEX_INITIALI << 60 << 61 // Mutex to lock saving of voxel statistics << 62 G4Mutex voxelStatsMutex = G4MUTEX_INITIALIZE << 63 << 64 // Mutex to provide Statistics Results << 65 G4Mutex statResultsMutex = G4MUTEX_INITIALIZ << 66 << 67 // Mutex to start wall clock (global) timer << 68 G4Mutex wallClockTimerMutex = G4MUTEX_INITIA << 69 << 70 // Mutex to write debug output << 71 G4Mutex outputDbgMutex = G4MUTEX_INITIALIZER << 72 } << 73 << 74 // ******************************************* 59 // *************************************************************************** 75 // Static class data << 60 // Static class variable: ptr to single instance of class 76 // ******************************************* 61 // *************************************************************************** 77 // 62 // 78 G4ThreadLocal G4GeometryManager* G4GeometryMan << 63 G4ThreadLocal G4GeometryManager* G4GeometryManager::fgInstance = 0; 79 << 80 // Static *global* class data << 81 G4bool G4GeometryManager::fParallelVoxelOptimi << 82 // Records User choice to use parallel voxel << 83 << 84 G4bool G4GeometryManager::fOptimiseInParallel << 85 // Configured = requested && available (ie i << 86 // Value calculated during each effort to op << 87 << 88 std::vector<G4LogicalVolume*> G4GeometryManage << 89 std::vector<G4LogicalVolume*>::const_iterator << 90 << 91 std::vector<G4SmartVoxelStat> G4GeometryManage << 92 // Container for statistics << 93 << 94 // Derived state << 95 G4bool G4GeometryManager::fVerboseParallel = f << 96 G4bool G4GeometryManager::fParallelVoxelOptimi << 97 G4bool G4GeometryManager::fParallelVoxelOptimi << 98 G4double G4GeometryManager::fSumVoxelTime = 0. << 99 << 100 G4int G4GeometryManager::fNumberThreadsReporti << 101 unsigned int G4GeometryManager::fTotalNumberVo << 102 << 103 // For Wall clock << 104 G4Timer* G4GeometryManager::fWallClockTimer = << 105 G4bool G4GeometryManager::fWallClockStarted = << 106 64 107 // ******************************************* 65 // *************************************************************************** 108 // Destructor << 66 // Constructor. Set the geometry to be open 109 // ******************************************* 67 // *************************************************************************** 110 // 68 // 111 G4GeometryManager::~G4GeometryManager() << 69 G4GeometryManager::G4GeometryManager() >> 70 : fIsClosed(false) 112 { 71 { 113 fgInstance = nullptr; << 114 fIsClosed = false; << 115 << 116 if( fWallClockTimer && G4Threading::IsMaster << 117 { << 118 delete fWallClockTimer; << 119 fWallClockTimer= nullptr; << 120 } << 121 } 72 } 122 73 123 // ******************************************* 74 // *************************************************************************** 124 // Closes geometry - performs sanity checks an 75 // Closes geometry - performs sanity checks and optionally builds optimisation 125 // for placed volumes (always built for replic 76 // for placed volumes (always built for replicas & parameterised). 126 // NOTE: Currently no sanity checks are perfor 77 // NOTE: Currently no sanity checks are performed. 127 // Applies to just a specific subtree if a phy 78 // Applies to just a specific subtree if a physical volume is specified. 128 // ******************************************* 79 // *************************************************************************** 129 // 80 // 130 G4bool G4GeometryManager::CloseGeometry(G4bool 81 G4bool G4GeometryManager::CloseGeometry(G4bool pOptimise, G4bool verbose, 131 G4VPhy 82 G4VPhysicalVolume* pVolume) 132 { 83 { 133 if (!fIsClosed && G4Threading::IsMasterThrea << 84 if (!fIsClosed) 134 { 85 { 135 if (pVolume != nullptr) << 86 if (pVolume) 136 { 87 { 137 BuildOptimisations(pOptimise, pVolume); 88 BuildOptimisations(pOptimise, pVolume); 138 } 89 } 139 else 90 else 140 { 91 { 141 BuildOptimisations(pOptimise, verbose); 92 BuildOptimisations(pOptimise, verbose); 142 } 93 } 143 fIsClosed = true; << 94 fIsClosed=true; 144 } 95 } 145 return true; 96 return true; 146 } 97 } 147 98 148 // ******************************************* 99 // *************************************************************************** 149 // Opens the geometry and removes optimisation 100 // Opens the geometry and removes optimisations (optionally, related to just 150 // the specified logical-volume). 101 // the specified logical-volume). 151 // Applies to just a specific subtree if a phy 102 // Applies to just a specific subtree if a physical volume is specified. 152 // ******************************************* 103 // *************************************************************************** 153 // 104 // 154 void G4GeometryManager::OpenGeometry(G4VPhysic 105 void G4GeometryManager::OpenGeometry(G4VPhysicalVolume* pVolume) 155 { 106 { 156 if (fIsClosed && G4Threading::IsMasterThread << 107 if (fIsClosed) 157 { 108 { 158 if (pVolume != nullptr) << 109 if (pVolume) 159 { 110 { 160 DeleteOptimisations(pVolume); 111 DeleteOptimisations(pVolume); 161 } 112 } 162 else 113 else 163 { 114 { 164 DeleteOptimisations(); 115 DeleteOptimisations(); 165 } 116 } 166 fIsClosed = false; << 117 fIsClosed=false; 167 } 118 } 168 } 119 } 169 120 170 // ******************************************* 121 // *************************************************************************** 171 // Returns the instance of the singleton. << 122 // Returns status of geometry 172 // Creates it in case it's called for the firs << 173 // ******************************************* 123 // *************************************************************************** 174 // 124 // 175 G4GeometryManager* G4GeometryManager::GetInsta << 125 G4bool G4GeometryManager::IsGeometryClosed() 176 { 126 { 177 if (fgInstance == nullptr) << 127 return fIsClosed; 178 { << 179 fgInstance = new G4GeometryManager; << 180 << 181 if( (fWallClockTimer == nullptr) && G4Thre << 182 { << 183 fWallClockTimer = new G4Timer; << 184 } << 185 } << 186 return fgInstance; << 187 } 128 } 188 129 189 // ******************************************* 130 // *************************************************************************** 190 // Returns the instance of the singleton. 131 // Returns the instance of the singleton. >> 132 // Creates it in case it's called for the first time. 191 // ******************************************* 133 // *************************************************************************** 192 // 134 // 193 G4GeometryManager* G4GeometryManager::GetInsta << 135 G4GeometryManager* G4GeometryManager::GetInstance() 194 { << 195 return fgInstance; << 196 } << 197 << 198 // ******************************************* << 199 // Simplest user method to request parallel op << 200 // ******************************************* << 201 // << 202 void G4GeometryManager::OptimiseInParallel( G4 << 203 { << 204 RequestParallelOptimisation(val); << 205 } << 206 << 207 // ******************************************* << 208 // Report about Voxel(isation) of a logical vo << 209 // ******************************************* << 210 // << 211 void << 212 G4GeometryManager::ReportVoxelInfo(G4LogicalVo << 213 { 136 { 214 G4SmartVoxelHeader* head = logVolume->GetVox << 137 if (!fgInstance) 215 if( head != nullptr ) << 216 { << 217 os << "** Created optimisations for logica << 218 << std::setw(50) << logVolume->GetName( << 219 << "- Result VoxelInfo - START: " << " << 220 << *head << 221 << "- Result VoxelInfo - END. " << G4 << 222 } << 223 else << 224 { 138 { 225 os << "** No optimisation for log-vol " << << 139 fgInstance = new G4GeometryManager; 226 } 140 } 227 os << "*** Report Voxel Info: END " << G4en << 141 return fgInstance; 228 } 142 } 229 143 230 // ******************************************* 144 // *************************************************************************** 231 // Creates optimisation info. Builds all voxel 145 // Creates optimisation info. Builds all voxels if allOpts=true 232 // otherwise it builds voxels only for replica 146 // otherwise it builds voxels only for replicated volumes. 233 // ******************************************* 147 // *************************************************************************** 234 // Returns whether optimisation is finished << 235 // 148 // 236 G4bool G4GeometryManager::BuildOptimisations(G << 149 void G4GeometryManager::BuildOptimisations(G4bool allOpts, G4bool verbose) 237 { 150 { 238 G4bool finishedOptimisation = false; << 151 G4Timer timer; 239 << 152 G4Timer allTimer; 240 fOptimiseInParallelConfigured = fParallelVox << 153 std::vector<G4SmartVoxelStat> stats; 241 && G4Threading: << 154 if (verbose) { allTimer.Start(); } 242 << 155 243 static unsigned int NumCallsBuildOptimisatio << 156 G4LogicalVolumeStore* Store = G4LogicalVolumeStore::GetInstance(); 244 if( fOptimiseInParallelConfigured && (NumCal << 157 G4LogicalVolume* volume; 245 { << 158 G4SmartVoxelHeader* head; 246 PrepareParallelOptimisation(allOpts, verbo << 159 247 ++NumCallsBuildOptimisations; << 160 for (size_t n=0; n<Store->size(); n++) 248 } << 161 { 249 else << 162 if (verbose) timer.Start(); 250 { << 163 volume=(*Store)[n]; 251 BuildOptimisationsSequential(allOpts, verb << 164 // For safety, check if there are any existing voxels and 252 finishedOptimisation= true; << 165 // delete before replacement 253 } << 166 // 254 << 167 head = volume->GetVoxelHeader(); 255 return finishedOptimisation; << 168 delete head; 256 } << 169 volume->SetVoxelHeader(0); 257 << 170 if ( ( (volume->IsToOptimise()) 258 // ******************************************* << 171 && (volume->GetNoDaughters()>=kMinVoxelVolumesLevel1&&allOpts) ) 259 // Creates optimisation info. Builds all voxel << 172 || ( (volume->GetNoDaughters()==1) 260 // otherwise it builds voxels only for replica << 173 && (volume->GetDaughter(0)->IsReplicated()==true) 261 // << 174 && (volume->GetDaughter(0)->GetRegularStructureId()!=1) ) ) 262 // This is the original sequential implementat << 175 { 263 // - at first initialisation to create voxels, << 264 // - at re-initialisation if the geometry has << 265 // ******************************************* << 266 // << 267 void G4GeometryManager::BuildOptimisationsSequ << 268 << 269 { << 270 G4Timer timer; << 271 G4Timer allTimer; << 272 std::vector<G4SmartVoxelStat> stats; << 273 << 274 if (verbose) { allTimer.Start(); } << 275 << 276 G4LogicalVolumeStore* Store = G4LogicalVolum << 277 G4LogicalVolume* volume; << 278 G4SmartVoxelHeader* head; << 279 << 280 #ifdef G4GEOMETRY_VOXELDEBUG << 281 G4cout << G4endl << 282 << "*** G4GeometryManager::BuildOptimisat << 283 << G4Threading::G4GetThreadId() << " all- << 284 #endif << 285 << 286 for (auto & n : *Store) << 287 { << 288 if (verbose) timer.Start(); << 289 volume=n; << 290 // For safety, check if there are any exis << 291 // delete before replacement << 292 // << 293 head = volume->GetVoxelHeader(); << 294 delete head; << 295 volume->SetVoxelHeader(nullptr); << 296 if ( ( (volume->IsToOptimise()) << 297 && (volume->GetNoDaughters()>=kMi << 298 || ( (volume->GetNoDaughters()==1) << 299 && (volume->GetDaughter(0)->IsRepl << 300 && (volume->GetDaughter(0)->GetReg << 301 { << 302 #ifdef G4GEOMETRY_VOXELDEBUG 176 #ifdef G4GEOMETRY_VOXELDEBUG 303 G4cout << "** G4GeometryManager::BuildOpti << 177 G4cout << "**** G4GeometryManager::BuildOptimisations" << G4endl 304 << " Examining logical volume nam << 178 << " Examining logical volume name = " 305 << "' #daughters= " << volume->Get << 179 << volume->GetName() << G4endl; 306 #endif 180 #endif 307 head = new G4SmartVoxelHeader(volume); << 181 head = new G4SmartVoxelHeader(volume); 308 << 182 if (head) 309 if (head != nullptr) << 183 { 310 { << 184 volume->SetVoxelHeader(head); 311 volume->SetVoxelHeader(head); << 185 } 312 } << 186 else 313 else << 187 { 314 { << 188 std::ostringstream message; 315 std::ostringstream message; << 189 message << "VoxelHeader allocation error." << G4endl 316 message << "VoxelHeader allocation err << 190 << "Allocation of new VoxelHeader" << G4endl 317 << "Allocation of new VoxelHea << 191 << " for volume " << volume->GetName() << " failed."; 318 << " for volume '" << v << 192 G4Exception("G4GeometryManager::BuildOptimisations()", "GeomMgt0003", 319 G4Exception("G4GeometryManager::BuildO << 193 FatalException, message); 320 FatalException, message); << 194 } 321 } << 195 if (verbose) 322 if (verbose) << 196 { 323 { << 197 timer.Stop(); 324 timer.Stop(); << 198 stats.push_back( G4SmartVoxelStat( volume, head, 325 stats.emplace_back( volume, head, << 199 timer.GetSystemElapsed(), 326 timer.GetSystemElap << 200 timer.GetUserElapsed() ) ); 327 timer.GetUserElapse << 201 } 328 } << 202 } 329 } << 203 else 330 else << 204 { 331 { << 205 // Don't create voxels for this node 332 // Don't create voxels for this node << 333 #ifdef G4GEOMETRY_VOXELDEBUG 206 #ifdef G4GEOMETRY_VOXELDEBUG 334 auto numDaughters = volume->GetNoDaughte << 207 G4cout << "**** G4GeometryManager::BuildOptimisations" << G4endl 335 G4cout << "- Skipping logical volume wit << 208 << " Skipping logical volume name = " << volume->GetName() 336 << " daughters and name = '" << v << 209 << G4endl; 337 if( numDaughters > 1 ) << 338 { << 339 G4cout << "[Placement]"; << 340 } << 341 else << 342 { << 343 if( numDaughters == 1 ) << 344 { << 345 G4cout << ( volume->GetDaughter(0)-> << 346 << 347 } << 348 } << 349 G4cout << G4endl; << 350 #endif 210 #endif 351 } << 211 } 352 } 212 } 353 if (verbose) 213 if (verbose) 354 { 214 { 355 allTimer.Stop(); << 215 allTimer.Stop(); 356 << 216 ReportVoxelStats( stats, allTimer.GetSystemElapsed() 357 ReportVoxelStats( stats, allTimer.GetSyste << 217 + allTimer.GetUserElapsed() ); 358 + allTimer.GetUserElapsed << 359 } 218 } 360 } 219 } 361 220 362 // ******************************************* 221 // *************************************************************************** 363 // Creates a list of logical volumes which wil << 222 // Creates optimisation info for the specified volumes subtree. 364 // if allOpts=true it lists all voxels << 365 // otherwise it lists only the voxels << 366 // This list will be used subsequently to buil << 367 // << 368 // Note: this method is NOT thread safe! << 369 // It expects to be called only once in eac << 370 // i.e. either by master thread or a select << 371 // ******************************************* << 372 // << 373 void << 374 G4GeometryManager::CreateListOfVolumesToOptimi << 375 { << 376 // Prepare the work - must be called only in << 377 << 378 G4LogicalVolumeStore* Store = G4LogicalVolum << 379 << 380 if( !fVolumesToOptimise.empty() ) << 381 { << 382 ResetListOfVolumesToOptimise(); << 383 } << 384 << 385 for (auto & n : *Store) << 386 { << 387 G4LogicalVolume* volume=n; << 388 << 389 if ( ( (volume->IsToOptimise()) << 390 && (volume->GetNoDaughters()>=kMi << 391 || ( (volume->GetNoDaughters()==1) << 392 && (volume->GetDaughter(0)->IsRepl << 393 && (volume->GetDaughter(0)->GetReg << 394 { << 395 fVolumesToOptimise.push_back(volume); << 396 << 397 // For safety, must check (later) if the << 398 // delete before replacement: << 399 // All 'clients' of this code must do th << 400 // delete volume->GetVoxelHeader(); << 401 // volume->SetVoxelHeader(nullptr); << 402 << 403 #ifdef G4GEOMETRY_VOXELDEBUG << 404 G4cout << "- Booking logical volume wit << 405 << " daughters and name = '" << volume-> << 406 << " -- for optimisation (ie voxels will << 407 #endif << 408 } << 409 else << 410 { << 411 #ifdef G4GEOMETRY_VOXELDEBUG << 412 G4cout << "- Skipping logical volume wit << 413 << " daughters and name = '" << volume-> << 414 #endif << 415 } << 416 } << 417 << 418 if(verbose) << 419 G4cout << "** G4GeometryManager::PrepareOp << 420 << " Number of volumes for voxelis << 421 << fVolumesToOptimise.size() << G4e << 422 << 423 fLogVolumeIterator = fVolumesToOptimise.cbeg << 424 } << 425 << 426 // ******************************************* << 427 // Obtain a logical volume from the list of vo << 428 // Must be thread-safe: its role is to be call << 429 // Critical method for parallel optimisation - << 430 // ******************************************* << 431 // << 432 G4LogicalVolume* G4GeometryManager::ObtainVolu << 433 { << 434 G4LogicalVolume* logVolume = nullptr; << 435 << 436 G4AutoLock lock(obtainVolumeMutex); << 437 << 438 if( fLogVolumeIterator != fVolumesToOptimise << 439 { << 440 logVolume = *fLogVolumeIterator; << 441 ++fLogVolumeIterator; << 442 } << 443 return logVolume; << 444 } << 445 << 446 // ******************************************* << 447 // Thread-safe method to clear the list of vol << 448 // ******************************************* << 449 // << 450 void G4GeometryManager::ResetListOfVolumesToOp << 451 { << 452 G4AutoLock lock(obtainVolumeMutex); << 453 << 454 std::vector<G4LogicalVolume*>().swap(fVolume << 455 // Swapping with an empty vector in order to << 456 // without calling destructors of logical vo << 457 // Must not call clear: i.e. fVolumesToOptim << 458 << 459 assert(fVolumesToOptimise.empty()); << 460 fLogVolumeIterator = fVolumesToOptimise.cbeg << 461 << 462 fGlobVoxelStats.clear(); << 463 // Reset also the statistics of volumes -- t << 464 } << 465 << 466 // ******************************************* << 467 // Method which user calls to ask for parallel << 468 // ******************************************* << 469 // << 470 void G4GeometryManager::RequestParallelOptimis << 471 { << 472 fParallelVoxelOptimisationRequested = flag; << 473 if( flag ) << 474 { << 475 ConfigureParallelOptimisation(verbose); << 476 } << 477 } << 478 << 479 // ******************************************* << 480 // Setup up state to enable parallel optimisat << 481 // ******************************************* << 482 // << 483 void G4GeometryManager::ConfigureParallelOptim << 484 { << 485 if(verbose) << 486 { << 487 G4cout << "** G4GeometryManager::Configure << 488 << " LEAVING all the work (of voxel optimi << 489 << G4endl; << 490 } << 491 fParallelVoxelOptimisationRequested = true; << 492 fParallelVoxelOptimisationUnderway = false; << 493 fParallelVoxelOptimisationFinished = false; << 494 << 495 // Keep values of options / verbosity for us << 496 fVerboseParallel = verbose; << 497 << 498 // New effort -- reset the total time -- and << 499 fSumVoxelTime = 0.0; << 500 fNumberThreadsReporting = 0; << 501 fTotalNumberVolumesOptimised = 0; // Numbe << 502 << 503 fWallClockStarted = false; // Will need to << 504 } << 505 << 506 // ******************************************* << 507 // Build voxel optimisation in parallel -- pre << 508 // ******************************************* << 509 // << 510 void << 511 G4GeometryManager::PrepareParallelOptimisation << 512 { << 513 if( verbose ) << 514 { << 515 G4cout << "** G4GeometryManager::PreparePa << 516 << G4endl; << 517 } << 518 CreateListOfVolumesToOptimise(allOpts, verbo << 519 ConfigureParallelOptimisation(verbose); << 520 } << 521 << 522 // ******************************************* << 523 // Method for a thread/task to contribute dyna << 524 // ******************************************* << 525 // << 526 void G4GeometryManager::UndertakeOptimisation( << 527 { << 528 G4bool verbose = fVerboseParallel; << 529 G4LogicalVolume* logVolume = nullptr; << 530 << 531 fParallelVoxelOptimisationUnderway = true; << 532 << 533 // Start timer - if not already done << 534 if( ( !fWallClockStarted ) && verbose ) << 535 { << 536 G4AutoLock startTimeLock(wallClockTimerMut << 537 if( !fWallClockStarted ) << 538 { << 539 fWallClockTimer->Start(); << 540 fWallClockStarted= true; << 541 } << 542 } << 543 << 544 G4Timer fetimer; << 545 unsigned int numVolumesOptimised = 0; << 546 << 547 while( (logVolume = ObtainVolumeToOptimise() << 548 { << 549 if (verbose) fetimer.Start(); << 550 << 551 G4SmartVoxelHeader* head = logVolume->GetV << 552 delete head; << 553 logVolume->SetVoxelHeader(nullptr); << 554 << 555 head = new G4SmartVoxelHeader(logVolume); << 556 // ********************************* << 557 logVolume->SetVoxelHeader(head); << 558 << 559 if (head != nullptr) << 560 { << 561 ++numVolumesOptimised; << 562 } << 563 else << 564 { << 565 G4ExceptionDescription message; << 566 message << "VoxelHeader allocation error << 567 << "Allocation of new VoxelHeade << 568 << "for logical volume " << logV << 569 G4Exception("G4GeometryManager::BuildOpt << 570 "GeomMgt0003", FatalExceptio << 571 } << 572 << 573 if(verbose) << 574 { << 575 fetimer.Stop(); << 576 auto feRealElapsed = fetimer.GetRealElap << 577 // Must use 'real' elapsed time -- canno << 578 // (it accounts for all threads) << 579 << 580 G4AutoLock lock(voxelStatsMutex); << 581 fGlobVoxelStats.emplace_back( logVolume, << 582 0.0, // << 583 feRealElapsed ); // << 584 fSumVoxelTime += feRealElapsed; << 585 } << 586 } << 587 << 588 G4bool allDone = false; << 589 G4int myCount= -1; << 590 << 591 myCount = ReportWorkerIsDoneOptimising(numVo << 592 allDone = IsParallelOptimisationFinished(); << 593 << 594 if( allDone && (myCount == G4Threading::GetN << 595 { << 596 G4int badVolumes = CheckOptimisation(); // << 597 if( badVolumes > 0 ) << 598 { << 599 G4ExceptionDescription errmsg; << 600 errmsg <<" Expected that all voxelisatio << 601 << "but found that voxels headers << 602 << badVolumes << " volumes."; << 603 G4Exception("G4GeometryManager::Undertak << 604 "GeomMng002", FatalException << 605 } << 606 << 607 // Create report << 608 << 609 if( verbose ) << 610 { << 611 fWallClockTimer->Stop(); << 612 << 613 std::ostream& report_stream = std::cout; << 614 report_stream << G4endl << 615 << "G4GeometryManager::UndertakeOptimi << 616 << " -- Timing for Voxel Optimisation" << 617 report_stream << " - Elapsed time (real << 618 << fWallClockTimer->GetRealElapsed() < << 619 << ", user " << fWallClockTimer->GetUs << 620 << ", system " << fWallClockTimer->Get << 621 << G4endl; << 622 report_stream << " - Sum voxel time (re << 623 << "s."; << 624 report_stream << std::setprecision(6) << << 625 << 626 ReportVoxelStats( fGlobVoxelStats, fSumV << 627 report_stream.flush(); << 628 } << 629 } << 630 else << 631 { << 632 WaitForVoxelisationFinish(false); << 633 } << 634 } << 635 << 636 // ******************************************* << 637 // Ensure that all the work of voxelisation is << 638 // Can be called in GeometryManager methods or << 639 // ******************************************* << 640 // << 641 void G4GeometryManager::WaitForVoxelisationFin << 642 { << 643 // Must wait until all workers are done ... << 644 using namespace std::chrono_literals; << 645 unsigned int trials = 0; << 646 auto tid = G4Threading::G4GetThreadId(); << 647 << 648 std::ostream& out_stream = std::cout; // G4c << 649 while( ! IsParallelOptimisationFinished() ) << 650 { << 651 // Each thread must wait until all are don << 652 std::this_thread::sleep_for(250ms); << 653 ++trials; << 654 } << 655 << 656 if( verbose ) << 657 { << 658 G4AutoLock lock(outputDbgMutex); << 659 out_stream << G4endl << 660 << "** UndertakeOptimisation do << 661 << " after waiting for " << tr << 662 out_stream.flush(); << 663 } << 664 } << 665 << 666 // ******************************************* << 667 // Ensure that all logical volumes in list hav << 668 // ******************************************* << 669 // << 670 G4int G4GeometryManager::CheckOptimisation() << 671 { << 672 unsigned int numErrors = 0; << 673 for ( const auto& logical : fVolumesToOptimi << 674 { << 675 if( logical->GetVoxelHeader() == nullptr ) << 676 } << 677 return numErrors; << 678 } << 679 << 680 // ******************************************* << 681 // Report that current thread/task is done opt << 682 // A thread call this method to reports that i << 683 // many volumes it optimised. The method: << 684 // - increments the count of workers that ha << 685 // - keeps count of number of volumes optimi << 686 // - if all works is done (ie all workers ha << 687 // in the 'Finished' state. << 688 // ******************************************* << 689 // << 690 G4int << 691 G4GeometryManager::ReportWorkerIsDoneOptimisin << 692 { << 693 // Check that all are done and, if so, signa << 694 G4int orderReporting; << 695 << 696 G4AutoLock lock(statResultsMutex); << 697 orderReporting = ++fNumberThreadsReporting; << 698 fTotalNumberVolumesOptimised += numVolumesOp << 699 << 700 if (fNumberThreadsReporting == G4Threading:: << 701 { << 702 InformOptimisationIsFinished(fVerboseParal << 703 } << 704 << 705 return orderReporting; << 706 } << 707 << 708 // ******************************************* << 709 // Inform that all work for parallel optimisat << 710 // ******************************************* << 711 // << 712 void G4GeometryManager::InformOptimisationIsFi << 713 { << 714 if(verbose) // G4cout does not work! << 715 { << 716 std::cout << "** G4GeometryManager: All vo << 717 << G4endl; << 718 std::cout << " Total number of volumes o << 719 << fTotalNumberVolumesOptimised << 720 << " of " << fVolumesToOptimise. << 721 std::cout << " Number of workers reporti << 722 << fNumberThreadsReporting << 723 << " of " << G4Threading::GetNum << 724 << " expected\n"; << 725 } << 726 assert ( fTotalNumberVolumesOptimised == fVo << 727 << 728 fParallelVoxelOptimisationFinished = true; << 729 // fParallelVoxelOptimisationRequested = fal << 730 fParallelVoxelOptimisationUnderway = false; << 731 } << 732 << 733 // ******************************************* << 734 // Creates Optimisation info for the specified << 735 // ******************************************* 223 // *************************************************************************** 736 // 224 // 737 void G4GeometryManager::BuildOptimisations(G4b 225 void G4GeometryManager::BuildOptimisations(G4bool allOpts, 738 G4V 226 G4VPhysicalVolume* pVolume) 739 { 227 { 740 if (pVolume == nullptr) { return; } << 228 if (!pVolume) { return; } 741 229 742 // Retrieve the mother logical volume, if no << 230 // Retrieve the mother logical volume, if not NULL, 743 // otherwise apply global optimisation for t << 231 // otherwise apply global optimisation for the world volume 744 // << 232 // 745 G4LogicalVolume* tVolume = pVolume->GetMothe << 233 G4LogicalVolume* tVolume = pVolume->GetMotherLogical(); 746 if (tVolume == nullptr) << 234 if (!tVolume) { return BuildOptimisations(allOpts, false); } 747 { << 235 748 BuildOptimisations(allOpts, false); << 236 G4SmartVoxelHeader* head = tVolume->GetVoxelHeader(); 749 return; << 237 delete head; 750 } << 238 tVolume->SetVoxelHeader(0); 751 << 239 if ( ( (tVolume->IsToOptimise()) 752 G4SmartVoxelHeader* head = tVolume->GetVoxel << 240 && (tVolume->GetNoDaughters()>=kMinVoxelVolumesLevel1&&allOpts) ) 753 delete head; << 241 || ( (tVolume->GetNoDaughters()==1) 754 tVolume->SetVoxelHeader(nullptr); << 242 && (tVolume->GetDaughter(0)->IsReplicated()==true) ) ) 755 if ( ( (tVolume->IsToOptimise()) << 243 { 756 && (tVolume->GetNoDaughters()>=kMinVo << 244 head = new G4SmartVoxelHeader(tVolume); 757 || ( (tVolume->GetNoDaughters()==1) << 245 if (head) 758 && (tVolume->GetDaughter(0)->IsReplic << 246 { 759 { << 247 tVolume->SetVoxelHeader(head); 760 head = new G4SmartVoxelHeader(tVolume); << 248 } 761 if (head != nullptr) << 249 else 762 { << 250 { 763 tVolume->SetVoxelHeader(head); << 251 std::ostringstream message; 764 } << 252 message << "VoxelHeader allocation error." << G4endl 765 else << 253 << "Allocation of new VoxelHeader" << G4endl 766 { << 254 << " for volume " << tVolume->GetName() << " failed."; 767 std::ostringstream message; << 255 G4Exception("G4GeometryManager::BuildOptimisations()", "GeomMgt0003", 768 message << "VoxelHeader allocation error << 256 FatalException, message); 769 << "Allocation of new VoxelHeade << 257 } 770 << " for volume " << tVol << 258 } 771 G4Exception("G4GeometryManager::BuildOpt << 259 else 772 FatalException, message); << 260 { 773 } << 261 // Don't create voxels for this node 774 } << 775 else << 776 { << 777 // Don't create voxels for this node << 778 #ifdef G4GEOMETRY_VOXELDEBUG 262 #ifdef G4GEOMETRY_VOXELDEBUG 779 G4cout << "** G4GeometryManager::BuildOpti << 263 G4cout << "**** G4GeometryManager::BuildOptimisations" << G4endl 780 << " Skipping logical volume na << 264 << " Skipping logical volume name = " << tVolume->GetName() 781 << G4endl; << 265 << G4endl; 782 #endif 266 #endif 783 } << 267 } 784 268 785 // Scan recursively the associated logical v << 269 // Scan recursively the associated logical volume tree 786 // << 270 // 787 tVolume = pVolume->GetLogicalVolume(); 271 tVolume = pVolume->GetLogicalVolume(); 788 if (tVolume->GetNoDaughters() != 0) << 272 if (tVolume->GetNoDaughters()) 789 { 273 { 790 BuildOptimisations(allOpts, tVolume->GetDa 274 BuildOptimisations(allOpts, tVolume->GetDaughter(0)); 791 } 275 } 792 } 276 } 793 277 794 // ******************************************* 278 // *************************************************************************** 795 // Removes all optimisation info. 279 // Removes all optimisation info. 796 // Loops over all logical volumes, deleting no << 280 // Loops over all logical volumes, deleting non-null voxels pointers, 797 // ******************************************* 281 // *************************************************************************** 798 // 282 // 799 void G4GeometryManager::DeleteOptimisations() 283 void G4GeometryManager::DeleteOptimisations() 800 { 284 { 801 G4LogicalVolume* tVolume = nullptr; << 285 G4LogicalVolume* tVolume = 0; 802 G4LogicalVolumeStore* Store = G4LogicalVolum 286 G4LogicalVolumeStore* Store = G4LogicalVolumeStore::GetInstance(); 803 for (auto & n : *Store) << 287 for (size_t n=0; n<Store->size(); n++) 804 { 288 { 805 tVolume=n; << 289 tVolume=(*Store)[n]; 806 delete tVolume->GetVoxelHeader(); 290 delete tVolume->GetVoxelHeader(); 807 tVolume->SetVoxelHeader(nullptr); << 291 tVolume->SetVoxelHeader(0); 808 } 292 } 809 } 293 } 810 294 811 // ******************************************* 295 // *************************************************************************** 812 // Removes optimisation info for the specified 296 // Removes optimisation info for the specified subtree. 813 // Scans recursively all daughter volumes, del 297 // Scans recursively all daughter volumes, deleting non-null voxels pointers. 814 // ******************************************* 298 // *************************************************************************** 815 // 299 // 816 void G4GeometryManager::DeleteOptimisations(G4 300 void G4GeometryManager::DeleteOptimisations(G4VPhysicalVolume* pVolume) 817 { 301 { 818 if (pVolume == nullptr) { return; } << 302 if (!pVolume) { return; } 819 303 820 // Retrieve the mother logical volume, if no 304 // Retrieve the mother logical volume, if not NULL, 821 // otherwise global deletion to world volume 305 // otherwise global deletion to world volume. 822 // 306 // 823 G4LogicalVolume* tVolume = pVolume->GetMothe 307 G4LogicalVolume* tVolume = pVolume->GetMotherLogical(); 824 if (tVolume == nullptr) { return DeleteOptim << 308 if (!tVolume) { return DeleteOptimisations(); } 825 delete tVolume->GetVoxelHeader(); 309 delete tVolume->GetVoxelHeader(); 826 tVolume->SetVoxelHeader(nullptr); << 310 tVolume->SetVoxelHeader(0); 827 311 828 // Scan recursively the associated logical v 312 // Scan recursively the associated logical volume tree 829 // 313 // 830 tVolume = pVolume->GetLogicalVolume(); 314 tVolume = pVolume->GetLogicalVolume(); 831 if (tVolume->GetNoDaughters() != 0) << 315 if (tVolume->GetNoDaughters()) 832 { 316 { 833 DeleteOptimisations(tVolume->GetDaughter(0 317 DeleteOptimisations(tVolume->GetDaughter(0)); 834 } 318 } 835 } 319 } 836 320 837 // ******************************************* 321 // *************************************************************************** 838 // Sets the maximum extent of the world volume 322 // Sets the maximum extent of the world volume. The operation is allowed only 839 // if NO solids have been created already. 323 // if NO solids have been created already. 840 // ******************************************* 324 // *************************************************************************** 841 // 325 // 842 void G4GeometryManager::SetWorldMaximumExtent( 326 void G4GeometryManager::SetWorldMaximumExtent(G4double extent) 843 { 327 { 844 if (!G4SolidStore::GetInstance()->empty()) << 328 if (G4SolidStore::GetInstance()->size()) 845 { 329 { 846 // Sanity check to assure that extent is 330 // Sanity check to assure that extent is fixed BEFORE creating 847 // any geometry object (solids in this ca 331 // any geometry object (solids in this case) 848 // 332 // 849 G4Exception("G4GeometryManager::SetMaximu 333 G4Exception("G4GeometryManager::SetMaximumExtent()", 850 "GeomMgt0003", FatalException 334 "GeomMgt0003", FatalException, 851 "Extent can be set only BEFOR 335 "Extent can be set only BEFORE creating any geometry object!"); 852 } 336 } 853 G4GeometryTolerance::GetInstance()->SetSurfa 337 G4GeometryTolerance::GetInstance()->SetSurfaceTolerance(extent); 854 } 338 } 855 339 856 // ******************************************* 340 // *************************************************************************** 857 // Reports statistics on voxel optimisation wh 341 // Reports statistics on voxel optimisation when closing geometry. 858 // ******************************************* 342 // *************************************************************************** 859 // 343 // 860 void 344 void 861 G4GeometryManager::ReportVoxelStats( std::vect 345 G4GeometryManager::ReportVoxelStats( std::vector<G4SmartVoxelStat> & stats, 862 G4double << 346 G4double totalCpuTime ) 863 std::ostr << 864 { 347 { 865 os << "------------------------------------- << 348 G4cout << "G4GeometryManager::ReportVoxelStats -- Voxel Statistics" 866 << G4endl; << 867 os << "G4GeometryManager::ReportVoxelStats - << 868 << G4endl << G4endl; 349 << G4endl << G4endl; 869 350 870 // 351 // 871 // Get total memory use 352 // Get total memory use 872 // 353 // 873 G4int i, nStat = (G4int)stats.size(); << 354 G4int i, nStat = stats.size(); 874 G4long totalMemory = 0; 355 G4long totalMemory = 0; 875 356 876 for( i=0; i<nStat; ++i ) { totalMemory += s << 357 for( i=0;i<nStat;++i ) { totalMemory += stats[i].GetMemoryUse(); } 877 358 878 os << " Total memory consumed for geometr << 359 G4cout << " Total memory consumed for geometry optimisation: " 879 << totalMemory/1024 << " kByte" << G4 360 << totalMemory/1024 << " kByte" << G4endl; 880 os << " Total CPU time elapsed for geomet << 361 G4cout << " Total CPU time elapsed for geometry optimisation: " 881 << std::setprecision(4) << totalCpuTi << 362 << std::setprecision(2) << totalCpuTime << " seconds" 882 << std::setprecision(6) << G4endl; 363 << std::setprecision(6) << G4endl; 883 364 884 // 365 // 885 // First list: sort by total CPU time 366 // First list: sort by total CPU time 886 // 367 // 887 std::sort( stats.begin(), stats.end(), << 368 std::sort( stats.begin(), stats.end(), G4SmartVoxelStat::ByCpu() ); 888 [](const G4SmartVoxelStat& a, const G4Smar << 889 { << 890 return a.GetTotalTime() > b.GetTotalTime() << 891 } ); << 892 369 893 const G4int maxPrint = 20; << 370 G4int nPrint = nStat > 10 ? 10 : nStat; 894 G4int nPrint = std::min ( nStat, maxPrint ); << 895 371 896 if (nPrint != 0) << 372 if (nPrint) 897 { 373 { 898 os << "\n Voxelisation: top CPU users:" << 374 G4cout << "\n Voxelisation: top CPU users:" << G4endl; 899 os << " Percent Total CPU System C << 375 G4cout << " Percent Total CPU System CPU Memory Volume\n" 900 << " ------- ---------- -------- << 376 << " ------- ---------- ---------- -------- ----------" 901 << G4endl; << 377 << G4endl; >> 378 // 12345678901.234567890123.234567890123.234567890123k . 902 } 379 } 903 380 904 for(i=0; i<nPrint; ++i) << 381 for(i=0;i<nPrint;++i) 905 { 382 { 906 G4double total = stats[i].GetTotalTime(); 383 G4double total = stats[i].GetTotalTime(); 907 G4double system = stats[i].GetSysTime(); 384 G4double system = stats[i].GetSysTime(); 908 G4double perc = 0.0; 385 G4double perc = 0.0; 909 386 910 if (system < 0) { system = 0.0; } 387 if (system < 0) { system = 0.0; } 911 if ((total < 0) || (totalCpuTime < perMill 388 if ((total < 0) || (totalCpuTime < perMillion)) 912 { total = 0; } 389 { total = 0; } 913 else 390 else 914 { perc = total*100/totalCpuTime; } 391 { perc = total*100/totalCpuTime; } 915 392 916 os << std::setprecision(2) << 393 G4cout << std::setprecision(2) 917 << std::setiosflags(std::ios::fixed 394 << std::setiosflags(std::ios::fixed|std::ios::right) 918 << std::setw(11) << perc 395 << std::setw(11) << perc 919 << std::setw(13) << total 396 << std::setw(13) << total 920 << std::setw(13) << system 397 << std::setw(13) << system 921 << std::setw(13) << (stats[i].GetMe 398 << std::setw(13) << (stats[i].GetMemoryUse()+512)/1024 922 << "k " << std::setiosflags(std::io 399 << "k " << std::setiosflags(std::ios::left) 923 << stats[i].GetVolume()->GetName() 400 << stats[i].GetVolume()->GetName() 924 << std::resetiosflags(std::ios::flo 401 << std::resetiosflags(std::ios::floatfield|std::ios::adjustfield) 925 << std::setprecision(6) 402 << std::setprecision(6) 926 << G4endl; 403 << G4endl; 927 } 404 } 928 405 929 // 406 // 930 // Second list: sort by memory use 407 // Second list: sort by memory use 931 // 408 // 932 std::sort( stats.begin(), stats.end(), << 409 std::sort( stats.begin(), stats.end(), G4SmartVoxelStat::ByMemory() ); 933 [](const G4SmartVoxelStat& a, const G4Smar << 934 { << 935 return a.GetMemoryUse() > b.GetMemoryUse() << 936 } ); << 937 410 938 if (nPrint != 0) << 411 if (nPrint) 939 { 412 { 940 os << "\n Voxelisation: top memory user << 413 G4cout << "\n Voxelisation: top memory users:" << G4endl; 941 os << " Percent Memory Heads << 414 G4cout << " Percent Memory Heads Nodes Pointers Total CPU Volume\n" 942 << " ------- -------- ------ << 415 << " ------- -------- ------ ------ -------- ---------- ----------" 943 << G4endl; << 416 << G4endl; >> 417 // 12345678901.2345678901k .23456789.23456789.2345678901.234567890123. . 944 } 418 } 945 419 946 for(i=0; i<nPrint; ++i) << 420 for(i=0;i<nPrint;++i) 947 { 421 { 948 G4long memory = stats[i].GetMemoryUse(); 422 G4long memory = stats[i].GetMemoryUse(); 949 G4double totTime = stats[i].GetTotalTime() 423 G4double totTime = stats[i].GetTotalTime(); 950 if (totTime < 0) { totTime = 0.0; } 424 if (totTime < 0) { totTime = 0.0; } 951 425 952 os << std::setprecision(2) << 426 G4cout << std::setprecision(2) 953 << std::setiosflags(std::ios::fixed|std << 427 << std::setiosflags(std::ios::fixed|std::ios::right) 954 << std::setw(11) << G4double(memory*100 << 428 << std::setw(11) << G4double(memory*100)/G4double(totalMemory) 955 << std::setw(11) << memory/1024 << "k " << 429 << std::setw(11) << memory/1024 << "k " 956 << std::setw( 9) << stats[i].GetNumberH << 430 << std::setw( 9) << stats[i].GetNumberHeads() 957 << std::setw( 9) << stats[i].GetNumberN << 431 << std::setw( 9) << stats[i].GetNumberNodes() 958 << std::setw(11) << stats[i].GetNumberP << 432 << std::setw(11) << stats[i].GetNumberPointers() 959 << std::setw(13) << totTime << " " << 433 << std::setw(13) << totTime << " " 960 << std::setiosflags(std::ios::left) << 434 << std::setiosflags(std::ios::left) 961 << stats[i].GetVolume()->GetName() << 435 << stats[i].GetVolume()->GetName() 962 << std::resetiosflags(std::ios::floatfi << 436 << std::resetiosflags(std::ios::floatfield|std::ios::adjustfield) 963 << std::setprecision(6) << 437 << std::setprecision(6) 964 << G4endl; << 438 << G4endl; 965 } 439 } 966 os << "------------------------------------- << 967 << G4endl << G4endl; << 968 } << 969 << 970 // ******************************************* << 971 // Check whether parallel optimisation was req << 972 // ******************************************* << 973 // << 974 G4bool G4GeometryManager::IsParallelOptimisati << 975 { << 976 return fOptimiseInParallelConfigured; << 977 } << 978 << 979 // ******************************************* << 980 // Report whether parallel optimisation is don << 981 // ******************************************* << 982 // << 983 G4bool G4GeometryManager::IsParallelOptimisati << 984 { << 985 return fParallelVoxelOptimisationFinished; << 986 } 440 } 987 441