Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4Scheduler Class Reference

#include <G4Scheduler.hh>

+ Inheritance diagram for G4Scheduler:

Public Member Functions

virtual G4bool Notify (G4ApplicationState requestedState)
 
virtual void RegisterModel (G4VITStepModel *, G4double)
 
void Initialize ()
 
void ForceReinitialization ()
 
G4bool IsInitialized ()
 
G4bool IsRunning ()
 
void Reset ()
 
void Process ()
 
void ClearList ()
 
void SetGun (G4ITGun *)
 
G4ITGunGetGun ()
 
void Stop ()
 
void Clear ()
 
void EndTracking ()
 
void SetEndTime (const G4double)
 
void SetTimeTolerance (G4double)
 
G4double GetTimeTolerance () const
 
void SetMaxZeroTimeAllowed (G4int)
 
G4int GetMaxZeroTimeAllowed () const
 
G4ITModelHandlerGetModelHandler ()
 
void SetTimeSteps (std::map< G4double, G4double > *)
 
void AddTimeStep (G4double, G4double)
 
void SetDefaultTimeStep (G4double)
 
G4double GetLimitingTimeStep () const
 
G4int GetNbSteps () const
 
void SetMaxNbSteps (G4int)
 
G4int GetMaxNbSteps () const
 
G4double GetStartTime () const
 
G4double GetEndTime () const
 
virtual G4double GetTimeStep () const
 
G4double GetPreviousTimeStep () const
 
G4double GetGlobalTime () const
 
void SetUserAction (G4UserTimeStepAction *)
 
G4UserTimeStepActionGetUserTimeStepAction () const
 
void UseDefaultTimeSteps (G4bool)
 
G4bool AreDefaultTimeStepsUsed ()
 
G4ITStepStatus GetStatus () const
 
void SetVerbose (G4int)
 
G4int GetVerbose () const
 
void WhyDoYouStop ()
 
void SetInteractivity (G4ITTrackingInteractivity *)
 
G4ITTrackingInteractivityGetInteractivity ()
 
virtual size_t GetNTracks ()
 
void GetCollisionType (G4String &interactionType)
 
void AddWatchedTime (G4double time)
 
G4double GetNextWatchedTime () const
 
void SetMaxTimeStep (G4double maxTimeStep)
 
G4double GetMaxTimeStep () const
 
G4VScavengerMaterialGetScavengerMaterial () const
 
void SetScavengerMaterial (std::unique_ptr< G4VScavengerMaterial > scavengerMaterial)
 
void ResetScavenger (bool)
 
virtual void Initialize ()
 
virtual void Reset ()
 
virtual void SetVerbose (int)
 
virtual void SetGun (G4ITGun *)
 
virtual void Process ()
 
virtual G4bool IsRunning ()
 
virtual G4ITModelHandlerGetModelHandler ()
 
virtual void RegisterModel (G4VITStepModel *, double)
 
virtual void SetEndTime (const double)
 
virtual void SetTimeTolerance (double)
 
virtual double GetTimeTolerance () const
 
virtual void SetMaxZeroTimeAllowed (int)
 
virtual int GetMaxZeroTimeAllowed () const
 
virtual void SetTimeSteps (std::map< double, double > *)
 
virtual void AddTimeStep (double, double)
 
virtual void SetDefaultTimeStep (double)
 
virtual double GetLimitingTimeStep () const
 
virtual G4int GetNbSteps () const
 
virtual void SetMaxNbSteps (G4int)
 
virtual G4int GetMaxNbSteps () const
 
virtual G4double GetStartTime () const
 
virtual G4double GetEndTime () const
 
virtual G4double GetTimeStep () const
 
virtual G4double GetPreviousTimeStep () const
 
virtual G4double GetGlobalTime () const
 
virtual void SetUserAction (G4UserTimeStepAction *)
 
virtual G4UserTimeStepActionGetUserTimeStepAction () const
 
virtual void SetInteractivity (G4ITTrackingInteractivity *)
 
virtual G4ITTrackingInteractivityGetInteractivity ()
 
- Public Member Functions inherited from G4VStateDependent
 G4VStateDependent (G4bool bottom=false)
 
virtual ~G4VStateDependent ()
 
G4bool operator== (const G4VStateDependent &right) const
 
G4bool operator!= (const G4VStateDependent &right) const
 
virtual G4bool Notify (G4ApplicationState requestedState)=0
 

Static Public Member Functions

static G4SchedulerInstance ()
 
static void DeleteInstance ()
 
- Static Public Member Functions inherited from G4VScheduler
static G4VSchedulerInstance ()
 

Protected Member Functions

virtual ~G4Scheduler ()
 
void DoProcess ()
 
void SynchronizeTracks ()
 
void Stepping ()
 
void FindUserPreDefinedTimeStep ()
 
G4bool CanICarryOn ()
 
void PrintWhyDoYouStop ()
 
- Protected Member Functions inherited from G4VScheduler
 G4VScheduler ()
 
virtual ~G4VScheduler ()
 

Detailed Description

G4Scheduler synchronizes (in time) track stepping

Definition at line 88 of file G4Scheduler.hh.

Constructor & Destructor Documentation

◆ ~G4Scheduler()

G4Scheduler::~G4Scheduler ( )
protectedvirtual

Definition at line 203 of file G4Scheduler.cc.

204{
205
206 if(fpMessenger) // is used as a flag to know whether the manager was cleared
207 {
208 Clear();
209 }
210
211 fgScheduler = 0;
212
213// if (fVerbose >= 1)
214// {
215// G4cout << "G4Scheduler is being deleted. Bye :) !" << G4endl;
216// }
217}
void Clear()
Definition: G4Scheduler.cc:219

Member Function Documentation

◆ AddTimeStep()

void G4Scheduler::AddTimeStep ( G4double  startingTime,
G4double  timeStep 
)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 318 of file G4Scheduler.hh.

319{
320 if (fpUserTimeSteps == 0)
321 {
322 fpUserTimeSteps = new std::map<G4double, G4double>();
323 fUsePreDefinedTimeSteps = true;
324 }
325
326 (*fpUserTimeSteps)[startingTime] = timeStep;
327}

◆ AddWatchedTime()

void G4Scheduler::AddWatchedTime ( G4double  time)
inline

Definition at line 177 of file G4Scheduler.hh.

178 {
179 fWatchedTimes.insert(time);
180 }

◆ AreDefaultTimeStepsUsed()

G4bool G4Scheduler::AreDefaultTimeStepsUsed ( )
inline

