Geant4-11
Data Structures | Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | Friends
G4MoleculeCounter Class Reference

#include <G4MoleculeCounter.hh>

Inheritance diagram for G4MoleculeCounter:
G4VMoleculeCounter

Data Structures

struct  Search
 

Public Types

using CounterMapType = std::map< Reactant *, NbMoleculeAgainstTime >
 
using Reactant = const G4MolecularConfiguration
 
using ReactantList = std::vector< Reactant * >
 
using RecordedMolecules = std::unique_ptr< ReactantList >
 

Public Member Functions

void CheckTimeForConsistency (G4bool flag)
 
void DontRegister (const G4MoleculeDefinition *) override
 
void Dump ()
 
const NbMoleculeAgainstTimeGetNbMoleculeAgainstTime (Reactant *molecule)
 
int GetNMoleculesAtTime (Reactant *molecule, double time)
 
RecordedMolecules GetRecordedMolecules ()
 
RecordedTimes GetRecordedTimes ()
 
G4int GetVerbose ()
 
void Initialize () override
 
bool IsRegistered (const G4MoleculeDefinition *) override
 
G4bool IsTimeCheckedForConsistency () const
 
void RegisterAll () override
 
void ResetCounter () override
 
void SetVerbose (G4int)
 

Static Public Member Functions

static void DeleteInstance ()
 
static void InitializeInstance ()
 
static G4MoleculeCounterInstance ()
 
static G4bool InUse ()
 
static void SetInstance (G4VMoleculeCounter *)
 
static void SetTimeSlice (double)
 
static void Use (G4bool flag=true)
 

Protected Member Functions

void AddAMoleculeAtTime (Reactant *, G4double time, const G4ThreeVector *position=nullptr, int number=1) override
 
 G4MoleculeCounter ()
 
void RemoveAMoleculeAtTime (Reactant *, G4double time, const G4ThreeVector *position=nullptr, int number=1) override
 
G4bool SearchTimeMap (Reactant *molecule)
 
int SearchUpperBoundTime (double time, bool sameTypeOfMolecule)
 
 ~G4MoleculeCounter () override
 

Protected Attributes

G4bool fCheckTimeIsConsistentWithScheduler
 
CounterMapType fCounterMap
 
std::map< const G4MoleculeDefinition *, G4boolfDontRegister
 
std::unique_ptr< SearchfpLastSearch
 
G4int fVerbose
 

Static Protected Attributes

static G4ThreadLocal G4VMoleculeCounterfpInstance = nullptr
 
static G4bool fUse = false
 

Friends

class G4Molecule
 
class G4VMoleculeCounter
 

Detailed Description

Definition at line 69 of file G4MoleculeCounter.hh.

Member Typedef Documentation

◆ CounterMapType

Definition at line 74 of file G4MoleculeCounter.hh.

◆ Reactant

Definition at line 58 of file G4VMoleculeCounter.hh.

◆ ReactantList

Definition at line 73 of file G4MoleculeCounter.hh.

◆ RecordedMolecules

Definition at line 75 of file G4MoleculeCounter.hh.

Constructor & Destructor Documentation

◆ G4MoleculeCounter()

G4MoleculeCounter::G4MoleculeCounter ( )
protected

Definition at line 71 of file G4MoleculeCounter.cc.

References fCheckTimeIsConsistentWithScheduler, and fVerbose.

Referenced by Instance().

◆ ~G4MoleculeCounter()

G4MoleculeCounter::~G4MoleculeCounter ( )
overrideprotecteddefault

Member Function Documentation

◆ AddAMoleculeAtTime()

void G4MoleculeCounter::AddAMoleculeAtTime ( Reactant molecule,
G4double  time,
const G4ThreeVector position = nullptr,
int  number = 1 
)
overrideprotectedvirtual

Implements G4VMoleculeCounter.

Definition at line 208 of file G4MoleculeCounter.cc.

