Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
GB06BOptnSplitAndKillByImportance.cc
이 파일의 문서화 페이지로 가기
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$
27 //
30 
32 #include "Randomize.hh"
33 
34 
38 #include "G4ProcessManager.hh"
39 
40 
41 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
42 
44 : G4VBiasingOperation (name),
45  fParallelWorldIndex ( -1 ),
46  fBiasingSharedData ( nullptr ),
47  fParticleChange(),
48  fDummyParticleChange()
49 {}
50 
51 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
52 
54 {}
55 
56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
57 
60  G4double,
62 {
63 
64  // -- Remember the touchable history (ie geometry state) at the beginning of the step:
65  // -- Start by getting the process handling the step limitation in parallel
66  // -- geometries for the generic biasing:
67  auto biasingLimiterProcess = fBiasingSharedData->GetParallelGeometriesLimiterProcess();
69  biasingLimiterProcess
70  ->GetNavigator( fParallelWorldIndex ) // -- get the navigator of the geometry...
71  ->CreateTouchableHistoryHandle(); // -- ...and collect the geometry state.
72 
73  // -- We request to be "forced" : meaning GenerateBiasingFinalState(...) will be called
74  // -- anyway at the PostStepDoIt(...) stage (ie, when the track will have been moved to
75  // -- its end point position) and, there, we will have to handle the decision to
76  // -- split/kill if we are on a volume boundary, or do nothing, if we're not:
77  *condition = Forced;
78 
79  // -- As we're forced, we make this returned length void:
80  return DBL_MAX;
81 }
82 
83 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
84 
87  const G4Step* )
88 {
89  // -- Given we used the "Forced" condition, this method is always called.
90  // -- We check the status of the step in the parallel geometry, and apply
91  // -- splitting/killing if the step has been limited by the geometry.
92 
93  // -- We first check if we limit the step in the parallel geometry:
94  G4bool isLimiting = fBiasingSharedData
97 
98  // -- if not limiting, we leave the track unchanged:
99  if ( !isLimiting )
100  {
102  return &fDummyParticleChange;
103  }
104 
105  // -- We collect the geometry state at the end point step:
106  // -- Note this is the same call than in the DistanceToApplyOperation(...) for the
107  // -- fPreStepTouchableHistory, but the navigator has pushed the track in the next
108  // -- volume since then (even if the track is still on the boundary), and hence the
109  // -- geometry state has changed.
110  auto biasingLimiterProcess = fBiasingSharedData->GetParallelGeometriesLimiterProcess();
112  biasingLimiterProcess
113  ->GetNavigator( fParallelWorldIndex )
114  ->CreateTouchableHistoryHandle();
115 
116 
117  // -- We verify we are still in the "same" physical volume:
118  // -- remember : using a replica, we have "one" physical volume
119  // -- but representing many different placements, distinguished
120  // -- by replica number. By checking we are in the same physical
121  // -- volume, we verify the end step point has not reached the
122  // -- world volume of the parallel geometry:
123  if ( fPreStepTouchableHistory ->GetVolume() !=
125  {
126  // -- the track is leaving the volumes in which biasing is applied,
127  // -- we leave this track unchanged:
129  return &fDummyParticleChange;
130  }
131 
132  // -------------------------------------------------------------------------------------
133  // -- At this stage, we know we have a particle crossing a boundary between two slices,
134  // -- each of this slice has a well defined importance : we apply the biasing.
135  // -- We will split if the track is entering a volume with higher importance, and
136  // -- kill (applying Russian roulette) in the other case.
137  // -------------------------------------------------------------------------------------
138 
139  // -- We start by getting the replica numbers:
140  G4int preReplicaNumber = fPreStepTouchableHistory->GetReplicaNumber();
141  G4int postReplicaNumber = fPostStepTouchableHistory->GetReplicaNumber();
142 
143  // -- and get the related importance we defined in the importance map:
144  G4int preImportance = (*fImportanceMap)[ preReplicaNumber];
145  G4int postImportance = (*fImportanceMap)[postReplicaNumber];
146 
147 
148  // -- Get track weight:
149  G4double initialWeight = track->GetWeight();
150 
151  // -- Initialize the "particle change" (it will communicate the new track state to
152  // -- the tracking):
153  fParticleChange.Initialize(*track);
154 
155  if ( postImportance > preImportance )
156  {
157  // -- We split :
158  G4int splittingFactor = postImportance/preImportance;
159 
160  // Define the tracks weight:
161  G4double weightOfTrack = initialWeight/splittingFactor;
162 
163  // Ask currect track weight to be changed to the new value:
164  fParticleChange.ProposeParentWeight( weightOfTrack );
165  // Now we clone this track (this is the actual splitting):
166  // we will then have the primary and clone of it, hence the
167  // splitting by a factor 2:
168  G4Track* clone = new G4Track( *track );
169  clone->SetWeight( weightOfTrack );
170  fParticleChange.AddSecondary( clone );
171  // -- Below's call added for safety & illustration : inform particle change to not
172  // -- modify the clone (ie : daughter) weight to make it that of the primary.
173  // -- Here this call is not mandatory as both tracks have same weights.
175  }
176  else
177  {
178  // -- We apply Russian roulette:
179  G4double survivingProbability = G4double(postImportance) / G4double(preImportance);
180 
181  // Shoot a random number (in ]0,1[ segment):
182  G4double random = G4UniformRand();
183  if ( random > survivingProbability )
184  {
185  // We ask for the the track to be killed:
187  }
188  else
189  {
190  // In this case, the track survives. We change its weight
191  // to conserve weight among killed and survival tracks:
192  fParticleChange.ProposeParentWeight( initialWeight/survivingProbability );
193  }
194  }
195 
196  return &fParticleChange;
197 
198 }
199 
200 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
const XML_Char * name
Definition: expat.h:151
const G4BiasingProcessSharedData * fBiasingSharedData
virtual G4double DistanceToApplyOperation(const G4Track *, G4double, G4ForceCondition *condition) final
void SetWeight(G4double aValue)
G4int GetReplicaNumber(G4int depth=0) const
const G4ParallelGeometriesLimiterProcess * GetParallelGeometriesLimiterProcess() const
void SetSecondaryWeightByProcess(G4bool)
G4double condition(const G4ErrorSymMatrix &m)
const std::vector< G4bool > & GetIsLimiting() const
void AddSecondary(G4Track *aSecondary)
G4double GetWeight() const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
virtual void Initialize(const G4Track &track)
virtual void Initialize(const G4Track &)
G4VPhysicalVolume * GetVolume(G4int depth=0) const
virtual G4VParticleChange * GenerateBiasingFinalState(const G4Track *, const G4Step *) final
#define G4UniformRand()
Definition: Randomize.hh:53
Definition: G4Step.hh:76
Definition of the GB06BOptnSplitAndKillByImportance class.
int G4int
Definition: G4Types.hh:78
void ProposeParentWeight(G4double finalWeight)
G4ForceCondition
#define DBL_MAX
Definition: templates.hh:83
void ProposeTrackStatus(G4TrackStatus status)