G4DNAIndependentReactionTimeStepper Class Reference

#include <G4DNAIndependentReactionTimeStepper.hh>

Inheritance diagram for G4DNAIndependentReactionTimeStepper:

Public Member Functions
	G4DNAIndependentReactionTimeStepper ()
	~G4DNAIndependentReactionTimeStepper () override=default
	G4DNAIndependentReactionTimeStepper (const G4DNAIndependentReactionTimeStepper &)=delete
G4DNAIndependentReactionTimeStepper &	operator= (const G4DNAIndependentReactionTimeStepper &)=delete
void	Prepare () override
G4double	CalculateStep (const G4Track &, const G4double &) override
G4double	CalculateMinTimeStep (G4double, G4double) override
void	SetReactionModel (G4VDNAReactionModel *)
G4VDNAReactionModel *	GetReactionModel ()
std::unique_ptr< G4ITReactionChange >	FindReaction (G4ITReactionSet *pReactionSet, const G4double &currentStepTime=0, const G4double &previousStepTime=0, const G4bool &reachedUserStepTimeLimit=false)
void	SetReactionProcess (G4VITReactionProcess *pReactionProcess)
void	SetVerbose (G4int)
Public Member Functions inherited from G4VITTimeStepComputer
	G4VITTimeStepComputer ()
virtual	~G4VITTimeStepComputer ()
	G4VITTimeStepComputer (const G4VITTimeStepComputer &)
G4VITTimeStepComputer &	operator= (const G4VITTimeStepComputer &other)
virtual void	Initialize ()
G4TrackVectorHandle	GetReactants ()
virtual void	ResetReactants ()
G4double	GetSampledMinTimeStep ()
void	SetReactionTable (const G4ITReactionTable *)
const G4ITReactionTable *	GetReactionTable ()

Additional Inherited Members
Static Public Member Functions inherited from G4VITTimeStepComputer
static void	SetTimes (const G4double &, const G4double &)
Protected Attributes inherited from G4VITTimeStepComputer
G4double	fSampledMinTimeStep
G4TrackVectorHandle	fReactants
const G4ITReactionTable *	fpReactionTable
Static Protected Attributes inherited from G4VITTimeStepComputer
static G4ThreadLocal G4double	fCurrentGlobalTime = -1
static G4ThreadLocal G4double	fUserMinTimeStep = -1

Detailed Description

Definition at line 51 of file G4DNAIndependentReactionTimeStepper.hh.

Constructor & Destructor Documentation

G4DNAIndependentReactionTimeStepper() [1/2]

G4DNAIndependentReactionTimeStepper::G4DNAIndependentReactionTimeStepper

(

)

Definition at line 54 of file G4DNAIndependentReactionTimeStepper.cc.

{
  fReactionSet->SortByTime();
}

Referenced by G4DNAIndependentReactionTimeStepper(), and operator=().

~G4DNAIndependentReactionTimeStepper()

G4DNAIndependentReactionTimeStepper::~G4DNAIndependentReactionTimeStepper ( )

overridedefault

G4DNAIndependentReactionTimeStepper() [2/2]

G4DNAIndependentReactionTimeStepper::G4DNAIndependentReactionTimeStepper ( const G4DNAIndependentReactionTimeStepper & )

delete

Member Function Documentation

CalculateMinTimeStep()

G4double G4DNAIndependentReactionTimeStepper::CalculateMinTimeStep ( G4double , G4double definedMinTimeStep )

overridevirtual

Implements G4VITTimeStepComputer.

Definition at line 357 of file G4DNAIndependentReactionTimeStepper.cc.