212{
213 if (fDontRegister[molecule->GetDefinition()])
214 {
215 return;
216 }
217
218 if (fVerbose)
219 {
220 G4cout << "G4MoleculeCounter::AddAMoleculeAtTime : " << molecule->GetName()
221 << " at time : " << G4BestUnit(time, "Time") << G4endl;
222 }
223
224 auto counterMap_i = fCounterMap.find(molecule);
225
226 if (counterMap_i == fCounterMap.end())
227 {
228 fCounterMap[molecule][time] = number;
229 }
230 else if (counterMap_i->second.empty())
231 {
232 counterMap_i->second[time] = number;
233 }
234 else
235 {
236 NbMoleculeAgainstTime::reverse_iterator end = counterMap_i->second.rbegin();
237
238 if (end->first <= time ||
239 fabs(end->first - time) <= G4::MoleculeCounter::TimePrecision::fPrecision)
240 // Case 1 = new time comes after last recorded data
241 // Case 2 = new time is about the same as the last recorded one
242 {
243 double newValue = end->second + number;
244 counterMap_i->second[time] = newValue;
245 }
246 else
247 {
248 // if(fabs(time - G4Scheduler::Instance()->GetGlobalTime()) >
249 // G4Scheduler::Instance()->GetTimeTolerance())
250 {
252 errMsg << "Time of species "
253 << molecule->GetName() << " is "
254 << G4BestUnit(time, "Time") << " while "
255 << " global time is "
256 << G4BestUnit(G4Scheduler::Instance()->GetGlobalTime(), "Time")
257 << G4endl;
258 G4Exception("G4MoleculeCounter::AddAMoleculeAtTime",
259 "TIME_DONT_MATCH",
260 FatalException, errMsg);
261 }
262 }
263 }
264}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
#define G4BestUnit(a, b)
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
std::map< const G4MoleculeDefinition *, G4bool > fDontRegister
CounterMapType fCounterMap
static G4Scheduler * Instance()
Definition: G4Scheduler.cc:101
static G4ThreadLocal double fPrecision

References FatalException, fCounterMap, fDontRegister, G4::MoleculeCounter::TimePrecision::fPrecision, fVerbose, G4BestUnit, G4cout, G4endl, G4Exception(), G4MolecularConfiguration::GetDefinition(), G4MolecularConfiguration::GetName(), and G4Scheduler::Instance().

◆ CheckTimeForConsistency()

void G4MoleculeCounter::CheckTimeForConsistency ( G4bool  flag)

Definition at line 509 of file G4MoleculeCounter.cc.

510{
512}

References fCheckTimeIsConsistentWithScheduler.

◆ DeleteInstance()

void G4VMoleculeCounter::DeleteInstance ( )
staticinherited

Definition at line 78 of file G4VMoleculeCounter.cc.

79{
80 if (fpInstance != nullptr)
81 {
82 delete fpInstance;
83 fpInstance = nullptr;
84 }
85}
static G4ThreadLocal G4VMoleculeCounter * fpInstance

References G4VMoleculeCounter::fpInstance.

Referenced by G4DNAChemistryManager::Clear().

◆ DontRegister()

void G4MoleculeCounter::DontRegister ( const G4MoleculeDefinition molDef)
overridevirtual

Reimplemented from G4VMoleculeCounter.

Definition at line 481 of file G4MoleculeCounter.cc.

482{
483 fDontRegister[molDef] = true;
484}

References fDontRegister.

◆ Dump()

void G4MoleculeCounter::Dump ( )

Definition at line 430 of file G4MoleculeCounter.cc.

431{
432 for (auto it : fCounterMap)
433 {
434 auto pReactant = it.first;
435
436 G4cout << " --- > For " << pReactant->GetName() << G4endl;
437
438 for (auto it2 : it.second)
439 {
440 G4cout << " " << G4BestUnit(it2.first, "Time")
441 << " " << it2.second << G4endl;
442 }
443 }
444}

References fCounterMap, G4BestUnit, G4cout, and G4endl.

Referenced by RemoveAMoleculeAtTime(), and G4DNAUpdateSystemModel::UpdateSystem().

