Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DNASmoluchowskiReactionModel.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4DNASmoluchowskiReactionModel.cc 95948 2016-03-03 10:40:33Z gcosmo $
27 //
29 #include "Randomize.hh"
30 #include "G4Track.hh"
32 #include "G4UnitsTable.hh"
33 #include "G4Molecule.hh"
34 //#include "G4Scheduler.hh"
35 #include "G4Exp.hh"
36 
39 {
40  fReactionData = 0;
41 }
42 
44  G4VDNAReactionModel(__right)
45 {
46  fReactionData = 0;
47 }
48 
49 G4DNASmoluchowskiReactionModel& G4DNASmoluchowskiReactionModel::operator=(const G4DNASmoluchowskiReactionModel& right)
50 {
51  if (this == &right) return *this;
52  fReactionData = 0;
53  return *this;
54 }
55 
57 {
58  fReactionData = 0;
59 }
60 
62  const G4Track&)
63 {
64  fReactionData = fReactionTable->GetReactionData(__molecule);
65 }
66 
67 void
70 {
71  fReactionData = fReactionTable->GetReactionData(__molecule);
72 }
73 
77 {
78  G4double __output = fReactionTable->GetReactionData(__mol1, __mol2)
80  return __output;
81 }
82 
84 {
85  G4double __output = (*fReactionData)[__i]->GetEffectiveReactionRadius();
86  return __output;
87 }
88 
90  const G4Track& __trackB,
91  const G4double __R,
92  G4double& __r,
93  const G4bool __alongStepReaction)
94 {
95  G4double postStepSeparation = 0;
96  bool do_break = false;
97  G4double R2 = __R * __R;
98  int k = 0;
99 
100  for (; k < 3; k++)
101  {
102  postStepSeparation += std::pow(
103  __trackA.GetPosition()[k] - __trackB.GetPosition()[k], 2);
104 
105  if (postStepSeparation > R2)
106  {
107  do_break = true;
108  break;
109  }
110  }
111 
112  if (do_break == false)
113  {
114  // The loop was not break
115  // => __r^2 < __R^2
116  __r = std::sqrt(postStepSeparation);
117  return true;
118  }
119  else if (__alongStepReaction == true)
120  {
121  //G4cout << "alongStepReaction==true" << G4endl;
122  //Along step cheack and
123  // the loop has break
124 
125  // Continue loop
126  for (; k < 3; k++)
127  {
128  postStepSeparation += std::pow(
129  __trackA.GetPosition()[k] - __trackB.GetPosition()[k], 2);
130  }
131  // Use Green approach : the Brownian bridge
132  __r = (postStepSeparation = std::sqrt(postStepSeparation));
133 
134  G4Molecule* __moleculeA = GetMolecule(__trackA);
135  G4Molecule* __moleculeB = GetMolecule(__trackB);
136 
137  G4double __D = __moleculeA->GetDiffusionCoefficient()
138  + __moleculeB->GetDiffusionCoefficient();
139 
140  G4ThreeVector __preStepPositionA = __trackA.GetStep()->GetPreStepPoint()
141  ->GetPosition();
142  G4ThreeVector __preStepPositionB = __trackB.GetStep()->GetPreStepPoint()
143  ->GetPosition();
144 
145  if (__preStepPositionA == __trackA.GetPosition())
146  {
147  G4ExceptionDescription exceptionDescription;
148  exceptionDescription << "The molecule : " << __moleculeA->GetName();
149  exceptionDescription << " with track ID :" << __trackA.GetTrackID();
150  exceptionDescription << " did not move since the previous step." << G4endl;
151  exceptionDescription << "Current position : "
152  << G4BestUnit(__trackA.GetPosition(), "Length")
153  << G4endl;
154  exceptionDescription << "Previous position : "
155  << G4BestUnit(__preStepPositionA, "Length") << G4endl;
156  G4Exception("G4DNASmoluchowskiReactionModel::FindReaction",
157  "G4DNASmoluchowskiReactionModel", FatalErrorInArgument,
158  exceptionDescription);
159  }
160 
161  G4double __preStepSeparation =
162  (__preStepPositionA - __preStepPositionB).mag();
163 
164  //===================================
165  // Brownian bridge
166 
167 
168 // if(G4Scheduler::Instance()->GetTimeStep() != __trackB.GetStep()->GetDeltaTime())
169 // {
170 // G4cout << G4Scheduler::Instance()->GetTimeStep() << G4endl;
171 // G4cout << __trackB.GetStep()->GetDeltaTime() << G4endl;
172 // assert(G4Scheduler::Instance()->GetTimeStep() == __trackB.GetStep()->GetDeltaTime());
173 // }
174 
175  G4double __probabiltyOfEncounter = G4Exp(
176  -(__preStepSeparation - __R) * (postStepSeparation - __R) / (__D
177  * (__trackB.GetStep()->GetDeltaTime())));
178  G4double __selectedPOE = G4UniformRand();
179 
180  if (__selectedPOE <= __probabiltyOfEncounter) return true;
181  //===================================
182  }
183 
184  return false;
185 }
virtual G4bool FindReaction(const G4Track &, const G4Track &, const G4double, G4double &, const G4bool)
const G4DNAMolecularReactionData * GetReactionData(G4MolecularConfiguration *, G4MolecularConfiguration *) const
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
const G4ThreeVector & GetPosition() const
G4double GetDiffusionCoefficient() const
Definition: G4Molecule.cc:560
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
const G4Step * GetStep() const
int G4int
Definition: G4Types.hh:78
const G4DNAMolecularReactionTable * fReactionTable
G4StepPoint * GetPreStepPoint() const
const G4String & GetName() const
Definition: G4Molecule.cc:356
#define G4UniformRand()
Definition: Randomize.hh:97
const G4ThreeVector & GetPosition() const
bool G4bool
Definition: G4Types.hh:79
virtual G4double GetReactionRadius(G4MolecularConfiguration *, G4MolecularConfiguration *)
G4double GetDeltaTime() const
G4int GetTrackID() const
G4Molecule * GetMolecule(const G4Track &track)
Definition: G4Molecule.cc:76
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
#define G4endl
Definition: G4ios.hh:61
virtual void Initialise(G4MolecularConfiguration *, const G4Track &)
virtual void InitialiseToPrint(G4MolecularConfiguration *)
double G4double
Definition: G4Types.hh:76