{
  G4double fTSTimeStep = DBL_MAX;
  fCheckedTracks.clear();
 
  for (auto pTrack : *fpTrackContainer->GetMainList()) {
    if (pTrack == nullptr) {
      G4ExceptionDescription msg;
      msg << "No track found.";
      G4Exception("G4DNAIndependentReactionTimeStepper::CalculateMinTimeStep",
                  "G4DNAIndependentReactionTimeStepper006", FatalErrorInArgument,
                  msg);
      continue;
    }
 
    G4TrackStatus trackStatus = pTrack->GetTrackStatus();
    if (trackStatus == fStopAndKill || trackStatus == fStopButAlive) {
      continue;
    }
 
    G4double sampledMinTimeStep = CalculateStep(*pTrack, definedMinTimeStep);
    G4TrackVectorHandle reactants = GetReactants();
 
    if (sampledMinTimeStep < fTSTimeStep) {
      fTSTimeStep = sampledMinTimeStep;
      if (reactants) {
        fReactionSet->AddReactions(fTSTimeStep, const_cast<G4Track*>(pTrack), std::move(reactants));
 
        fSampledPositions[pTrack->GetTrackID()] = pTrack->GetPosition();
        for (const auto& it : *fReactants) {
          auto pTrackB = it;
          fSampledPositions[pTrackB->GetTrackID()] = pTrackB->GetPosition();
        }
        ResetReactants();
      }
    }
    else if (fTSTimeStep == sampledMinTimeStep && G4bool(reactants)) {
      fReactionSet->AddReactions(fTSTimeStep, const_cast<G4Track*>(pTrack), std::move(reactants));
 
      fSampledPositions[pTrack->GetTrackID()] = pTrack->GetPosition();
      for (const auto& it : *fReactants) {
        auto pTrackB = it;
        fSampledPositions[pTrackB->GetTrackID()] = pTrackB->GetPosition();
      }
      ResetReactants();
    }
    else if (reactants) {
      ResetReactants();
    }
  }
  return fTSTimeStep;
}

CalculateStep()

G4double G4DNAIndependentReactionTimeStepper::CalculateStep ( const G4Track & trackA, const G4double & userMinTimeStep )

overridevirtual

Implements G4VITTimeStepComputer.

Definition at line 75 of file G4DNAIndependentReactionTimeStepper.cc.