Definition at line 451 of file G4Scheduler.hh.

452{
453 return (fUseDefaultTimeSteps == false && fUsePreDefinedTimeSteps == false);
454}

Referenced by G4SchedulerMessenger::GetCurrentValue().

◆ CanICarryOn()

G4bool G4Scheduler::CanICarryOn ( )
protected

Definition at line 564 of file G4Scheduler.cc.

565{
566 return fGlobalTime < fEndTime && (fMaxSteps == -1 ? true : fNbSteps < fMaxSteps)
567 && fContinue == true;
568}

Referenced by SynchronizeTracks().

◆ Clear()

void G4Scheduler::Clear ( )

Definition at line 219 of file G4Scheduler.cc.

220{
221// if (fVerbose) G4cout << "*** G4Scheduler is being cleared ***" << G4endl;
222
223 if(fpMessenger)
224 {
225 delete fpMessenger;
226 fpMessenger = 0;
227 }
228 if(fpStepProcessor)
229 {
230 delete fpStepProcessor;
231 fpStepProcessor = 0;
232 }
233 if(fpModelProcessor)
234 {
235 delete fpModelProcessor;
236 fpModelProcessor = 0;
237 }
238
240 ClearList();
241 if(fpTrackingManager)
242 {
243 delete fpTrackingManager;
244 fpTrackingManager = 0;
245 }
246
247 if(fReactionSet)
248 {
249 delete fReactionSet;
250 fReactionSet = 0;
251 }
252
253 if(fpModelHandler)
254 {
255 delete fpModelHandler;
256 fpModelHandler = 0;
257 }
258
259 //* DEBUG
260 //* assert(G4StateManager::GetStateManager()->
261 //* DeregisterDependent(this) == true);
262
263}
static G4ITTypeManager * Instance()
Definition: G4ITType.cc:57
void ReleaseRessource()
Definition: G4ITType.cc:82
void ClearList()
Definition: G4Scheduler.cc:267

Referenced by Notify(), and ~G4Scheduler().

◆ ClearList()

void G4Scheduler::ClearList ( )

Definition at line 267 of file G4Scheduler.cc.

268{
269// if (fNbTracks == 0) return;
270
271 fTrackContainer.Clear();
272
274}
static void DeleteInstance()

Referenced by Clear(), and Process().

◆ DeleteInstance()

void G4Scheduler::DeleteInstance ( )
static

DeleteInstance should be used instead of the destructor

Definition at line 123 of file G4Scheduler.cc.

124{
125 if(fgScheduler)
126 {
127 delete fgScheduler;
128 }
129}

◆ DoProcess()

void G4Scheduler::DoProcess ( )
protected

Definition at line 616 of file G4Scheduler.cc.

618{
619 if(fpUserTimeStepAction) fpUserTimeStepAction->NewStage();
620
621#if defined (DEBUG_MEM) && defined (DEBUG_MEM_STEPPING)
622 MemStat mem_first, mem_second, mem_diff;
623#endif
624
625#if defined (DEBUG_MEM) && defined (DEBUG_MEM_STEPPING)
626 mem_first = MemoryUsage();
627#endif
628
629 while (fGlobalTime < fStopTime
630 && fTrackContainer.MainListsNOTEmpty()
631 && (fMaxSteps == -1 ? true : fNbSteps < fMaxSteps)
632 && fContinue == true)
633 {
634// G4cout << "Mainlist size : " << fTrackContainer.GetMainList()->size()
635// << G4endl;
636
637 Stepping();
638
639#if defined (DEBUG_MEM) && defined (DEBUG_MEM_STEPPING)
640 mem_second = MemoryUsage();
641 mem_diff = mem_second-mem_first;
642 G4cout << "\t || MEM || After step " << fNbSteps << ", diff is : "
643 << mem_diff << G4endl;
644#endif
645 }
646
648
649#if defined (DEBUG_MEM) && defined (DEBUG_MEM_STEPPING)
650 mem_second = MemoryUsage();
651 mem_diff = mem_second-mem_first;
652 G4cout << "\t || MEM || After stepping, diff is : " << mem_diff << G4endl;
653#endif
654
655#ifdef G4VERBOSE
656 if(fVerbose > 2)
657 G4cout << "*** G4Scheduler has finished processing a track list at time : "
658 << G4BestUnit(fGlobalTime, "Time") << G4endl;
659#endif
660}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
void Stepping()
Definition: G4Scheduler.cc:663
void PrintWhyDoYouStop()
Definition: G4Scheduler.cc:572
MemStat MemoryUsage()
Definition: G4MemStat.cc:55

Referenced by SynchronizeTracks().

◆ EndTracking()

void G4Scheduler::EndTracking ( )

Definition at line 1074 of file G4Scheduler.cc.

1075{
1076 if(fRunning)
1077 {
1078 G4ExceptionDescription exceptionDescription;
1079 exceptionDescription
1080 << "End tracking is called while G4Scheduler is still running."
1081 << G4endl;
1082
1083 G4Exception("G4Scheduler::EndTracking",
1084 "Scheduler017",
1086 exceptionDescription);
1087 }
1088
1089 fTrackContainer.MergeSecondariesWithMainList();
1090
1091 if (fTrackContainer.MainListsNOTEmpty())
1092 {
1093 G4TrackManyList* mainList = fTrackContainer.GetMainList();
1094 G4TrackManyList::iterator it = mainList->begin();
1095 G4TrackManyList::iterator end = mainList->end();
1096 for (; it != end; ++it)
1097 {
1098 fpTrackingManager->EndTrackingWOKill(*it);
1099 }
1100 }
1101
1102 if (fTrackContainer.SecondaryListsNOTEmpty()) // should be empty
1103 {
1104 G4TrackManyList* secondaries = fTrackContainer.GetSecondariesList();
1105 G4TrackManyList::iterator it = secondaries->begin();
1106 G4TrackManyList::iterator end = secondaries->end();
1107
1108 for (; it != end; ++it)
1109 {
1110 fpTrackingManager->EndTrackingWOKill(*it);
1111 }
1112 }
1113}
@ FatalErrorInArgument
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
G4TrackList * GetMainList(Key)
void MergeSecondariesWithMainList()
G4TrackManyList * GetSecondariesList()
bool SecondaryListsNOTEmpty()
void EndTrackingWOKill(G4Track *)

Referenced by Process().

◆ FindUserPreDefinedTimeStep()

void G4Scheduler::FindUserPreDefinedTimeStep ( )
protected

