This section presents the generic biasing classes which have been introduced since release 10.0. Theses classes are meant to virtually allow any type of process-level based biasing.
In 10.0 and 10.1, only PostStep
biasing actions are possible,
which allows:
PostStepGetPhysicalInteractionLength(...)
level;PostStepDoIt(...)
We have decided to split the actual biasing actions (change probability occurence of a process, change its final state generation, split, kill, etc...) from the decisions about these actions to be taken. The reason for this is that several biasing actions are often needed to build one biasing technique, these actions having to be selected along some specific logic of the technique.
For example, a technique like the "force collision" of MCNP involves a splitting, a force interaction of one copy in the volume and a force free flight (under zero weight) of the other copy in the volume.
The classes which provides the interfaces for these biasing actions and decisions are respectively:
G4VBiasingOperation
G4VBiasingOperator
A third class,
G4BiasingProcessInterface
provides the interface between these classes and the stepping.
G4VBiasingOperation
defines two types of methods:
virtual const G4VBiasingInteractionLaw* ProvideOccurenceBiasingInteractionLaw( const G4BiasingProcessInterface* /* callingProcess */,
G4ForceCondition& /* proposeForceCondition */ ) = 0;
virtual G4VParticleChange* ApplyFinalStateBiasing( const G4BiasingProcessInterface* /* callingProcess */,
const G4Track* /* track */,
const G4Step* /* step */,
G4bool& /* forceBiasedFinalState */) = 0;
G4BiasingProcessInterface callingProcess
) and has to return a biased particle
change. This one should take care of providing the proper weight values for the biasing applied.
virtual G4double DistanceToApplyOperation( const G4Track* /* track */,
G4double /* previousStepSize */,
G4ForceCondition* /* condition */) = 0;
G4BiasingProcessInterface callingProcess
to make the biasing opetation to compete for limiting the step.
virtual G4VParticleChange* GenerateBiasingFinalState( const G4Track* /* track */,
const G4Step* /* step */) = 0;
The G4VBiasingOperator
class is meant to define the interface to pilot biasing operations. It selects biasing operations or sequences of biasing
operations to build up the logic of a specific biasing technique. It is messaged by each G4BiasingProcessInterface callingProcess
instance, which
pass their "identity" through their own pointer values. The operator has hence all the necessary information for taking decisions.
In addition, since 10.1, the G4BiasingProcessInterface callingProcess
instances do "anticipated" calls to their underneath wrapped physics process (if any) :
the first of the G4BiasingProcessInterface callingProcess
instance, trigger calls
to all PostStepGetPhysicalInteractionLength(...)
methods of biased processes. In this way, all process cross-sections have been updated before the
first call to the biasing operator. This one has hence all ready-to-use physics information the first it is messaged in the step by the first
G4BiasingProcessInterface callingProcess
instance.
The G4VBiasingOperator
class defines the following interface:
virtual G4VBiasingOperation* ProposeNonPhysicsBiasingOperation( const G4Track* track, const G4BiasingProcessInterface* callingProcess ) = 0;
G4BiasingProcessInterface
instance which does not hold/wrap a physics process.
G4VBiasingOperation
returned will have its DistanceToApplyOperation(...)
and
GenerateBiasingFinalState(...)
methods called at proper times (post step GPIL and post step DoIt times, respectively)
by the G4BiasingProcessInterface callingProcess
.
virtual G4VBiasingOperation* ProposeOccurenceBiasingOperation( const G4Track* track, const G4BiasingProcessInterface* callingProcess ) = 0;
G4VBiasingOperation
returned will have its ProvideOccurenceBiasingInteractionLaw(...)
called to get the biased
interaction law to be used.
PostStepGetPhysicalInteractionLength(...)
level.
virtual G4VBiasingOperation* ProposeFinalStateBiasingOperation( const G4Track* track, const G4BiasingProcessInterface* callingProcess ) = 0;
G4VBiasingOperation
returned will have its ApplyFinalStateBiasing(...)
called to generated the physics-based
biased final state.
PostStepDoIt(...)
level.
The occurence biasing is the biasing which affects the probability for a physics process to occur. For now, we discuss only the case of neutral particles, ie, having no continuous energy loss along a step. The case of charged particles is expected to be treated in later releases.
The process occurence is driven by the exponential law, which parameter in Geant4 is the process mean free path λ, which is also the inverse of the cross-section σ. The (analog or natural) probability density function (pdf) of interactions is given by pa(ℓ)=σa·exp(-σaℓ), where ℓ is the distance at which the interaction occurs, and where we have relabelled σ as σa. In Geant4, volumes are made of a single material, meaning that σa does not depend on ℓ : at some position, starting point of a step, the track "sees" the same cross-section at all positions ℓ in the volume.
The occurence biasing consists in substituting pa(ℓ) by some arbitrary biased interaction law pb(ℓ).
More details about the formalism will be provided in the physics reference manual biasing related part (to come).
An arbitrary interaction law pb(ℓ) can be recasted in term of an "effective" cross-section σb(ℓ), which depends on ℓ in general, as
which is the formalism that corresponds to non-constant cross-sections.
When applying this formalism to the analog pdf pa(ℓ), we simply have σb(ℓ) = σa and PNI,b(ℓ) = exp(-σaℓ), as this must.
In contrast with the analog case, the "effective" cross-section may depend on ℓ, even in a volume made of single material. An example of this is the forced interaction: if it is decided that the interaction of a process, with analog cross-section σa, will be forced and will happen somewhere on the segment [0,L], then we have
In the occurence biasing, two weights have to be taken into account:
G4VBiasingInteractionLaw
is the interface for implementing "interaction laws", it defines the pure virtual methods
virtual G4double ComputeNonInteractionProbabilityAt(G4double length) const = 0;
virtual G4double ComputeEffectiveCrossSectionAt(G4double length) const = 0;
virtual G4double SampleInteractionLength() = 0;
In the case of the occurence biasing, the dedicated virtual const G4VBiasingOperation::ProvideOccurenceBiasingInteractionLaw(const G4BiasingProcessInterface*
callingProcess, G4ForceCondition& proposeForceCondition ) = 0;
is meant to return the biased law. The G4BiasingProcessInterface
will not change the state of the law. It will collect the sampled interaction length (that the biasing operation must have asked to law to sample) and will
used the non-interaction probability method in its AlongStepDoIt(...)
to compute the weight for non-interaction,
these weight being multiplied among biased processes,
and it will use the effective cross-section of process "i", if process "i" wins the interaction length race, in its PostStepDoIt(...)
to compute the weight for interaction.
To compute these weights, the G4BiasingProcessInterface
holds a private interaction law, to which it sets the analog process cross-section that
it collects at the beginning of the step.
As occurence biasing and final state biasing are independent operations, the weight correction for interaction due to the occurence biasing is applied on top of the final state generated by the process (this final state being biased or not).