{
  auto pMoleculeA = GetMolecule(trackA);
  InitializeForNewTrack();
  fUserMinTimeStep = userMinTimeStep;
  fCheckedTracks.insert(trackA.GetTrackID());
 
#ifdef G4VERBOSE
  if (fVerbose != 0) {
    G4cout << "________________________________________________________________"
              "_______"
           << G4endl;
    G4cout << "G4DNAIndependentReactionTimeStepper::CalculateStep" << G4endl;
    G4cout << "Check done for molecule : " << pMoleculeA->GetName() << " (" << trackA.GetTrackID()
           << ") " << G4endl;
  }
#endif
 
  auto pMolConfA = pMoleculeA->GetMolecularConfiguration();
 
  const auto pReactantList = fMolecularReactionTable->CanReactWith(pMolConfA);
 
  if (pReactantList == nullptr) {
    if(fVerbose > 1) {
      G4ExceptionDescription msg;
      msg << "G4DNAIndependentReactionTimeStepper::CalculateStep will return infinity "
             "for the reaction because the molecule "
          << pMoleculeA->GetName() << " does not have any reactants given in the reaction table."
          << G4endl;
      G4Exception("G4DNAIndependentReactionTimeStepper::CalculateStep",
                  "G4DNAIndependentReactionTimeStepper03", JustWarning, msg);
    }
    return DBL_MAX;
  }
 
  auto nbReactives = (G4int)pReactantList->size();
 
  if (nbReactives == 0) {
    if(fVerbose != 0){
      G4ExceptionDescription msg;
      msg << "G4DNAIndependentReactionTimeStepper::CalculateStep will "
             "return infinity "
             "for the reaction because the molecule "
          << pMoleculeA->GetName() << " does not have any reactants given in the reaction table."
          << "This message can also result from a wrong implementation of "
             "the reaction table."
          << G4endl;
      G4Exception("G4DNAIndependentReactionTimeStepper::CalculateStep",
                  "G4DNAIndependentReactionTimeStepper04", JustWarning, msg);
    }
    return DBL_MAX;
  }
  fReactants = std::make_shared<vector<G4Track*>>();
  fReactionModel->Initialise(pMolConfA, trackA);
  for (G4int i = 0; i < nbReactives; ++i) {
    auto pMoleculeB = (*pReactantList)[i];
    G4int key = pMoleculeB->GetMoleculeID();
 
    // fRCutOff = G4IRTUtils::GetRCutOff(1 * ps);
    fRCutOff = G4IRTUtils::GetRCutOff();
    //______________________________________________________________
    // Retrieve reaction range
    const G4double Reff = fReactionModel->GetReactionRadius(i);
    std::vector<std::pair<G4TrackList::iterator, G4double>> resultIndices;
    resultIndices.clear();
    G4ChemicalMoleculeFinder::Instance()->FindNearestInRange(trackA, key, fRCutOff, resultIndices);
 
    if (resultIndices.empty()) {
      continue;
    }
    for (auto& it : resultIndices) {
      G4Track* pTrackB = *(std::get<0>(it));
 
      if (pTrackB == &trackA) {
        continue;
      }
      if (pTrackB == nullptr) {
        G4ExceptionDescription exceptionDescription;
        exceptionDescription << "No trackB no valid";
        G4Exception(
          "G4DNAIndependentReactionTimeStepper"
          "::CalculateStep()",
          "G4DNAIndependentReactionTimeStepper007", FatalException, exceptionDescription);
      }
      else {
        if (fCheckedTracks.find(pTrackB->GetTrackID()) != fCheckedTracks.end()) {
          continue;
        }
 
        Utils utils(trackA, *pTrackB);
 
        auto pMolB = GetMolecule(pTrackB);
        auto pMolConfB = pMolB->GetMolecularConfiguration();
        G4double distance = (trackA.GetPosition() - pTrackB->GetPosition()).mag();
        if (distance * distance < Reff * Reff) {
          auto reactionData = fMolecularReactionTable->GetReactionData(pMolConfA, pMolConfB);
          if (G4Scheduler::Instance()->GetGlobalTime() == G4Scheduler::Instance()->GetStartTime()) {
            if (reactionData->GetProbability() > G4UniformRand()) {
              if (!fHasAlreadyReachedNullTime) {
                fReactants->clear();
                fHasAlreadyReachedNullTime = true;
              }
              fSampledMinTimeStep = 0.;
              CheckAndRecordResults(utils);
            }
          }
        }
        else {
          G4double tempMinET = GetTimeToEncounter(trackA, *pTrackB);
          if (tempMinET < 0 || tempMinET > G4Scheduler::Instance()->GetEndTime()) {
            continue;
          }
          if (tempMinET >= fSampledMinTimeStep) {
            continue;
          }
          fSampledMinTimeStep = tempMinET;
          fReactants->clear();
          CheckAndRecordResults(utils);
        }
      }
    }
  }
 
#ifdef G4VERBOSE
  if (fVerbose != 0) {
    G4cout << "G4DNAIndependentReactionTimeStepper::CalculateStep will finally "
              "return :"
           << G4BestUnit(fSampledMinTimeStep, "Time") << G4endl;
 
    if (fVerbose > 1) {
      G4cout << "Selected reactants for trackA: " << pMoleculeA->GetName() << " ("
             << trackA.GetTrackID() << ") are: ";
 
      vector<G4Track*>::iterator it;
      for (it = fReactants->begin(); it != fReactants->end(); it++) {
        G4Track* trackB = *it;
        G4cout << GetMolecule(trackB)->GetName() << " (" << trackB->GetTrackID() << ") \t ";
      }
      G4cout << G4endl;
    }
  }
#endif
  return fSampledMinTimeStep;
}

Referenced by CalculateMinTimeStep().

FindReaction()

std::unique_ptr< G4ITReactionChange > G4DNAIndependentReactionTimeStepper::FindReaction	(	G4ITReactionSet *	pReactionSet,
		const G4double &	currentStepTime = 0,
		const G4double &	previousStepTime = 0,
		const G4bool &	reachedUserStepTimeLimit = false )

Definition at line 266 of file G4DNAIndependentReactionTimeStepper.cc.

