G4DNAIRTMoleculeEncounterStepper Class Reference

#include <G4DNAIRTMoleculeEncounterStepper.hh>

Inheritance diagram for G4DNAIRTMoleculeEncounterStepper:

Public Member Functions
	G4DNAIRTMoleculeEncounterStepper ()
virtual	~G4DNAIRTMoleculeEncounterStepper ()
	G4DNAIRTMoleculeEncounterStepper (const G4DNAIRTMoleculeEncounterStepper &)=delete
G4DNAIRTMoleculeEncounterStepper &	operator= (const G4DNAIRTMoleculeEncounterStepper &)=delete
virtual void	Prepare ()
virtual G4double	CalculateStep (const G4Track &, const G4double &)
virtual G4double	CalculateMinTimeStep (G4double, G4double)
void	SetReactionModel (G4VDNAReactionModel *)
G4VDNAReactionModel *	GetReactionModel ()
void	SetVerbose (int)
Public Member Functions inherited from G4VITTimeStepComputer
	G4VITTimeStepComputer ()
virtual	~G4VITTimeStepComputer ()
	G4VITTimeStepComputer (const G4VITTimeStepComputer &)
G4VITTimeStepComputer &	operator= (const G4VITTimeStepComputer &other)
virtual void	Initialize ()
virtual void	Prepare ()
virtual G4double	CalculateStep (const G4Track &, const G4double &)=0
virtual G4double	CalculateMinTimeStep (G4double, G4double)=0
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

Given a molecule G4DNAIRTMoleculeEncounterStepper will calculate for its possible reactants what will be the minimum encounter time and the associated molecules.*

This model includes dynamical time steps as explained in "Computer-Aided Stochastic Modeling of the Radiolysis of Liquid Water", V. Michalik, M. Begusová, E. A. Bigildeev, Radiation Research, Vol. 149, No. 3 (Mar., 1998), pp. 224-236

Definition at line 59 of file G4DNAIRTMoleculeEncounterStepper.hh.

Constructor & Destructor Documentation

G4DNAIRTMoleculeEncounterStepper() [1/2]

G4DNAIRTMoleculeEncounterStepper::G4DNAIRTMoleculeEncounterStepper

(

)

Definition at line 67 of file G4DNAIRTMoleculeEncounterStepper.cc.

    : G4VITTimeStepComputer()
    , fHasAlreadyReachedNullTime(false)
    , fMolecularReactionTable(reference_cast<const G4DNAMolecularReactionTable*>(fpReactionTable))
    , fReactionModel(nullptr)
    , fVerbose(0)
{
    fpTrackContainer = G4ITTrackHolder::Instance();
    fReactionSet = G4ITReactionSet::Instance();
}

~G4DNAIRTMoleculeEncounterStepper()

G4DNAIRTMoleculeEncounterStepper::~G4DNAIRTMoleculeEncounterStepper ( )

virtualdefault

G4DNAIRTMoleculeEncounterStepper() [2/2]

G4DNAIRTMoleculeEncounterStepper::G4DNAIRTMoleculeEncounterStepper ( const G4DNAIRTMoleculeEncounterStepper &  )

delete

Member Function Documentation

CalculateMinTimeStep()

G4double G4DNAIRTMoleculeEncounterStepper::CalculateMinTimeStep ( G4double  currentGlobalTime, G4double  definedMinTimeStep  )

virtual

Implements G4VITTimeStepComputer.