Definition at line 1019 of file G4Scheduler.cc.

1020{
1021
1022 if(fpUserTimeSteps == 0)
1023 {
1024 G4ExceptionDescription exceptionDescription;
1025 exceptionDescription
1026 << "You are asking to use user defined steps but you did not give any.";
1027 G4Exception("G4Scheduler::FindUserPreDefinedTimeStep",
1028 "Scheduler004",
1030 exceptionDescription);
1031 return; // makes coverity happy
1032 }
1033 map<G4double, G4double>::iterator fpUserTimeSteps_i =
1034 fpUserTimeSteps->upper_bound(fGlobalTime);
1035 map<G4double, G4double>::iterator fpUserTimeSteps_low = fpUserTimeSteps
1036 ->lower_bound(fGlobalTime);
1037
1038 // DEBUG
1039 // G4cout << "fGlobalTime : " << G4BestUnit(fGlobalTime,"Time") << G4endl;
1040 // G4cout << "fpUserTimeSteps_i : "
1041 // <<"<"<<G4BestUnit(fpUserTimeSteps_i->first,"Time")<<", "
1042 // << G4BestUnit(fpUserTimeSteps_i->second,"Time")<<">"
1043 // << "\t fpUserTimeSteps_low : "
1044 // <<"<"<<G4BestUnit(fpUserTimeSteps_low->first,"Time")<<", "
1045 // << G4BestUnit(fpUserTimeSteps_low->second,"Time")<<">"
1046 // << G4endl;
1047
1048 if(fpUserTimeSteps_i == fpUserTimeSteps->end())
1049 {
1050 fpUserTimeSteps_i--;
1051 }
1052 else if(fabs(fGlobalTime - fpUserTimeSteps_low->first) < fTimeTolerance)
1053 {
1054 // Case : fGlobalTime = X picosecond
1055 // and fpUserTimeSteps_low->first = X picosecond
1056 // but the precision is not good enough
1057 fpUserTimeSteps_i = fpUserTimeSteps_low;
1058 }
1059 else if(fpUserTimeSteps_i == fpUserTimeSteps_low)
1060 {
1061 // "Normal" cases
1062 fpUserTimeSteps_i--;
1063 }
1064 else
1065 {
1066 fpUserTimeSteps_i = fpUserTimeSteps_low;
1067 }
1068
1069 fDefinedMinTimeStep = fpUserTimeSteps_i->second;
1070}

◆ ForceReinitialization()

void G4Scheduler::ForceReinitialization ( )

Definition at line 1128 of file G4Scheduler.cc.

1129{
1130 fInitialized = false;
1131 Initialize();
1132}
void Initialize()
Definition: G4Scheduler.cc:284

◆ GetCollisionType()

void G4Scheduler::GetCollisionType ( G4String interactionType)

Definition at line 1159 of file G4Scheduler.cc.

1160{
1161 switch(fITStepStatus)
1162 {
1164 interactionType = "eInteractionWithMedium";
1165 break;
1167 interactionType = "eCollisionBetweenTracks";
1168 break;
1169 default:
1170 interactionType = "eCollisionBetweenTracks";
1171 break;
1172 }
1173}
@ eInteractionWithMedium
@ eCollisionBetweenTracks

Referenced by Stepping().

◆ GetEndTime()

G4double G4Scheduler::GetEndTime ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 349 of file G4Scheduler.hh.

350{
351 return fEndTime;
352}

Referenced by G4DNAIRT::G4DNAIRT(), G4SchedulerMessenger::GetCurrentValue(), and G4DNAIRT::Initialize().

◆ GetGlobalTime()

G4double G4Scheduler::GetGlobalTime ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 364 of file G4Scheduler.hh.

365{
366 return fGlobalTime;
367}

Referenced by G4ITTrackHolder::_PushTrack(), G4DNAIndependentReactionTimeStepper::CalculateStep(), G4DNAIRT::MakeReaction(), and G4DNAIRT::Sampling().

◆ GetGun()

G4ITGun * G4Scheduler::GetGun ( )
inline

Definition at line 436 of file G4Scheduler.hh.

437{
438 return fpGun;
439}

◆ GetInteractivity()

G4ITTrackingInteractivity * G4Scheduler::GetInteractivity ( )
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 426 of file G4Scheduler.hh.

427{
428 return fpTrackingInteractivity;
429}

◆ GetLimitingTimeStep()

G4double G4Scheduler::GetLimitingTimeStep ( ) const
virtual

Reimplemented from G4VScheduler.

Definition at line 953 of file G4Scheduler.cc.

954{
955 if (fpUserTimeSteps == 0) return fDefaultMinTimeStep;
956 if (fabs(fGlobalTime - fUserUpperTimeLimit) < fTimeTolerance)
957 return fDefinedMinTimeStep;
958
959 map<G4double, G4double>::const_iterator it_fpUserTimeSteps_i = fpUserTimeSteps
960 ->upper_bound(fGlobalTime);
961 map<G4double, G4double>::const_iterator it_fpUserTimeSteps_low = fpUserTimeSteps
962 ->lower_bound(fGlobalTime);
963
964 // DEBUG
965 // G4cout << "fGlobalTime : " << G4BestUnit(fGlobalTime,"Time")
966 // << G4endl;
967 // G4cout << "fpUserTimeSteps_i : "
968 // <<"<"<<G4BestUnit(it_fpUserTimeSteps->first,"Time")
969 // <<", "<< G4BestUnit(it_fpUserTimeSteps->second,"Time")<<">"
970 // << "\t fpUserTimeSteps_low : "
971 // <<"<"<<G4BestUnit(fpUserTimeSteps_low->first,"Time")<<", "*
972 // << G4BestUnit(fpUserTimeSteps_low->second,"Time")<<">"
973 // << G4endl;
974
975 if (it_fpUserTimeSteps_i == fpUserTimeSteps->end())
976 {
977 it_fpUserTimeSteps_i--;
978 fUserUpperTimeLimit = fStopTime;
979 }
980 else if (fabs(fGlobalTime - it_fpUserTimeSteps_low->first) < fTimeTolerance)
981 {
982 // Case : fGlobalTime = X picosecond
983 // and fpUserTimeSteps_low->first = X picosecond
984 // but the precision is not good enough
985 it_fpUserTimeSteps_i = it_fpUserTimeSteps_low;
986 map<G4double, G4double>::const_iterator tmp_it = it_fpUserTimeSteps_low;
987 ++tmp_it;
988 if (tmp_it == fpUserTimeSteps->end())
989 {
990 fUserUpperTimeLimit = fStopTime;
991 }
992 else
993 {
994 fUserUpperTimeLimit = tmp_it->first;
995 }
996 }
997 else if (it_fpUserTimeSteps_i == it_fpUserTimeSteps_low)
998 {
999 // "Normal" cases
1000 fUserUpperTimeLimit = it_fpUserTimeSteps_i->first;
1001// it_fpUserTimeSteps_i++;
1002// G4cout << "Global time = " << fGlobalTime << G4endl;
1003// G4cout << "Is begin = "
1004// << (it_fpUserTimeSteps_i == fpUserTimeSteps->begin())<< G4endl;
1005
1006 if(it_fpUserTimeSteps_i != fpUserTimeSteps->begin()) it_fpUserTimeSteps_i--;
1007 }
1008 else
1009 {
1010 fUserUpperTimeLimit = it_fpUserTimeSteps_i->first;
1011 it_fpUserTimeSteps_i = it_fpUserTimeSteps_low;
1012 }
1013
1014 return it_fpUserTimeSteps_i->second;
1015}