◆ GetNbMoleculeAgainstTime()

const NbMoleculeAgainstTime & G4MoleculeCounter::GetNbMoleculeAgainstTime ( Reactant molecule)

Definition at line 460 of file G4MoleculeCounter.cc.

461{
462 return fCounterMap[molecule];
463}

References fCounterMap.

◆ GetNMoleculesAtTime()

int G4MoleculeCounter::GetNMoleculesAtTime ( Reactant molecule,
double  time 
)

Definition at line 199 of file G4MoleculeCounter.cc.

201{
202 G4bool sameTypeOfMolecule = SearchTimeMap(molecule);
203 return SearchUpperBoundTime(time, sameTypeOfMolecule);
204}
bool G4bool
Definition: G4Types.hh:86
G4bool SearchTimeMap(Reactant *molecule)
int SearchUpperBoundTime(double time, bool sameTypeOfMolecule)

References SearchTimeMap(), and SearchUpperBoundTime().

Referenced by G4DNAEventScheduler::ClearAndReChargeCounter().

◆ GetRecordedMolecules()

G4MoleculeCounter::RecordedMolecules G4MoleculeCounter::GetRecordedMolecules ( )

Definition at line 370 of file G4MoleculeCounter.cc.

371{
372 if (fVerbose > 1)
373 {
374 G4cout << "Entering in G4MoleculeCounter::RecordMolecules" << G4endl;
375 }
376
377 RecordedMolecules output(new ReactantList());
378
379 for (auto it : fCounterMap)
380 {
381 output->push_back(it.first);
382 }
383 return output;
384}
std::unique_ptr< ReactantList > RecordedMolecules
std::vector< Reactant * > ReactantList

References fCounterMap, fVerbose, G4cout, and G4endl.

Referenced by G4DNAEventScheduler::ClearAndReChargeCounter().

◆ GetRecordedTimes()

RecordedTimes G4MoleculeCounter::GetRecordedTimes ( )

Definition at line 388 of file G4MoleculeCounter.cc.

389{
390 RecordedTimes output(new std::set<G4double>);
391
392 for(const auto& it : fCounterMap)
393 {
394 for(const auto& it2 : it.second)
395 {
396 //time = it2->first;
397 output->insert(it2.first);
398 }
399 }
400
401 return output;
402}
std::unique_ptr< std::set< G4double > > RecordedTimes

References fCounterMap.

◆ GetVerbose()

G4int G4MoleculeCounter::GetVerbose ( )

Definition at line 474 of file G4MoleculeCounter.cc.

475{
476 return fVerbose;
477}

References fVerbose.

◆ Initialize()

void G4MoleculeCounter::Initialize ( )
overridevirtual

Implements G4VMoleculeCounter.

Definition at line 83 of file G4MoleculeCounter.cc.

84{
86 while ((mol_iterator)())
87 {
88 if (IsRegistered(mol_iterator.value()->GetDefinition()) == false)
89 {
90 continue;
91 }
92
93 fCounterMap[mol_iterator.value()]; // initialize the second map
94 }
95}
bool IsRegistered(const G4MoleculeDefinition *) override
static G4MoleculeTable * Instance()
G4ConfigurationIterator GetConfigurationIterator()

References fCounterMap, G4MoleculeTable::GetConfigurationIterator(), G4MoleculeTable::Instance(), and IsRegistered().

◆ InitializeInstance()

void G4VMoleculeCounter::InitializeInstance ( )
staticinherited

Definition at line 89 of file G4VMoleculeCounter.cc.

90{
91 if (fpInstance)
92 {
94 }
95}
virtual void Initialize()=0

References G4VMoleculeCounter::fpInstance, and G4VMoleculeCounter::Initialize().

Referenced by G4DNAChemistryManager::InitializeThread().

◆ Instance()

G4MoleculeCounter * G4MoleculeCounter::Instance ( )
static

◆ InUse()

G4bool G4VMoleculeCounter::InUse ( )
staticinherited