Definition at line 431 of file G4DNAIRTMoleculeEncounterStepper.cc.

                                                                                                                      {
 
    G4bool start = true;
    G4bool active = false;
 
    fUserMinTimeStep = definedMinTimeStep;
 
    if(fReactionSet->Empty()){
        if(currentGlobalTime == G4Scheduler::Instance()->GetStartTime()){
 
            for (auto pTrack : *fpTrackContainer->GetMainList())
            {
                if (pTrack == nullptr)
                {
                    G4ExceptionDescription exceptionDescription;
                    exceptionDescription << "No track found.";
                    G4Exception("G4Scheduler::CalculateMinStep", "ITScheduler006",
                                FatalErrorInArgument, exceptionDescription);
                    continue;
                }
 
                G4TrackStatus trackStatus = pTrack->GetTrackStatus();
                if (trackStatus == fStopAndKill || trackStatus == fStopButAlive)
                {
                    start = false;
                    continue;
                }
                active = true;
            }
 
            if(start == true){
                return -1;
            }else if(active == false){
                G4Scheduler::Instance()->Stop();
                return fSampledMinTimeStep;
            }else{
                return fSampledMinTimeStep;
            }
 
        }else{
            for (auto pTrack : *fpTrackContainer->GetMainList())
            {
                pTrack->SetGlobalTime(G4Scheduler::Instance()->GetEndTime());
            }
            return fSampledMinTimeStep;
        }
    }
 
    auto fReactionSetInTime = fReactionSet->GetReactionsPerTime();
    fSampledMinTimeStep = fReactionSetInTime.begin()->get()->GetTime() - currentGlobalTime;
 
    return fSampledMinTimeStep;
}

CalculateStep()

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

virtual

Implements G4VITTimeStepComputer.

Definition at line 107 of file G4DNAIRTMoleculeEncounterStepper.cc.