Referenced by Stepping().

◆ GetMaxNbSteps()

G4int G4Scheduler::GetMaxNbSteps ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 339 of file G4Scheduler.hh.

340{
341 return fMaxSteps;
342}

Referenced by G4SchedulerMessenger::GetCurrentValue().

◆ GetMaxTimeStep()

G4double G4Scheduler::GetMaxTimeStep ( ) const
inline

Definition at line 189 of file G4Scheduler.hh.

190 {
191 return fMaxTimeStep;
192 }

◆ GetMaxZeroTimeAllowed()

G4int G4Scheduler::GetMaxZeroTimeAllowed ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 396 of file G4Scheduler.hh.

397{
398 return fMaxNZeroTimeStepsAllowed;
399}

Referenced by G4SchedulerMessenger::GetCurrentValue().

◆ GetModelHandler()

G4ITModelHandler * G4Scheduler::GetModelHandler ( )
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 301 of file G4Scheduler.hh.

302{
303 return fpModelHandler;
304}

◆ GetNbSteps()

G4int G4Scheduler::GetNbSteps ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 329 of file G4Scheduler.hh.

330{
331 return fNbSteps;
332}

◆ GetNextWatchedTime()

G4double G4Scheduler::GetNextWatchedTime ( ) const

Definition at line 502 of file G4Scheduler.cc.

503{
504 std::set<G4double>::const_iterator up = fWatchedTimes.upper_bound(fGlobalTime);
505 if(up == fWatchedTimes.end()) return DBL_MAX;
506 return *up;
507}
#define DBL_MAX
Definition: templates.hh:62

Referenced by SynchronizeTracks().

◆ GetNTracks()

size_t G4Scheduler::GetNTracks ( )
virtual

Definition at line 1153 of file G4Scheduler.cc.

1154{
1155 return fTrackContainer.GetNTracks();
1156}
virtual size_t GetNTracks()

◆ GetPreviousTimeStep()

G4double G4Scheduler::GetPreviousTimeStep ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 411 of file G4Scheduler.hh.

412{
413 return fPreviousTimeStep;
414}

◆ GetScavengerMaterial()

G4VScavengerMaterial * G4Scheduler::GetScavengerMaterial ( ) const
inline

◆ GetStartTime()

G4double G4Scheduler::GetStartTime ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 344 of file G4Scheduler.hh.

345{
346 return fStartTime;
347}

Referenced by G4DNAIndependentReactionTimeStepper::CalculateStep(), G4DNAIRT::G4DNAIRT(), and G4DNAIRT::Initialize().

◆ GetStatus()

G4ITStepStatus G4Scheduler::GetStatus ( ) const
inline

Definition at line 416 of file G4Scheduler.hh.

417{
418 return fITStepStatus;
419}

◆ GetTimeStep()

G4double G4Scheduler::GetTimeStep ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 354 of file G4Scheduler.hh.

355{
356 return fTimeStep;
357}

◆ GetTimeTolerance()

G4double G4Scheduler::GetTimeTolerance ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 406 of file G4Scheduler.hh.

407{
408 return fTimeTolerance;
409}

Referenced by G4ITTrackHolder::_PushTrack(), and G4SchedulerMessenger::GetCurrentValue().

◆ GetUserTimeStepAction()

G4UserTimeStepAction * G4Scheduler::GetUserTimeStepAction ( ) const
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 375 of file G4Scheduler.hh.

376{
377 return fpUserTimeStepAction;
378}

◆ GetVerbose()

G4int G4Scheduler::GetVerbose ( ) const
inline

Definition at line 385 of file G4Scheduler.hh.

386{
387 return fVerbose;
388}

Referenced by G4SchedulerMessenger::GetCurrentValue().

◆ Initialize()

void G4Scheduler::Initialize ( )
virtual

Reimplemented from G4VScheduler.

Definition at line 284 of file G4Scheduler.cc.

285{
286 if(fpStepProcessor)
287 {
288 delete fpStepProcessor;
289 }
290 if(fpModelProcessor)
291 {
292 delete fpModelProcessor;
293 }
294 // if(fpMasterModelProcessor)
295 // {
296 // delete fpMasterModelProcessor;
297 // }
298
299 //______________________________________________________________
300
301 fpModelProcessor = new G4ITModelProcessor();
302 fpModelProcessor->SetModelHandler(fpModelHandler);
303 fpModelProcessor->SetTrackingManager(fpTrackingManager);
304
305 // fpMasterModelProcessor = new G4ITModelProcessor();
306 // fpMasterModelProcessor->SetModelHandler(fpModelHandler);
307 // fpModelProcessor = fpMasterModelProcessor;
308
309 //______________________________________________________________
310
311 fpStepProcessor = new G4ITStepProcessor();
312 fpStepProcessor->SetTrackingManager(fpTrackingManager);
313
314 fpTrackingManager->SetInteractivity(fpTrackingInteractivity);
315
316 // fpMasterStepProcessor = new G4ITStepProcessor();
317 // fpMasterStepProcessor->SetTrackingManager(fpTrackingManager);
318 // fpStepProcessor = fpMasterStepProcessor ;
319 //______________________________________________________________
320
321 if(fUsePreDefinedTimeSteps)
322 {
323 if(fpUserTimeSteps == 0) // Extra checking, is it necessary ?
324 {
325 G4ExceptionDescription exceptionDescription;
326 exceptionDescription
327 << "You are asking to use user defined steps but you did not give any.";
328 G4Exception("G4Scheduler::FindUserPreDefinedTimeStep",
329 "Scheduler004",
331 exceptionDescription);
332 return; // makes coverity happy
333 }
334 }
335
336// if (fComputeTimeStep)
337// {
338// if (fpModelProcessor == 0)
339// {
340// G4ExceptionDescription exceptionDescription;
341// exceptionDescription
342// << "There is no G4ITModelProcessor to handle IT reaction. ";
343// exceptionDescription
344// << "You probably did not initialize the G4Scheduler. ";
345// exceptionDescription
346// << "Just do G4Scheduler::Instance()->Initialize(); ";
347// exceptionDescription << " but only after initializing the run manager.";
348// G4Exception("G4Scheduler::CalculateMinStep", "ITScheduler005",
349// FatalErrorInArgument, exceptionDescription);
350// return; // makes coverity happy
351// }
352// }
353
354 fInitialized = true;
355}
void SetModelHandler(G4ITModelHandler *)
void SetTrackingManager(G4ITTrackingManager *trackingManager)
void SetTrackingManager(G4ITTrackingManager *trackMan)
void SetInteractivity(G4ITTrackingInteractivity *)