◆ IsRegistered()

bool G4MoleculeCounter::IsRegistered ( const G4MoleculeDefinition molDef)
overridevirtual

Reimplemented from G4VMoleculeCounter.

Definition at line 488 of file G4MoleculeCounter.cc.

489{
490 if (fDontRegister.find(molDef) == fDontRegister.end())
491 {
492 return true;
493 }
494 return false;
495}

References fDontRegister.

Referenced by Initialize().

◆ IsTimeCheckedForConsistency()

G4bool G4MoleculeCounter::IsTimeCheckedForConsistency ( ) const

Definition at line 504 of file G4MoleculeCounter.cc.

505{
507}

References fCheckTimeIsConsistentWithScheduler.

◆ RegisterAll()

void G4MoleculeCounter::RegisterAll ( )
overridevirtual

Reimplemented from G4VMoleculeCounter.

Definition at line 499 of file G4MoleculeCounter.cc.

500{
501 fDontRegister.clear();
502}

References fDontRegister.

◆ RemoveAMoleculeAtTime()

void G4MoleculeCounter::RemoveAMoleculeAtTime ( Reactant ,
G4double  time,
const G4ThreeVector position = nullptr,
int  number = 1 
)
overrideprotectedvirtual

Implements G4VMoleculeCounter.

Definition at line 268 of file G4MoleculeCounter.cc.

272{
273 if (fDontRegister[pMolecule->GetDefinition()])
274 {
275 return;
276 }
277
278 if (fVerbose)
279 {
280 G4cout << "G4MoleculeCounter::RemoveAMoleculeAtTime : "
281 << pMolecule->GetName() << " at time : " << G4BestUnit(time, "Time")
282 << G4endl;
283 }
284
286 {
287 if (fabs(time - G4Scheduler::Instance()->GetGlobalTime()) >
288 G4Scheduler::Instance()->GetTimeTolerance())
289 {
291 errMsg << "Time of species "
292 << pMolecule->GetName() << " is "
293 << G4BestUnit(time, "Time") << " while "
294 << " global time is "
295 << G4BestUnit(G4Scheduler::Instance()->GetGlobalTime(), "Time")
296 << G4endl;
297 G4Exception("G4MoleculeCounter::RemoveAMoleculeAtTime",
298 "TIME_DONT_MATCH",
299 FatalException, errMsg);
300 }
301 }
302
303 NbMoleculeAgainstTime& nbMolPerTime = fCounterMap[pMolecule];
304
305 if (nbMolPerTime.empty())
306 {
307 pMolecule->PrintState();
308 Dump();
309 G4String errMsg =
310 "You are trying to remove molecule " + pMolecule->GetName() +
311 " from the counter while this kind of molecules has not been registered yet";
312 G4Exception("G4MoleculeCounter::RemoveAMoleculeAtTime", "",
313 FatalErrorInArgument, errMsg);
314
315 return;
316 }
317 else
318 {
319 NbMoleculeAgainstTime::reverse_iterator it = nbMolPerTime.rbegin();
320
321 if (it == nbMolPerTime.rend())
322 {
323 it--;
324
325 G4String errMsg =
326 "There was no " + pMolecule->GetName() + " recorded at the time or even before the time asked";
327 G4Exception("G4MoleculeCounter::RemoveAMoleculeAtTime", "",
328 FatalErrorInArgument, errMsg);
329 }
330
332 {
333 Dump();
335 errMsg << "Is time going back?? " << pMolecule->GetName()
336 << " is being removed at time " << G4BestUnit(time, "Time")
337 << " while last recorded time was "
338 << G4BestUnit(it->first, "Time") << ".";
339 G4Exception("G4MoleculeCounter::RemoveAMoleculeAtTime",
340 "RETURN_TO_THE_FUTUR",
342 errMsg);
343 }
344
345 double finalN = it->second - number;
346
347 if (finalN < 0)
348 {
349 Dump();
351 errMsg << "After removal of " << number << " species of "
352 << pMolecule->GetName() << " the final number at time "
353 << G4BestUnit(time, "Time") << " is less than zero and so not valid."
354 << " Global time is "
355 << G4BestUnit(G4Scheduler::Instance()->GetGlobalTime(), "Time")
356 << ". Previous selected time is "
357 << G4BestUnit(it->first, "Time")
358 << G4endl;
359 G4Exception("G4MoleculeCounter::RemoveAMoleculeAtTime",
360 "N_INF_0",
361 FatalException, errMsg);
362 }
363
364 nbMolPerTime[time] = finalN;
365 }
366}
@ FatalErrorInArgument
std::map< G4double, G4int, G4::MoleculeCounter::TimePrecision > NbMoleculeAgainstTime