{
 
    auto pMoleculeA = GetMolecule(trackA);
    InitializeForNewTrack();
    fUserMinTimeStep = userMinTimeStep;
 
#ifdef G4VERBOSE
    if (fVerbose)
    {
        G4cout
            << "_______________________________________________________________________"
            << G4endl;
        G4cout << "G4DNAMoleculeEncounterStepper::CalculateStep" << G4endl;
        G4cout << "Check done for molecule : " << pMoleculeA->GetName()
            << " (" << trackA.GetTrackID() << ") "
            << G4endl;
    }
#endif
 
    //__________________________________________________________________
    // Retrieve general informations for making reactions
    auto pMolConfA = pMoleculeA->GetMolecularConfiguration();
 
    const auto pReactantList = fMolecularReactionTable->CanReactWith(pMolConfA);
 
    if (!pReactantList)
    {
#ifdef G4VERBOSE
        //    DEBUG
        if (fVerbose > 1)
        {
            G4cout << "!!!!!!!!!!!!!!!!!!!!" << G4endl;
            G4cout << "!!! WARNING" << G4endl;
            G4cout << "G4MoleculeEncounterStepper::CalculateStep will return infinity "
                "for the reaction because the molecule "
                << pMoleculeA->GetName()
                << " does not have any reactants given in the reaction table."
                << G4endl;
            G4cout << "!!!!!!!!!!!!!!!!!!!!" << G4endl;
        }
#endif
        return DBL_MAX;
    }
 
    G4int nbReactives = pReactantList->size();
 
    if (nbReactives == 0)
    {
#ifdef G4VERBOSE
        //    DEBUG
        if (fVerbose)
        {
            // TODO replace with the warning mode of G4Exception
            G4cout << "!!!!!!!!!!!!!!!!!!!!" << G4endl;
            G4cout << "!!! WARNING" << G4endl;
            G4cout << "G4MoleculeEncounterStepper::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;
            G4cout << "!!!!!!!!!!!!!!!!!!!!" << G4endl;
        }
#endif
        return DBL_MAX;
    }
 
    fReactants.reset(new vector<G4Track*>());
    fReactionModel->Initialise(pMolConfA, trackA);
 
    //__________________________________________________________________
    // Start looping on possible reactants
    for (G4int i = 0; i < nbReactives; i++)
    {
        auto pMoleculeB = (*pReactantList)[i];
 
        //______________________________________________________________
        // Retrieve reaction range
        const G4double R = fReactionModel->GetReactionRadius(i);
 
        //______________________________________________________________
        // Use KdTree algorithm to find closest reactants
        G4KDTreeResultHandle resultsNearest(
            G4MoleculeFinder::Instance()->FindNearest(pMoleculeA,
                                                      pMoleculeB->GetMoleculeID()));
 
        if (resultsNearest == 0) continue;
 
        G4double r2 = resultsNearest->GetDistanceSqr();
        Utils utils(trackA, pMoleculeB);
 
        if (r2 <= R * R) // ==> Record in range
        {
            // Entering in this condition may due to the fact that molecules are very close
            // to each other
            // Therefore, if we only take the nearby reactant into account, it might have already
            // reacted. Instead, we will take all possible reactants that satisfy the condition r<R
 
            if (fHasAlreadyReachedNullTime == false)
            {
                fReactants->clear();
                fHasAlreadyReachedNullTime = true;
            }
 
            fSampledMinTimeStep = 0.;
            G4KDTreeResultHandle resultsInRange(
                G4MoleculeFinder::Instance()->FindNearestInRange(pMoleculeA,
                                                                 pMoleculeB->GetMoleculeID(),
                                                                 R));
            CheckAndRecordResults(utils,
#ifdef G4VERBOSE
                                  R,
#endif
                                  resultsInRange);
        }
        else
        {
            G4double r = sqrt(r2);
            G4double tempMinET = pow(r - R, 2) / utils.fConstant;
            // constant = 16 * (fDA + fDB + 2*sqrt(fDA*fDB))
 
            if (tempMinET <= fSampledMinTimeStep)
            {
                if (fUserMinTimeStep < DBL_MAX/*IsInf(fUserMinTimeStep) == false*/
                    && tempMinET <= fUserMinTimeStep) // ==> Record in range
                {
                    if (fSampledMinTimeStep > fUserMinTimeStep)
                    {
                        fReactants->clear();
                    }
 
                    fSampledMinTimeStep = fUserMinTimeStep;
 
                    G4double range = R + sqrt(fUserMinTimeStep*utils.fConstant);
 
                    G4KDTreeResultHandle resultsInRange(
                        G4MoleculeFinder::Instance()->
                        FindNearestInRange(pMoleculeA,
                                           pMoleculeB->GetMoleculeID(),
                                           range));
 
                    CheckAndRecordResults(utils,
#ifdef G4VERBOSE
                                          range,
#endif
                                          resultsInRange);
                }
                else // ==> Record nearest
                {
                    if (tempMinET < fSampledMinTimeStep)
                        // to avoid cases where fSampledMinTimeStep == tempMinET
                    {
                        fSampledMinTimeStep = tempMinET;
                        fReactants->clear();
                    }
 
                    CheckAndRecordResults(utils,
#ifdef G4VERBOSE
                                          R,
#endif
                                          resultsNearest);
                }
            }
        }
    }
 
#ifdef G4VERBOSE
    if (fVerbose)
    {
        G4cout << "G4MoleculeEncounterStepper::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;
}

GetReactionModel()

G4VDNAReactionModel * G4DNAIRTMoleculeEncounterStepper::GetReactionModel

(

)

Definition at line 421 of file G4DNAIRTMoleculeEncounterStepper.cc.

{
    return fReactionModel;
}

operator=()

G4DNAIRTMoleculeEncounterStepper & G4DNAIRTMoleculeEncounterStepper::operator= ( const G4DNAIRTMoleculeEncounterStepper &  )

delete

Prepare()

void G4DNAIRTMoleculeEncounterStepper::Prepare ( )

virtual

Reimplemented from G4VITTimeStepComputer.

Definition at line 80 of file G4DNAIRTMoleculeEncounterStepper.cc.

{
    fSampledMinTimeStep = DBL_MAX;
    if(G4Scheduler::Instance()->GetGlobalTime() == G4Scheduler::Instance()->GetStartTime()){
        G4VITTimeStepComputer::Prepare();
        G4MoleculeFinder::Instance()->UpdatePositionMap();
    }
}

SetReactionModel()

void G4DNAIRTMoleculeEncounterStepper::SetReactionModel

(

G4VDNAReactionModel *

pReactionModel

)

Definition at line 416 of file G4DNAIRTMoleculeEncounterStepper.cc.

{
    fReactionModel = pReactionModel;
}

SetVerbose()

void G4DNAIRTMoleculeEncounterStepper::SetVerbose

(

int

flag

)

Definition at line 426 of file G4DNAIRTMoleculeEncounterStepper.cc.

{
    fVerbose = flag;
}

---

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

• G4DNAIRTMoleculeEncounterStepper.hh
• G4DNAIRTMoleculeEncounterStepper.cc