{
  if (pReactionSet == nullptr) {
    return nullptr;
  }
 
  G4ITReactionPerTime& reactionPerTime = pReactionSet->GetReactionsPerTime();
  if (reactionPerTime.empty()) {
    return nullptr;
  }
  for (auto reaction_i = reactionPerTime.begin(); reaction_i != reactionPerTime.end();
       reaction_i = reactionPerTime.begin())
  {
    if ((*reaction_i)->GetTime() > currentStepTime) {
      fReactionSet->CleanAllReaction();
      return nullptr;
    }
 
    G4Track* pTrackA = (*reaction_i)->GetReactants().first;
    if (pTrackA->GetTrackStatus() == fStopAndKill) {
      continue;
    }
    G4Track* pTrackB = (*reaction_i)->GetReactant(pTrackA);
    if (pTrackB->GetTrackStatus() == fStopAndKill) {
      continue;
    }
 
    if (pTrackB == pTrackA) {
      G4ExceptionDescription msg;
      msg << "The IT reaction process sent back a reaction "
                              "between trackA and trackB. ";
      msg << "The problem is trackA == trackB";
      G4Exception("G4DNAIndependentReactionTimeStepper::FindReaction",
                  "G4DNAIndependentReactionTimeStepper02", FatalErrorInArgument,
                  msg);
    }
    pReactionSet->SelectThisReaction(*reaction_i);
    if (fpReactionProcess != nullptr
        && fpReactionProcess->TestReactibility(*pTrackA, *pTrackB, currentStepTime, false))
    {
      if ((fSampledPositions.find(pTrackA->GetTrackID()) == fSampledPositions.end()
           && (fSampledPositions.find(pTrackB->GetTrackID()) == fSampledPositions.end())))
      {
        G4ExceptionDescription msg;
        msg << "The positions of trackA and trackB have no counted ";
        G4Exception("G4DNAIndependentReactionTimeStepper::FindReaction",
                    "G4DNAIndependentReactionTimeStepper0001", FatalErrorInArgument,
                    msg);
      }
 
      pTrackA->SetPosition(fSampledPositions[pTrackA->GetTrackID()]);
      pTrackB->SetPosition(fSampledPositions[pTrackB->GetTrackID()]);
      auto pReactionChange = fpReactionProcess->MakeReaction(*pTrackA, *pTrackB);
      if (pReactionChange == nullptr) {
        return nullptr;
      }
      return pReactionChange;
    }
  }
  return nullptr;
}

GetReactionModel()

G4VDNAReactionModel * G4DNAIndependentReactionTimeStepper::GetReactionModel

(

)

Definition at line 335 of file G4DNAIndependentReactionTimeStepper.cc.

{
  return fReactionModel;
}

operator=()

G4DNAIndependentReactionTimeStepper & G4DNAIndependentReactionTimeStepper::operator= ( const G4DNAIndependentReactionTimeStepper & )

delete

Prepare()

void G4DNAIndependentReactionTimeStepper::Prepare ( )

overridevirtual

Reimplemented from G4VITTimeStepComputer.

Definition at line 59 of file G4DNAIndependentReactionTimeStepper.cc.

{
  G4VITTimeStepComputer::Prepare();
  fSampledPositions.clear();
  BuildChemicalMoleculeFinder()
}

SetReactionModel()

void G4DNAIndependentReactionTimeStepper::SetReactionModel

(

G4VDNAReactionModel *

pReactionModel

)

Definition at line 330 of file G4DNAIndependentReactionTimeStepper.cc.

{
  fReactionModel = pReactionModel;
}

SetReactionProcess()

void G4DNAIndependentReactionTimeStepper::SetReactionProcess

(

G4VITReactionProcess *

pReactionProcess

)

Definition at line 353 of file G4DNAIndependentReactionTimeStepper.cc.

{
  fpReactionProcess = pReactionProcess;
}

SetVerbose()

void G4DNAIndependentReactionTimeStepper::SetVerbose

(

G4int

flag

)

Definition at line 340 of file G4DNAIndependentReactionTimeStepper.cc.

{
  fVerbose = flag;
}

---

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

• G4DNAIndependentReactionTimeStepper.hh
• G4DNAIndependentReactionTimeStepper.cc