References Dump(), FatalErrorInArgument, FatalException, fCheckTimeIsConsistentWithScheduler, fCounterMap, fDontRegister, G4::MoleculeCounter::TimePrecision::fPrecision, fVerbose, G4BestUnit, G4cout, G4endl, G4Exception(), G4MolecularConfiguration::GetDefinition(), G4MolecularConfiguration::GetName(), G4Scheduler::Instance(), and G4MolecularConfiguration::PrintState().

◆ ResetCounter()

void G4MoleculeCounter::ResetCounter ( )
overridevirtual

Implements G4VMoleculeCounter.

Definition at line 448 of file G4MoleculeCounter.cc.

449{
450 if (fVerbose)
451 {
452 G4cout << " ---> G4MoleculeCounter::ResetCounter" << G4endl;
453 }
454 fCounterMap.clear();
455 fpLastSearch.reset(0);
456}
std::unique_ptr< Search > fpLastSearch

References fCounterMap, fpLastSearch, fVerbose, G4cout, and G4endl.

Referenced by G4DNAEventScheduler::ClearAndReChargeCounter().

◆ SearchTimeMap()

G4bool G4MoleculeCounter::SearchTimeMap ( Reactant molecule)
protected

Definition at line 106 of file G4MoleculeCounter.cc.

107{
108 if (fpLastSearch == nullptr)
109 {
110 fpLastSearch.reset(new Search());
111 }
112 else
113 {
114 if (fpLastSearch->fLowerBoundSet &&
115 fpLastSearch->fLastMoleculeSearched->first == molecule)
116 {
117 return true;
118 }
119 }
120
121 auto mol_it = fCounterMap.find(molecule);
122 fpLastSearch->fLastMoleculeSearched = mol_it;
123
124 if (mol_it != fCounterMap.end())
125 {
126 fpLastSearch->fLowerBoundTime = fpLastSearch->fLastMoleculeSearched->second
127 .end();
128 fpLastSearch->fLowerBoundSet = true;
129 }
130 else
131 {
132 fpLastSearch->fLowerBoundSet = false;
133 }
134
135 return false;
136}

References fCounterMap, and fpLastSearch.

Referenced by GetNMoleculesAtTime().

◆ SearchUpperBoundTime()

int G4MoleculeCounter::SearchUpperBoundTime ( double  time,
bool  sameTypeOfMolecule 
)
protected

Definition at line 140 of file G4MoleculeCounter.cc.

142{
143 auto mol_it = fpLastSearch->fLastMoleculeSearched;
144 if (mol_it == fCounterMap.end())
145 {
146 return 0;
147 }
148
149 NbMoleculeAgainstTime& timeMap = mol_it->second;
150 if (timeMap.empty())
151 {
152 return 0;
153 }
154
155 if (sameTypeOfMolecule == true)
156 {
157 if (fpLastSearch->fLowerBoundSet && fpLastSearch->fLowerBoundTime != timeMap.end())
158 {
159 if (fpLastSearch->fLowerBoundTime->first < time)
160 {
161 auto upperToLast = fpLastSearch->fLowerBoundTime;
162 upperToLast++;
163
164 if (upperToLast == timeMap.end())
165 {
166 return fpLastSearch->fLowerBoundTime->second;
167 }
168
169 if (upperToLast->first > time)
170 {
171 return fpLastSearch->fLowerBoundTime->second;
172 }
173 }
174 }
175 }
176
177 auto up_time_it = timeMap.upper_bound(time);
178
179 if (up_time_it == timeMap.end())
180 {
181 NbMoleculeAgainstTime::reverse_iterator last_time = timeMap.rbegin();
182 return last_time->second;
183 }
184 if (up_time_it == timeMap.begin())
185 {
186 return 0;
187 }
188
189 up_time_it--;
190
191 fpLastSearch->fLowerBoundTime = up_time_it;
192 fpLastSearch->fLowerBoundSet = true;
193
194 return fpLastSearch->fLowerBoundTime->second;
195}

