Geant4  10.02.p03
G4BOptnForceFreeFlight Class Reference

#include <G4BOptnForceFreeFlight.hh>

Inheritance diagram for G4BOptnForceFreeFlight:
Collaboration diagram for G4BOptnForceFreeFlight:

Public Member Functions

 G4BOptnForceFreeFlight (G4String name)
 
virtual ~G4BOptnForceFreeFlight ()
 
virtual const G4VBiasingInteractionLawProvideOccurenceBiasingInteractionLaw (const G4BiasingProcessInterface *, G4ForceCondition &)
 
virtual G4bool DenyProcessPostStepDoIt (const G4BiasingProcessInterface *, const G4Track *, const G4Step *, G4double &)
 
virtual void AlongMoveBy (const G4BiasingProcessInterface *, const G4Step *, G4double)
 
virtual G4VParticleChange * ApplyFinalStateBiasing (const G4BiasingProcessInterface *, const G4Track *, const G4Step *, G4bool &)
 
virtual G4double DistanceToApplyOperation (const G4Track *, G4double, G4ForceCondition *)
 
virtual G4VParticleChange * GenerateBiasingFinalState (const G4Track *, const G4Step *)
 
G4ILawForceFreeFlightGetForceFreeFlightLaw ()
 
void ResetInitialTrackWeight (G4double w)
 
G4bool OperationComplete () const
 
- Public Member Functions inherited from G4VBiasingOperation
 G4VBiasingOperation (G4String name)
 
virtual ~G4VBiasingOperation ()
 
virtual G4double ProposeAlongStepLimit (const G4BiasingProcessInterface *)
 
virtual G4GPILSelection ProposeGPILSelection (const G4GPILSelection wrappedProcessSelection)
 
const G4StringGetName () const
 
std::size_t GetUniqueID () const
 

Private Attributes

G4ILawForceFreeFlightfForceFreeFlightInteractionLaw
 
G4double fCumulatedWeightChange
 
G4double fInitialTrackWeight
 
G4ParticleChange fParticleChange
 
G4bool fOperationComplete
 

Detailed Description

Definition at line 56 of file G4BOptnForceFreeFlight.hh.

Constructor & Destructor Documentation

◆ G4BOptnForceFreeFlight()

G4BOptnForceFreeFlight::G4BOptnForceFreeFlight ( G4String  name)

Definition at line 33 of file G4BOptnForceFreeFlight.cc.

34  : G4VBiasingOperation ( name ),
35  fCumulatedWeightChange ( -1.0 ),
36  fInitialTrackWeight ( -1.0 ),
37  fOperationComplete ( true )
38 {
39  fForceFreeFlightInteractionLaw = new G4ILawForceFreeFlight("LawForOperation"+name);
40 }
G4ILawForceFreeFlight * fForceFreeFlightInteractionLaw
G4VBiasingOperation(G4String name)

◆ ~G4BOptnForceFreeFlight()

G4BOptnForceFreeFlight::~G4BOptnForceFreeFlight ( )
virtual

Definition at line 42 of file G4BOptnForceFreeFlight.cc.

43 {
45 }
G4ILawForceFreeFlight * fForceFreeFlightInteractionLaw

Member Function Documentation

◆ AlongMoveBy()

void G4BOptnForceFreeFlight::AlongMoveBy ( const G4BiasingProcessInterface ,
const G4Step *  ,
G4double  weightChange 
)
virtual

Reimplemented from G4VBiasingOperation.

Definition at line 131 of file G4BOptnForceFreeFlight.cc.

132 {
133  fCumulatedWeightChange *= weightChange;
134 }

◆ ApplyFinalStateBiasing()

G4VParticleChange * G4BOptnForceFreeFlight::ApplyFinalStateBiasing ( const G4BiasingProcessInterface callingProcess,
const G4Track *  track,
const G4Step *  step,
G4bool forceFinalState 
)
virtual

Implements G4VBiasingOperation.

Definition at line 85 of file G4BOptnForceFreeFlight.cc.