Referenced by ForceReinitialization(), G4DNAChemistryManager::InitializeThread(), Process(), and G4SchedulerMessenger::SetNewValue().

◆ Instance()

◆ IsInitialized()

G4bool G4Scheduler::IsInitialized ( )
inline

Definition at line 296 of file G4Scheduler.hh.

297{
298 return fInitialized;
299}

Referenced by G4SchedulerMessenger::GetCurrentValue().

◆ IsRunning()

G4bool G4Scheduler::IsRunning ( )
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 108 of file G4Scheduler.hh.

108{return fRunning;}

◆ Notify()

G4bool G4Scheduler::Notify ( G4ApplicationState  requestedState)
virtual

Implements G4VStateDependent.

Definition at line 108 of file G4Scheduler.cc.

109{
110 if(requestedState == G4State_Quit)
111 {
112 if(fVerbose >= 4)
113 {
114 G4cout << "G4Scheduler received G4State_Quit" << G4endl;
115 }
116 Clear();
117 //DeleteInstance();
118 }
119 return true;
120}
@ G4State_Quit

◆ PrintWhyDoYouStop()

void G4Scheduler::PrintWhyDoYouStop ( )
protected

Definition at line 572 of file G4Scheduler.cc.

573{
574#ifdef G4VERBOSE
575 if(fWhyDoYouStop)
576 {
577 G4cout << "G4Scheduler has reached a stage: it might be"
578 " a transition or the end"
579 << G4endl;
580
581 G4bool normalStop = false;
582
583 if(fGlobalTime >= fStopTime)
584 {
585 G4cout << "== G4Scheduler: I stop because I reached the stop time : "
586 << G4BestUnit(fStopTime,"Time") << " =="<< G4endl;
587 normalStop = true;
588 }
589 if(fTrackContainer.MainListsNOTEmpty() == false) // is empty
590 {
591 G4cout << "G4Scheduler: I stop because the current main list of tracks "
592 "is empty"
593 << G4endl;
594 normalStop = true;
595 }
596 if(fMaxSteps == -1 ? false : fNbSteps >= fMaxSteps)
597 {
598 G4cout << "G4Scheduler: I stop because I reached the maximum allowed "
599 "number of steps=" << fMaxSteps
600 << G4endl;
601 normalStop = true;
602 }
603 if(fContinue && normalStop == false)
604 {
605 G4cout << "G4Scheduler: It might be that I stop because "
606 "I have been told so. You may check "
607 "member fContinue and usage of the method G4Scheduler::Stop()."
608 << G4endl;
609 }
610 }
611#endif
612}
bool G4bool
Definition: G4Types.hh:86

Referenced by DoProcess().

◆ Process()

void G4Scheduler::Process ( )
virtual

Reimplemented from G4VScheduler.

Definition at line 379 of file G4Scheduler.cc.

380{
381
382#ifdef G4VERBOSE
383 if(fVerbose)
384 {
385 G4cout << "*** G4Scheduler starts processing " << G4endl;
386 if(fVerbose > 2)
387 G4cout << "___________________________________________"
388 "___________________________" << G4endl;
389 }
390#endif
391
392 if (fInitialized == false) Initialize();
393
394 // fpTrackingManager->Initialize();
395 fpModelProcessor->Initialize();
396 fpStepProcessor->Initialize();
397
398 // TODO
399 // fpMasterModelProcessor->Initialize();
400 // fpMasterStepProcessor->Initialize();
401
402 if (fpGun) fpGun->DefineTracks();
403
404 if (fpTrackingInteractivity) fpTrackingInteractivity->Initialize();
405
406 // ___________________
407 fRunning = true;
408 Reset();
409
410 //hoang 7/12/2020
411
412 if(fResetScavenger) {
413 if (fpUserScavenger != nullptr) {
414 fpUserScavenger->Reset();
415 }
416 }
417
418 if (fpUserTimeStepAction)
419 {
420 fpUserTimeStepAction->StartProcessing();
421 }
422
423#ifdef G4VERBOSE
424 G4bool trackFound = false;
425 G4IosFlagsSaver iosfs(G4cout);
426 G4cout.precision(5);
427#endif
428
429 //===========================================================================
430 // By default, before the G4Scheduler is launched, the tracks are pushed to
431 // the delayed lists
432 //===========================================================================
433
434 if(fTrackContainer.DelayListsNOTEmpty())
435 {
436 fStartTime = fTrackContainer.GetNextTime();
437#ifdef G4VERBOSE
438 trackFound = true;
439 G4Timer localtimer;
440 if(fVerbose>1)
441 {
442 localtimer.Start();
443 }
444#endif
446#ifdef G4VERBOSE
447 if(fVerbose>1)
448 {
449 localtimer.Stop();
450 G4cout << "G4Scheduler: process time= "<< localtimer << G4endl;
451 }
452#endif
453 }
454
455// //---------------------------------
456// // TODO: This could be removed ?
457// if(fTrackContainer.MainListsNOTEmpty()
458// && (fMaxSteps == -1 ? true : fNbSteps < fMaxSteps)
459// && fGlobalTime < fEndTime
460// && fContinue == true)
461//{
462//#ifdef G4VERBOSE
463// trackFound = true;
464//#endif
465// DoProcess();
466//}
467// //---------------------------------
468
469#ifdef G4VERBOSE
470 if(fVerbose)
471 {
472 if(trackFound)
473 {
474 G4cout << "*** G4Scheduler ends at time : "
475 << G4BestUnit(fGlobalTime,"Time") << G4endl;
476 G4cout << "___________________________________" << G4endl;
477 }
478 else
479 {
480 G4cout << "*** G4Scheduler did not start because no "
481 "track was found to be processed"<< G4endl;
482 G4cout << "___________________________________" << G4endl;
483 }
484 }
485#endif
486
487 fRunning = false;
488
489 if (fpUserTimeStepAction) fpUserTimeStepAction->EndProcessing();
490
491 // ___________________
492 EndTracking();
493 ClearList();
494
495 Reset();
496
497 if (fpTrackingInteractivity) fpTrackingInteractivity->Finalize();
498}
virtual void DefineTracks()
Definition: G4ITGun.hh:56
virtual void Initialize()
void SynchronizeTracks()
Definition: G4Scheduler.cc:511
void EndTracking()
void Reset()
Definition: G4Scheduler.cc:359
void Stop()
void Start()
virtual void StartProcessing()

