| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // G4BOptnForceCommonTruncatedExp <<
27 // ------------------------------------------- <<
28 #include "G4BOptnForceCommonTruncatedExp.hh" 26 #include "G4BOptnForceCommonTruncatedExp.hh"
29 #include "G4ILawCommonTruncatedExp.hh" 27 #include "G4ILawCommonTruncatedExp.hh"
30 #include "G4ILawForceFreeFlight.hh" <<
31 #include "G4TransportationManager.hh" 28 #include "G4TransportationManager.hh"
32 29
33 #include "Randomize.hh" 30 #include "Randomize.hh"
34 #include "G4BiasingProcessInterface.hh" 31 #include "G4BiasingProcessInterface.hh"
35 32
36 G4BOptnForceCommonTruncatedExp:: << 33 G4BOptnForceCommonTruncatedExp::G4BOptnForceCommonTruncatedExp(G4String name)
37 G4BOptnForceCommonTruncatedExp(const G4String& << 34 : G4VBiasingOperation(name),
38 : G4VBiasingOperation(name) << 35 fInteractionOccured(false)
39 { 36 {
40 fCommonTruncatedExpLaw = new G4ILawCommonTru 37 fCommonTruncatedExpLaw = new G4ILawCommonTruncatedExp("ExpLawForOperation"+name);
41 fForceFreeFlightLaw = new G4ILawForceFree << 38 fCommonTruncatedExpLaw->SetOperation( this );
>> 39 fTotalCrossSection = 0.0;
42 } 40 }
43 41
44 G4BOptnForceCommonTruncatedExp::~G4BOptnForceC 42 G4BOptnForceCommonTruncatedExp::~G4BOptnForceCommonTruncatedExp()
>> 43 {}
>> 44
>> 45 const G4VBiasingInteractionLaw* G4BOptnForceCommonTruncatedExp::ProvideOccurenceBiasingInteractionLaw( const G4BiasingProcessInterface* )
>> 46 {
>> 47 return fCommonTruncatedExpLaw;
>> 48 }
>> 49
>> 50 G4ForceCondition G4BOptnForceCommonTruncatedExp::ProposeForceCondition( const G4ForceCondition )
>> 51 {
>> 52 return Forced;
>> 53 }
>> 54
>> 55 G4double G4BOptnForceCommonTruncatedExp::ProposeAlongStepLimit( const G4BiasingProcessInterface* callingProcess )
45 { 56 {
46 delete fCommonTruncatedExpLaw; << 57 if ( callingProcess->GetWrappedProcess() == fProcessToApply )
47 delete fForceFreeFlightLaw; << 58 return fCommonTruncatedExpLaw->GetInteractionDistance();
>> 59 else return DBL_MAX;
48 } 60 }
49 61
50 const G4VBiasingInteractionLaw* G4BOptnForceCo << 62 G4GPILSelection G4BOptnForceCommonTruncatedExp::ProposeGPILSelection( const G4GPILSelection )
51 ProvideOccurenceBiasingInteractionLaw( const G << 63 {
52 G4Force << 64 return CandidateForSelection;
53 { << 65 }
54 if ( callingProcess->GetWrappedProcess() == << 66
55 { << 67 G4bool G4BOptnForceCommonTruncatedExp::DenyProcessPostStepDoIt( const G4BiasingProcessInterface* callingProcess,
56 proposeForceCondition = Forced; << 68 const G4Track*, const G4Step* step, G4double& )
57 return fCommonTruncatedExpLaw; << 69 {
58 } << 70 G4double processGPIL = callingProcess->GetPostStepGPIL() < callingProcess->GetAlongStepGPIL() ?
59 else << 71 callingProcess->GetPostStepGPIL() : callingProcess->GetAlongStepGPIL() ;
60 { <<
61 proposeForceCondition = Forced; <<
62 return fForceFreeFlightLaw; <<
63 } <<
64 } <<
65 <<
66 G4GPILSelection G4BOptnForceCommonTruncatedExp <<
67 ProposeGPILSelection( const G4GPILSelection ) <<
68 { <<
69 return NotCandidateForSelection; <<
70 } <<
71 <<
72 G4VParticleChange* G4BOptnForceCommonTruncated <<
73 ApplyFinalStateBiasing( const G4BiasingProcess <<
74 const G4Track* track, <<
75 const G4Step* step, <<
76 G4bool& forceFin <<
77 { <<
78 if ( callingProcess->GetWrappedProcess() != <<
79 { <<
80 forceFinalState = true; <<
81 fDummyParticleChange.Initialize( *track ); <<
82 return &fDummyParticleChange; <<
83 } <<
84 if ( fInteractionOccured ) <<
85 { <<
86 forceFinalState = true; <<
87 fDummyParticleChange.Initialize( *track ); <<
88 return &fDummyParticleChange; <<
89 } <<
90 <<
91 // -- checks if process won the GPIL race: <<
92 G4double processGPIL = callingProcess->GetPo <<
93 < callingProcess->GetAl <<
94 ? callingProcess->GetPo <<
95 : callingProcess->GetAl <<
96 if ( processGPIL <= step->GetStepLength() ) 72 if ( processGPIL <= step->GetStepLength() )
97 { << 73 {
98 // -- if process won, wrapped process prod << 74 G4bool denyProcess = (callingProcess->GetWrappedProcess() != fProcessToApply );
99 // -- In this case, the weight for occurre << 75 fInteractionOccured = (fInteractionOccured || !denyProcess);
100 // -- by the callingProcess, at exit of pr << 76 return denyProcess;
101 // -- selected by "forceFinalState = false << 77 }
102 forceFinalState = false; << 78 else return true;
103 fInteractionOccured = true; << 79 }
104 return callingProcess->GetWrappedProcess() << 80
105 } << 81 G4double G4BOptnForceCommonTruncatedExp::GetTriggeredProcessXSfraction()
106 else << 82 {
107 { << 83 return fCrossSections[fProcessToApply] / fTotalCrossSection;
108 forceFinalState = true; << 84 }
109 fDummyParticleChange.Initialize( *track ); << 85
110 return &fDummyParticleChange; << 86 void G4BOptnForceCommonTruncatedExp::PostStepInteractionOccured( const G4VProcess* )
111 } << 87 {
112 } << 88 }
113 << 89
114 void G4BOptnForceCommonTruncatedExp:: << 90 void G4BOptnForceCommonTruncatedExp::AddCrossSection( const G4VProcess* process, G4double crossSection )
115 AddCrossSection( const G4VProcess* process, G4 << 91 {
116 { << 92 fTotalCrossSection += crossSection;
117 fTotalCrossSection += crossSection; << 93 fCrossSections[process] = crossSection;
118 fCrossSections[process] = crossSection; << 94 fNumberOfSharing = fCrossSections.size();
119 fNumberOfSharing = fCrossSections.si << 95 fCommonTruncatedExpLaw->SetNumberOfSharing ( fNumberOfSharing );
>> 96 fCommonTruncatedExpLaw->SetForceCrossSection( fTotalCrossSection );
120 } 97 }
121 98
>> 99
122 void G4BOptnForceCommonTruncatedExp::Initializ 100 void G4BOptnForceCommonTruncatedExp::Initialize( const G4Track* track )
123 { 101 {
124 fCrossSections.clear(); 102 fCrossSections.clear();
125 fTotalCrossSection = 0.0; 103 fTotalCrossSection = 0.0;
126 fNumberOfSharing = 0; 104 fNumberOfSharing = 0;
127 fProcessToApply = 0; 105 fProcessToApply = 0;
128 fInteractionOccured = false; <<
129 fInitialMomentum = track->GetMomentum(); 106 fInitialMomentum = track->GetMomentum();
130 107
>> 108 fCommonTruncatedExpLaw->reset();
131 G4VSolid* currentSolid = track->GetVolume()- 109 G4VSolid* currentSolid = track->GetVolume()->GetLogicalVolume()->GetSolid();
132 G4ThreeVector localPosition = (G4Transporta 110 G4ThreeVector localPosition = (G4TransportationManager::GetTransportationManager()->
133 GetNavigator << 111 GetNavigatorForTracking()->
134 GetGlobalToL << 112 GetGlobalToLocalTransform()).TransformPoint(track->GetPosition());
135 G4ThreeVector localDirection = (G4Transporta 113 G4ThreeVector localDirection = (G4TransportationManager::GetTransportationManager()->
136 GetNavigator << 114 GetNavigatorForTracking()->
137 GetGlobalToL << 115 GetGlobalToLocalTransform()).TransformAxis(track->GetMomentumDirection());
138 fMaximumDistance = currentSolid->DistanceToO 116 fMaximumDistance = currentSolid->DistanceToOut(localPosition, localDirection);
139 if ( fMaximumDistance <= DBL_MIN ) { fMaximu << 117 if ( fMaximumDistance <= DBL_MIN ) fMaximumDistance = 0.0;
140 fCommonTruncatedExpLaw->SetMaximumDistance( 118 fCommonTruncatedExpLaw->SetMaximumDistance( fMaximumDistance );
141 } 119 }
142 120
>> 121
143 void G4BOptnForceCommonTruncatedExp::UpdateFor 122 void G4BOptnForceCommonTruncatedExp::UpdateForStep( const G4Step* step )
144 { 123 {
145 fCrossSections.clear(); 124 fCrossSections.clear();
146 fTotalCrossSection = 0.0; 125 fTotalCrossSection = 0.0;
147 fNumberOfSharing = 0; 126 fNumberOfSharing = 0;
148 fProcessToApply = 0; 127 fProcessToApply = 0;
149 << 128
150 fCommonTruncatedExpLaw->UpdateForStep( step- 129 fCommonTruncatedExpLaw->UpdateForStep( step->GetStepLength() );
151 fMaximumDistance = fCommonTruncatedExpLaw->G 130 fMaximumDistance = fCommonTruncatedExpLaw->GetMaximumDistance();
152 } 131 }
153 132
154 void G4BOptnForceCommonTruncatedExp::Sample() 133 void G4BOptnForceCommonTruncatedExp::Sample()
155 { 134 {
156 fCommonTruncatedExpLaw->SetForceCrossSection << 135 fCommonTruncatedExpLaw->reset();
157 fCommonTruncatedExpLaw->Sample(); 136 fCommonTruncatedExpLaw->Sample();
158 ChooseProcessToApply(); << 137 ChooseProcessToApply();
159 fCommonTruncatedExpLaw->SetSelectedProcessXS <<
160 } 138 }
161 139
>> 140
162 void G4BOptnForceCommonTruncatedExp::ChoosePro 141 void G4BOptnForceCommonTruncatedExp::ChooseProcessToApply()
163 { 142 {
164 G4double sigmaRand = G4UniformRand() * fTo << 143 G4double sigmaRand = G4UniformRand() * fTotalCrossSection;
165 G4double sigmaSelect = 0.0; 144 G4double sigmaSelect = 0.0;
166 for ( auto it = fCrossSections.cbegin(); it << 145 for ( std::map< const G4VProcess*, G4double>::const_iterator it = fCrossSections.begin();
167 { << 146 it != fCrossSections.end();
168 sigmaSelect += (*it).second; << 147 it++)
169 if ( sigmaRand <= sigmaSelect ) <<
170 { 148 {
171 fProcessToApply = (*it).first; << 149 sigmaSelect += (*it).second;
172 break; << 150 if ( sigmaRand <= sigmaSelect )
>> 151 {
>> 152 fProcessToApply = (*it).first;
>> 153 break;
>> 154 }
173 } 155 }
174 } <<
175 } 156 }
176 157