89 {
90  // -- If the track is reaching the volume boundary, its free flight ends. In this case, its zero
91  // -- weight is brought back to non-zero value: its initial weight is restored by the first
92  // -- ApplyFinalStateBiasing operation called, and the weight for force free flight is applied
93  // -- is applied by each operation.
94  // -- If the track is not reaching the volume boundary, it zero weight flight continues.
95 
96  fParticleChange.Initialize( *track );
97  forceFinalState = true;
98  if ( step->GetPostStepPoint()->GetStepStatus() == fGeomBoundary )
99  {
100  // -- Sanity checks:
101  if ( fInitialTrackWeight <= DBL_MIN )
102  {
104  ed << " Initial track weight is null ! " << G4endl;
105  G4Exception(" G4BOptnForceFreeFlight::ApplyFinalStateBiasing(...)",
106  "BIAS.GEN.05",
107  JustWarning,
108  ed);
109  }
111  {
113  ed << " Cumulated weight is null ! " << G4endl;
114  G4Exception(" G4BOptnForceFreeFlight::ApplyFinalStateBiasing(...)",
115  "BIAS.GEN.06",
116  JustWarning,
117  ed);
118  }
119 
120  G4double proposedWeight = track->GetWeight();
121  if ( callingProcess->GetIsFirstPostStepDoItInterface() ) proposedWeight = fCumulatedWeightChange * fInitialTrackWeight;
122  else proposedWeight *= fCumulatedWeightChange;
123  fParticleChange.ProposeWeight(proposedWeight);
124  fOperationComplete = true;
125  }
126 
127  return &fParticleChange;
128 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4bool GetIsFirstPostStepDoItInterface(G4bool physOnly=true) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define DBL_MIN
Definition: templates.hh:75
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
Here is the call graph for this function:

◆ DenyProcessPostStepDoIt()

G4bool G4BOptnForceFreeFlight::DenyProcessPostStepDoIt ( const G4BiasingProcessInterface ,
const G4Track *  ,
const G4Step *  step,
G4double proposedWeight 
)
virtual

Definition at line 54 of file G4BOptnForceFreeFlight.cc.

55 {
56  // -- force free flight always deny process to apply its doit.
57  // -- if reaching boundary, track is restored with non-zero weight
59  {
61  ed << " Initial track weight is null ! " << G4endl;
62  G4Exception(" G4BOptnForceFreeFlight::DenyProcessPostStepDoIt(...)",
63  "BIAS.GEN.05",
65  ed);
66  }
68  {
70  ed << " Cumulated weight is null ! " << G4endl;
71  G4Exception(" G4BOptnForceFreeFlight::DenyProcessPostStepDoIt(...)",
72  "BIAS.GEN.06",
74  ed);
75  }
76  if ( step->GetPostStepPoint()->GetStepStatus() == fGeomBoundary )
77  {
78  if ( proposedWeight <= DBL_MIN ) proposedWeight = fCumulatedWeightChange * fInitialTrackWeight;
79  else proposedWeight *= fCumulatedWeightChange;
80  }
81 
82  return true;
83 }
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define DBL_MIN
Definition: templates.hh:75
#define G4endl
Definition: G4ios.hh:61
Here is the call graph for this function:

◆ DistanceToApplyOperation()

virtual G4double G4BOptnForceFreeFlight::DistanceToApplyOperation ( const G4Track *  ,
G4double  ,
G4ForceCondition *   
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 73 of file G4BOptnForceFreeFlight.hh.

75  {return DBL_MAX;}
#define DBL_MAX
Definition: templates.hh:83

◆ GenerateBiasingFinalState()

virtual G4VParticleChange* G4BOptnForceFreeFlight::GenerateBiasingFinalState ( const G4Track *  ,
const G4Step *   
)
inlinevirtual

Implements G4VBiasingOperation.

Definition at line 76 of file G4BOptnForceFreeFlight.hh.

77  {return 0;}

◆ GetForceFreeFlightLaw()

G4ILawForceFreeFlight* G4BOptnForceFreeFlight::GetForceFreeFlightLaw ( )
inline

Definition at line 84 of file G4BOptnForceFreeFlight.hh.

84  {
86  }
G4ILawForceFreeFlight * fForceFreeFlightInteractionLaw

◆ OperationComplete()

G4bool G4BOptnForceFreeFlight::OperationComplete ( ) const
inline

Definition at line 89 of file G4BOptnForceFreeFlight.hh.

◆ ProvideOccurenceBiasingInteractionLaw()

const G4VBiasingInteractionLaw * G4BOptnForceFreeFlight::ProvideOccurenceBiasingInteractionLaw ( const G4BiasingProcessInterface ,
G4ForceCondition &  proposeForceCondition 
)
virtual

Implements G4VBiasingOperation.

Definition at line 47 of file G4BOptnForceFreeFlight.cc.

48 {
49  fOperationComplete = false;
50  proposeForceCondition = Forced;
52 }
G4ILawForceFreeFlight * fForceFreeFlightInteractionLaw

◆ ResetInitialTrackWeight()

void G4BOptnForceFreeFlight::ResetInitialTrackWeight ( G4double  w)
inline

Definition at line 88 of file G4BOptnForceFreeFlight.hh.

Here is the caller graph for this function:

Member Data Documentation

◆ fCumulatedWeightChange

G4double G4BOptnForceFreeFlight::fCumulatedWeightChange
private

Definition at line 93 of file G4BOptnForceFreeFlight.hh.

◆ fForceFreeFlightInteractionLaw

G4ILawForceFreeFlight* G4BOptnForceFreeFlight::fForceFreeFlightInteractionLaw
private

Definition at line 92 of file G4BOptnForceFreeFlight.hh.

◆ fInitialTrackWeight

G4double G4BOptnForceFreeFlight::fInitialTrackWeight
private

Definition at line 93 of file G4BOptnForceFreeFlight.hh.

◆ fOperationComplete

G4bool G4BOptnForceFreeFlight::fOperationComplete
private

Definition at line 96 of file G4BOptnForceFreeFlight.hh.

◆ fParticleChange

G4ParticleChange G4BOptnForceFreeFlight::fParticleChange
private

Definition at line 95 of file G4BOptnForceFreeFlight.hh.


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