Referenced by G4DNAChemistryManager::Run(), and G4SchedulerMessenger::SetNewValue().

◆ RegisterModel()

void G4Scheduler::RegisterModel ( G4VITStepModel model,
G4double  time 
)
virtual

Reimplemented from G4VScheduler.

Definition at line 277 of file G4Scheduler.cc.

278{
279 fpModelHandler->RegisterModel(model, time);
280}
void RegisterModel(G4VITStepModel *pModel, G4double globalTime)

◆ Reset()

void G4Scheduler::Reset ( )
virtual

Reimplemented from G4VScheduler.

Definition at line 359 of file G4Scheduler.cc.

360{
361 fStartTime = 0;
362 fUserUpperTimeLimit = -1;
363 fTimeStep = DBL_MAX;
364 fTSTimeStep = DBL_MAX;
365 fILTimeStep = DBL_MAX;
366 fPreviousTimeStep = DBL_MAX;
367 fGlobalTime = -1;
368 fInteractionStep = true;
369 fITStepStatus = eUndefined;
370 fZeroTimeCount = 0;
371
372 fNbSteps = 0;
373 fContinue = true;
374 // fReactingTracks.clear();
375 fReactionSet->CleanAllReaction();
376}
@ eUndefined
void CleanAllReaction()

Referenced by Process().

◆ ResetScavenger()

void G4Scheduler::ResetScavenger ( bool  value)
inline

Definition at line 456 of file G4Scheduler.hh.

457{
458 fResetScavenger = value;
459}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ SetDefaultTimeStep()

void G4Scheduler::SetDefaultTimeStep ( G4double  timeStep)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 359 of file G4Scheduler.hh.

360{
361 fDefaultMinTimeStep = timeStep;
362}

◆ SetEndTime()

void G4Scheduler::SetEndTime ( const G4double  __endtime)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 306 of file G4Scheduler.hh.

307{
308 fEndTime = __endtime;
309}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ SetGun()

void G4Scheduler::SetGun ( G4ITGun gun)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 431 of file G4Scheduler.hh.

432{
433 fpGun = gun;
434}

Referenced by G4DNAChemistryManager::SetGun().

◆ SetInteractivity()

void G4Scheduler::SetInteractivity ( G4ITTrackingInteractivity interactivity)
virtual

Reimplemented from G4VScheduler.

Definition at line 1116 of file G4Scheduler.cc.

1117{
1118 fpTrackingInteractivity = interactivity;
1119 if(fpTrackingManager)
1120 {
1121 fpTrackingManager->SetInteractivity(fpTrackingInteractivity);
1122 }
1123
1124 //G4MIWorkspace::GetWorldWorkspace()->SetTrackingInteractivity(interactivity);
1125}

◆ SetMaxNbSteps()

void G4Scheduler::SetMaxNbSteps ( G4int  maxSteps)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 334 of file G4Scheduler.hh.

335{
336 fMaxSteps = maxSteps;
337}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ SetMaxTimeStep()

void G4Scheduler::SetMaxTimeStep ( G4double  maxTimeStep)
inline

Definition at line 184 of file G4Scheduler.hh.

185 {
186 fMaxTimeStep = maxTimeStep;
187 }

◆ SetMaxZeroTimeAllowed()

void G4Scheduler::SetMaxZeroTimeAllowed ( G4int  maxTimeStepAllowed)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 391 of file G4Scheduler.hh.

392{
393 fMaxNZeroTimeStepsAllowed = maxTimeStepAllowed;
394}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ SetScavengerMaterial()

void G4Scheduler::SetScavengerMaterial ( std::unique_ptr< G4VScavengerMaterial scavengerMaterial)
inline

Definition at line 198 of file G4Scheduler.hh.

199 {
200 fpUserScavenger = std::move(scavengerMaterial);
201 }

◆ SetTimeSteps()

void G4Scheduler::SetTimeSteps ( std::map< G4double, G4double > *  steps)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 312 of file G4Scheduler.hh.

313{
314 fUsePreDefinedTimeSteps = true;
315 fpUserTimeSteps = steps;
316}

◆ SetTimeTolerance()

void G4Scheduler::SetTimeTolerance ( G4double  time)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 401 of file G4Scheduler.hh.

402{
403 fTimeTolerance = time;
404}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ SetUserAction()

void G4Scheduler::SetUserAction ( G4UserTimeStepAction userITAction)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 370 of file G4Scheduler.hh.

371{
372 fpUserTimeStepAction = userITAction;
373}

◆ SetVerbose()

void G4Scheduler::SetVerbose ( G4int  verbose)
inlinevirtual

Reimplemented from G4VScheduler.

Definition at line 380 of file G4Scheduler.hh.

381{
382 fVerbose = verbose;
383}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ Stepping()

void G4Scheduler::Stepping ( )
protected

Definition at line 663 of file G4Scheduler.cc.