References fCounterMap, and fpLastSearch.

Referenced by GetNMoleculesAtTime().

◆ SetInstance()

void G4VMoleculeCounter::SetInstance ( G4VMoleculeCounter pCounterInstance)
staticinherited

Definition at line 43 of file G4VMoleculeCounter.cc.

44{
45 if (fpInstance != nullptr)
46 {
48 errMsg << "The G4MoleculeCounter was already initialized." << G4endl
49 << "The previous instance will be deleted in order to use yours." << G4endl
50 << "However this can generate conflicts. Make sure you call G4MoleculeCounter::SetInstance"
51 "at the beginning of your application."
52 << "A good place would be ActionInitialization::Build & BuildForMaster"
53 << G4endl;
54
55 G4Exception("G4MoleculeCounter::SetInstance",
56 "SINGLETON_ALREADY_INITIALIZED",
57 JustWarning, errMsg);
58 delete fpInstance;
59 fpInstance = nullptr;
60 }
61
62 fpInstance = pCounterInstance;
63}
@ JustWarning

References G4VMoleculeCounter::fpInstance, G4endl, G4Exception(), and JustWarning.

◆ SetTimeSlice()

void G4MoleculeCounter::SetTimeSlice ( double  timeSlice)
static

◆ SetVerbose()

void G4MoleculeCounter::SetVerbose ( G4int  level)

Definition at line 467 of file G4MoleculeCounter.cc.

468{
469 fVerbose = level;
470}

References fVerbose.

◆ Use()

void G4VMoleculeCounter::Use ( G4bool  flag = true)
staticinherited

Definition at line 99 of file G4VMoleculeCounter.cc.

100{
101 fUse = flag;
102}

References G4VMoleculeCounter::fUse.

Referenced by G4DNAEventScheduler::ClearAndReChargeCounter().

Friends And Related Function Documentation

◆ G4Molecule

friend class G4Molecule
friend

Definition at line 148 of file G4MoleculeCounter.hh.

◆ G4VMoleculeCounter

friend class G4VMoleculeCounter
friend

Definition at line 149 of file G4MoleculeCounter.hh.

Field Documentation

◆ fCheckTimeIsConsistentWithScheduler

G4bool G4MoleculeCounter::fCheckTimeIsConsistentWithScheduler
protected

◆ fCounterMap

CounterMapType G4MoleculeCounter::fCounterMap
protected

◆ fDontRegister

std::map<const G4MoleculeDefinition*, G4bool> G4MoleculeCounter::fDontRegister
protected

◆ fpInstance

G4ThreadLocal G4VMoleculeCounter * G4VMoleculeCounter::fpInstance = nullptr
staticprotectedinherited

◆ fpLastSearch

std::unique_ptr<Search> G4MoleculeCounter::fpLastSearch
protected

Definition at line 146 of file G4MoleculeCounter.hh.

Referenced by ResetCounter(), SearchTimeMap(), and SearchUpperBoundTime().

◆ fUse

G4bool G4VMoleculeCounter::fUse = false
staticprotectedinherited

Definition at line 47 of file G4VMoleculeCounter.hh.

Referenced by G4VMoleculeCounter::InUse(), and G4VMoleculeCounter::Use().

◆ fVerbose

G4int G4MoleculeCounter::fVerbose
protected

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