664{
665 fTimeStep = fMaxTimeStep;
666
667 fTSTimeStep = DBL_MAX;
668 fILTimeStep = DBL_MAX;
669
670 fInteractionStep = false;
671 fReachedUserTimeLimit = false;
672
673 fITStepStatus = eUndefined;
674
675 // Start of step
676#ifdef G4VERBOSE
677 if (fVerbose > 2)
678 {
679#ifdef USE_COLOR
680 G4cout << LIGHT_RED;
681#endif
682 G4cout << "*** Start Of Step N°" << fNbSteps + 1 << " ***" << G4endl;
683 G4cout << "Current Global time : " << G4BestUnit(fGlobalTime, "Time")
684 <<G4endl;
685#ifdef USE_COLOR
687#endif
688 }
689#endif
690
691#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
692 MemStat mem_first, mem_second, mem_diff;
693#endif
694
695#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
696 mem_first = MemoryUsage();
697#endif
698
699 fDefinedMinTimeStep = GetLimitingTimeStep();
700
701 if (fUsePreDefinedTimeSteps)
702 {
703#ifdef G4VERBOSE
704 if (fVerbose > 2)
705 {
706#ifdef USE_COLOR
707 G4cout << LIGHT_RED;
708#endif
709 G4cout << "*** At time : " << G4BestUnit(fGlobalTime, "Time")
710 << " the chosen user time step is : "
711 << G4BestUnit(fDefinedMinTimeStep, "Time") << " ***" << G4endl;
712#ifdef USE_COLOR
714#endif
715 }
716#endif
717 }
718
719 if (fpModelProcessor->GetComputeTimeStep()) // fComputeTimeStep)
720 {
721 fTSTimeStep = fpModelProcessor->CalculateMinTimeStep(fGlobalTime,
722 fDefinedMinTimeStep);
723 // => at least N (N = nb of tracks) loops
724 }
725 else if(fUseDefaultTimeSteps)
726 {
727 fTSTimeStep = fDefinedMinTimeStep;
728 }
729
730#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
731 mem_second = MemoryUsage();
732 mem_diff = mem_second-mem_first;
733 G4cout << "|| MEM || After computing TS, diff is : " << mem_diff << G4endl;
734#endif
735
736#ifdef G4VERBOSE
737 if (fVerbose > 2)
738 {
739#ifdef USE_COLOR
740 G4cout << LIGHT_RED;
741#endif
742 G4cout << "*** Time stepper returned : " << G4BestUnit(fTSTimeStep, "Time")
743 << " ***" << G4endl;
744#ifdef USE_COLOR
746#endif
747 }
748#endif
749
750#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
751 mem_first = MemoryUsage();
752#endif
753
754 // Call IL even if fTSTimeStep == 0
755 // if fILTimeStep == 0 give the priority to DoIt processes
756
757 fILTimeStep = fpStepProcessor->ComputeInteractionLength(fPreviousTimeStep);
758 // => at least N loops
759 // All process returns the physical step of interaction
760 // The transportation process calculates the corresponding
761 // time step
762
763#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
764 mem_second = MemoryUsage();
765 mem_diff = mem_second-mem_first;
766 G4cout << "|| MEM || After IL, diff is : " << mem_diff << G4endl;
767#endif
768
769#ifdef G4VERBOSE
770 if (fVerbose > 2)
771 {
772#ifdef USE_COLOR
773 G4cout << LIGHT_RED;
774#endif
775 G4cout << "*** The minimum time returned by the processes is : "
776 << G4BestUnit(fILTimeStep, "Time") << " ***" << G4endl;
777#ifdef USE_COLOR
779#endif
780 }
781#endif
782
783#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
784 mem_first = MemoryUsage();
785#endif
786
787 if (fILTimeStep <= fTSTimeStep)
788 // Give the priority to the IL
789 {
790 fInteractionStep = true;
791 fReactionSet->CleanAllReaction();
792 fTimeStep = fILTimeStep;
793 fITStepStatus = eInteractionWithMedium;
794 fpStepProcessor->PrepareLeadingTracks();
795 }
796 else
797 {
798 fInteractionStep = false;
799 fpStepProcessor->ResetLeadingTracks();
800 fTimeStep = fTSTimeStep;
801 fITStepStatus = eCollisionBetweenTracks;
802 }
803
804 if (fGlobalTime + fTimeStep > fStopTime)
805 // This check is done at every time step
806 {
807 fTimeStep = fStopTime - fGlobalTime;
808 fITStepStatus = eInteractionWithMedium; // ie: transportation
809 fInteractionStep = true;
810 fReactionSet->CleanAllReaction();
811 fpStepProcessor->ResetLeadingTracks();
812 }
813
814 if (fTimeStep == 0) // < fTimeTolerance)
815 {
816 ++fZeroTimeCount;
817 if (fZeroTimeCount >= fMaxNZeroTimeStepsAllowed)
818 {
819 G4ExceptionDescription exceptionDescription;
820
821 exceptionDescription << "Too many zero time steps were detected. ";
822 exceptionDescription << "The simulation is probably stuck. ";
823 exceptionDescription
824 << "The maximum number of zero time steps is currently : "
825 << fMaxNZeroTimeStepsAllowed;
826 exceptionDescription << ".";
827
828 G4Exception("G4Scheduler::Stepping",
829 "SchedulerNullTimeSteps",
831 exceptionDescription);
832 }
833 }
834 else
835 {
836 fZeroTimeCount = 0;
837 }
838
839 fReachedUserTimeLimit =
840 ((fTimeStep <= fDefinedMinTimeStep) || ((fTimeStep > fDefinedMinTimeStep)
841 && fabs(fTimeStep - fDefinedMinTimeStep) < fTimeTolerance)) ?
842 true : false;
843
844 if (fpUserTimeStepAction) fpUserTimeStepAction->UserPreTimeStepAction();
845 // TODO: pre/post
846
847#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
848 mem_second = MemoryUsage();
849 mem_diff = mem_second-mem_first;
850 G4cout << "|| MEM || After LeadingTracks and UserPreTimeStepAction: "
851 << mem_diff << G4endl;
852#endif
853
854#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
855 mem_first = MemoryUsage();
856#endif
857
858
859 fGlobalTime += fTimeStep;
860
861 // if fTSTimeStep > 0 => still need to call the transportation process
862 // if fILTimeStep < fTSTimeStep => call only DoIt processes, no reactions
863 // if fILTimeStep == fTSTimeStep => give the priority to the DoIt processes
864 if (fTSTimeStep > 0 || fILTimeStep <= fTSTimeStep)
865 {
866 // G4cout << "Will call DoIT" << G4endl;
867 fpStepProcessor->DoIt(fTimeStep);
868 // fTrackContainer.MergeSecondariesWithMainList();
869 // fTrackContainer.KillTracks(); // remove ?
870 }
871 // else
872 // {
873 // G4cout << "fTSTimeStep : " << fTSTimeStep << G4endl;
874 // G4cout << "fILTimeStep : " << fILTimeStep << G4endl;
875 // }
876
877#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
878 mem_second = MemoryUsage();
879 mem_diff = mem_second-mem_first;
880 G4cout << "|| MEM || After DoIT, diff is : " << mem_diff << G4endl;
881#endif
882
883#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
884 mem_first = MemoryUsage();
885#endif
886
887 fpModelProcessor->ComputeTrackReaction(fITStepStatus,
888 fGlobalTime,
889 fTimeStep,
890 fPreviousTimeStep,
891 fReachedUserTimeLimit,
892 fTimeTolerance,
893 fpUserTimeStepAction,
894 fVerbose);
895
896 ++fNbSteps;
897
898 if (fpUserTimeStepAction)
899 {
900 fpUserTimeStepAction->UserPostTimeStepAction();
901 }
902
903 fPreviousTimeStep = fTimeStep;
904
905#if defined (DEBUG_MEM) && defined (DEBUG_MEM_DETAILED_STEPPING)
906 mem_second = MemoryUsage();
907 mem_diff = mem_second-mem_first;
908 G4cout << "|| MEM || After computing reactions + UserPostTimeStepAction, "
909 "diff is : " << mem_diff << G4endl;
910#endif
911
912 // End of step
913#ifdef G4VERBOSE
914 if (fVerbose >= 2)
915 {
916#ifdef USE_COLOR
917 G4cout << LIGHT_RED;
918#endif
919
920 G4String interactionType;
921 GetCollisionType(interactionType);
922
923 std::stringstream finalOutput;
924
925 finalOutput << "*** End of step N°" << fNbSteps
926 << "\t T_i= " << G4BestUnit(fGlobalTime-fTimeStep, "Time")
927 << "\t dt= " << G4BestUnit(fTimeStep, "Time")
928 << "\t T_f= " << G4BestUnit(fGlobalTime, "Time")
929 << "\t " << interactionType
930 << G4endl;
931
932 if(fVerbose>2)
933 {
934 if(fReachedUserTimeLimit)
935 {
936 finalOutput << "It has also reached the user time limit" << G4endl;
937 }
938 finalOutput << "_______________________________________________________________"
939 "_______"<< G4endl;
940 }
941
942 G4cout << finalOutput.str();
943
944#ifdef USE_COLOR
946#endif
947 }
948#endif
949
950}
#define LIGHT_RED
Definition: G4Scheduler.cc:83
#define RESET_COLOR
Definition: G4Scheduler.cc:86
G4double CalculateMinTimeStep(G4double currentGlobalTime, G4double definedMinTimeStep)
void ComputeTrackReaction(G4ITStepStatus fITStepStatus, G4double fGlobalTime, G4double currentTimeStep, G4double previousTimeStep, G4bool reachedUserTimeLimit, G4double fTimeTolerance, G4UserTimeStepAction *fpUserTimeStepAction, G4int fVerbose)
bool GetComputeTimeStep() const
G4double ComputeInteractionLength(double previousTimeStep)
void DoIt(double timeStep)
G4double GetLimitingTimeStep() const
Definition: G4Scheduler.cc:953
void GetCollisionType(G4String &interactionType)
virtual void UserPostTimeStepAction()
virtual void UserPreTimeStepAction()

Referenced by DoProcess().

◆ Stop()

void G4Scheduler::Stop ( )
inline

Definition at line 421 of file G4Scheduler.hh.

422{
423 fContinue = false;
424}

Referenced by G4DNAIRTMoleculeEncounterStepper::CalculateMinTimeStep().

◆ SynchronizeTracks()

void G4Scheduler::SynchronizeTracks ( )
protected

Definition at line 511 of file G4Scheduler.cc.

512{
513// if(fTrackContainer.WaitingListsNOTEmpty())
514// {
515// G4ExceptionDescription exceptionDescription;
516// exceptionDescription
517// << "There is a waiting track list (fpWaitingList != 0).";
518// exceptionDescription
519// << " When G4Scheduler::SynchronizeTracks() is called, ";
520// exceptionDescription
521// << "no more tracks should remain in the fpWaitingList.";
522// G4Exception("G4Scheduler::SynchronizeTracks", "ITScheduler002",
523// FatalErrorInArgument, exceptionDescription);
524// }
525
526 // Backup main list and time feature
527 // Reminder : the main list here, should
528 // have the biggest global time
529 // fTrackContainer.MoveMainToWaitingList();
530 // TODO: not yet supported
531
532 fTmpGlobalTime = fGlobalTime;
533 //fTmpEndTime = fEndTime;
534
535 fGlobalTime = fTrackContainer.GetNextTime();
536 G4double tmpGlobalTime = fGlobalTime;
537
538 G4double nextWatchedTime = -1;
539 G4bool carryOn = true;
540
541 while(fTrackContainer.MergeNextTimeToMainList(tmpGlobalTime) && carryOn)
542 {
543// assert(tmpGlobalTime == fGlobalTime);
544 fStopTime = min(fTrackContainer.GetNextTime(), fEndTime);
545 while((nextWatchedTime = GetNextWatchedTime()) < fTrackContainer.GetNextTime()
546 && (carryOn = CanICarryOn()))
547 {
548 fStopTime = min(nextWatchedTime, fEndTime);
549 DoProcess();
550 }
551
552 carryOn = CanICarryOn();
553
554 if(nextWatchedTime > fEndTime && carryOn)
555 {
556 fStopTime = min(fTrackContainer.GetNextTime(), fEndTime);
557 DoProcess();
558 }
559 }
560}
double G4double
Definition: G4Types.hh:83
bool MergeNextTimeToMainList(double &time)
G4bool CanICarryOn()
Definition: G4Scheduler.cc:564
G4double GetNextWatchedTime() const
Definition: G4Scheduler.cc:502
void DoProcess()
Definition: G4Scheduler.cc:616
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

Referenced by Process().

◆ UseDefaultTimeSteps()

void G4Scheduler::UseDefaultTimeSteps ( G4bool  flag)
inline

Definition at line 446 of file G4Scheduler.hh.

447{
448 fUseDefaultTimeSteps = flag;
449}

Referenced by G4SchedulerMessenger::SetNewValue().

◆ WhyDoYouStop()

void G4Scheduler::WhyDoYouStop ( )
inline

Definition at line 441 of file G4Scheduler.hh.

442{
443 fWhyDoYouStop = true;
444}

Referenced by G4SchedulerMessenger::SetNewValue().


The documentation for this class was generated from the following files: