Geant4 11.4.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4MicroElecSurface Class Reference

#include <G4MicroElecSurface.hh>

Inheritance diagram for G4MicroElecSurface:

Public Member Functions

 G4MicroElecSurface (const G4String &processName="MicroElecSurface", G4ProcessType type=fElectromagnetic)
 ~G4MicroElecSurface () override
G4bool IsApplicable (const G4ParticleDefinition &aParticleType) override
void SetFlagFranchissement ()
G4double GetMeanFreePath (const G4Track &, G4double, G4ForceCondition *condition) override
G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep) override
void BuildPhysicsTable (const G4ParticleDefinition &) override
G4MicroElecSurfaceStatus GetStatus () const
 G4MicroElecSurface (const G4MicroElecSurface &right)=delete
G4MicroElecSurfaceoperator= (const G4MicroElecSurface &right)=delete
void Initialise ()
Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &aName, G4ProcessType aType=fNotDefined)
 G4VDiscreteProcess (G4VDiscreteProcess &)
virtual ~G4VDiscreteProcess ()
G4VDiscreteProcessoperator= (const G4VDiscreteProcess &)=delete
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
virtual G4double GetCrossSection (const G4double, const G4MaterialCutsCouple *)
virtual G4double MinPrimaryEnergy (const G4ParticleDefinition *, const G4Material *)
Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 G4VProcess (const G4VProcess &right)
virtual ~G4VProcess ()
G4VProcessoperator= (const G4VProcess &)=delete
G4bool operator== (const G4VProcess &right) const
G4bool operator!= (const G4VProcess &right) const
G4double GetCurrentInteractionLength () const
void SetPILfactor (G4double value)
G4double GetPILfactor () const
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
const G4StringGetProcessName () const
G4ProcessType GetProcessType () const
void SetProcessType (G4ProcessType)
G4int GetProcessSubType () const
void SetProcessSubType (G4int)
virtual const G4VProcessGetCreatorProcess () const
virtual void StartTracking (G4Track *)
virtual void EndTracking ()
virtual void SetProcessManager (const G4ProcessManager *)
virtual const G4ProcessManagerGetProcessManager ()
virtual void ResetNumberOfInteractionLengthLeft ()
G4double GetNumberOfInteractionLengthLeft () const
G4double GetTotalNumberOfInteractionLengthTraversed () const
G4bool isAtRestDoItIsEnabled () const
G4bool isAlongStepDoItIsEnabled () const
G4bool isPostStepDoItIsEnabled () const
virtual void DumpInfo () const
virtual void ProcessDescription (std::ostream &outfile) const
void SetVerboseLevel (G4int value)
G4int GetVerboseLevel () const
virtual void SetMasterProcess (G4VProcess *masterP)
const G4VProcessGetMasterProcess () const
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)

Additional Inherited Members

Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double prevStepSize)
void ClearNumberOfInteractionLengthLeft ()
Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager = nullptr
G4VParticleChangepParticleChange = nullptr
G4ParticleChange aParticleChange
G4double theNumberOfInteractionLengthLeft = -1.0
G4double currentInteractionLength = -1.0
G4double theInitialNumberOfInteractionLength = -1.0
G4String theProcessName
G4String thePhysicsTableFileName
G4ProcessType theProcessType = fNotDefined
G4int theProcessSubType = -1
G4double thePILfactor = 1.0
G4int verboseLevel = 0
G4bool enableAtRestDoIt = true
G4bool enableAlongStepDoIt = true
G4bool enablePostStepDoIt = true

Detailed Description

Definition at line 90 of file G4MicroElecSurface.hh.

Constructor & Destructor Documentation

◆ G4MicroElecSurface() [1/2]

G4MicroElecSurface::G4MicroElecSurface ( const G4String & processName = "MicroElecSurface",
G4ProcessType type = fElectromagnetic )
explicit

Definition at line 59 of file G4MicroElecSurface.cc.

60 : G4VDiscreteProcess(processName, type),
61 oldMomentum(0.,0.,0.), previousMomentum(0.,0.,0.),
62 theGlobalNormal(0.,0.,0.), theFacetNormal(0.,0.,0.)
63{
64 if ( verboseLevel > 0)
65 {
66 G4cout << GetProcessName() << " is created " << G4endl;
67 }
68
69 isInitialised=false;
71
72 theStatus = UndefinedSurf;
73 material1 = nullptr;
74 material2 = nullptr;
75
77 theParticleMomentum = 0.;
78
79 flag_franchissement_surface = false;
80 flag_normal = false;
81 flag_reflexion = false;
82 teleportToDo = teleportDone = false;
83
84 ekint = thetat = thetaft = energyThreshold = crossingProbability = 0.0;
85}
@ fSurfaceReflection
@ UndefinedSurf
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout
G4double GetSurfaceTolerance() const
static G4GeometryTolerance * GetInstance()
G4VDiscreteProcess(const G4String &aName, G4ProcessType aType=fNotDefined)
G4int verboseLevel
void SetProcessSubType(G4int)
const G4String & GetProcessName() const

Referenced by G4MicroElecSurface(), and operator=().

◆ ~G4MicroElecSurface()

G4MicroElecSurface::~G4MicroElecSurface ( )
override

Definition at line 89 of file G4MicroElecSurface.cc.

90{}

◆ G4MicroElecSurface() [2/2]

G4MicroElecSurface::G4MicroElecSurface ( const G4MicroElecSurface & right)
delete

Member Function Documentation

◆ BuildPhysicsTable()

void G4MicroElecSurface::BuildPhysicsTable ( const G4ParticleDefinition & )
overridevirtual

Reimplemented from G4VProcess.

Definition at line 132 of file G4MicroElecSurface.cc.

133{
134
135 G4ProductionCutsTable* theCoupleTable =
137 G4int numOfCouples = (G4int)theCoupleTable->GetTableSize();
138 G4cout << "G4MicroElecSurface::Initialise: Ncouples= "
139 << numOfCouples << G4endl;
140
141 for (G4int i = 0; i < numOfCouples; ++i)
142 {
143 const G4Material* material =
144 theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial();
145
146 G4cout << "G4Surface, Material " << i + 1 << " / "
147 << numOfCouples << " : " << material->GetName() << G4endl;
148 if (material->GetName() == "Vacuum")
149 {
150 tableWF[material->GetName()] = 0;
151 continue;
152 }
153 G4String mat = material->GetName();
154 G4MicroElecMaterialStructure str = G4MicroElecMaterialStructure(mat);
155 tableWF[mat] = str.GetWorkFunction();
156 }
157 isInitialised = true;
158}
int G4int
Definition G4Types.hh:85
const G4Material * GetMaterial() const
const G4String & GetName() const
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
std::size_t GetTableSize() const
static G4ProductionCutsTable * GetProductionCutsTable()

◆ GetMeanFreePath()

G4double G4MicroElecSurface::GetMeanFreePath ( const G4Track & ,
G4double ,
G4ForceCondition * condition )
overridevirtual

Implements G4VDiscreteProcess.

Definition at line 446 of file G4MicroElecSurface.cc.

448{
449 *condition = Forced;
450 return DBL_MAX;
451}
G4double condition(const G4ErrorSymMatrix &m)
@ Forced
#define DBL_MAX
Definition templates.hh:62

◆ GetStatus()

G4MicroElecSurfaceStatus G4MicroElecSurface::GetStatus ( ) const

Definition at line 531 of file G4MicroElecSurface.cc.

532{
533 return theStatus;
534}

◆ Initialise()

void G4MicroElecSurface::Initialise ( )

Definition at line 102 of file G4MicroElecSurface.cc.

103{
104 if (isInitialised) { return; }
105
106 G4ProductionCutsTable* theCoupleTable =
108 G4int numOfCouples = (G4int)theCoupleTable->GetTableSize();
109 G4cout << numOfCouples << G4endl;
110
111 for (G4int i = 0; i < numOfCouples; ++i)
112 {
113 const G4Material* material =
114 theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial();
115
116 G4cout << this->GetProcessName() << ", Material " << i + 1
117 << " / " << numOfCouples << " : " << material->GetName() << G4endl;
118 if (material->GetName() == "Vacuum")
119 {
120 tableWF[material->GetName()] = 0; continue;
121 }
122 G4String mat = material->GetName();
123 G4MicroElecMaterialStructure str = G4MicroElecMaterialStructure(mat);
124 tableWF[mat] = str.GetWorkFunction();
125 }
126
127 isInitialised = true;
128}

Referenced by PostStepDoIt().

◆ IsApplicable()

G4bool G4MicroElecSurface::IsApplicable ( const G4ParticleDefinition & aParticleType)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 95 of file G4MicroElecSurface.cc.

96{
97 return ( aParticleType.GetPDGEncoding() == 11 );
98}

◆ operator=()

G4MicroElecSurface & G4MicroElecSurface::operator= ( const G4MicroElecSurface & right)
delete

◆ PostStepDoIt()

G4VParticleChange * G4MicroElecSurface::PostStepDoIt ( const G4Track & aTrack,
const G4Step & aStep )
overridevirtual

Reimplemented from G4VDiscreteProcess.

Definition at line 162 of file G4MicroElecSurface.cc.

163{
164
165 if (!isInitialised) { Initialise(); }
166 theStatus = UndefinedSurf;
167
168 // Definition of the parameters for the particle
169 aParticleChange.Initialize(aTrack);
170 aParticleChange.ProposeVelocity(aTrack.GetVelocity());
171
172 G4StepPoint* pPreStepPoint = aStep.GetPreStepPoint();
173 G4StepPoint* pPostStepPoint = aStep.GetPostStepPoint();
174
175 material1 = pPreStepPoint -> GetMaterial();
176 material2 = pPostStepPoint -> GetMaterial();
177
178 theStatus = UndefinedSurf;
179
180 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
181
182 theParticleMomentum = aParticle->GetTotalMomentum();
183 previousMomentum = oldMomentum;
184 oldMomentum = aParticle->GetMomentumDirection();
185
186 // Fisrt case: not a boundary
187 if (pPostStepPoint->GetStepStatus() != fGeomBoundary
188 || pPostStepPoint->GetPhysicalVolume()->GetName() == pPreStepPoint->GetPhysicalVolume()->GetName())
189 {
190 theStatus = NotAtBoundarySurf;
191 flag_franchissement_surface = false;
192 flag_reflexion = false;
193 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
194 }
195 theStatus = UndefinedSurf;
196
197 // Third case: same material
198 if (material1 == material2)
199 {
200 theStatus = SameMaterialSurf;
201 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
202 }
203 if (verboseLevel > 3)
204 {
205 G4cout << G4endl << " Electron at Boundary! " << G4endl;
206 G4VPhysicalVolume* thePrePV = pPreStepPoint->GetPhysicalVolume();
207 G4VPhysicalVolume* thePostPV = pPostStepPoint->GetPhysicalVolume();
208 if (thePrePV) G4cout << " thePrePV: " << thePrePV->GetName() << G4endl;
209 if (thePostPV) G4cout << " thePostPV: " << thePostPV->GetName() << G4endl;
210 G4cout << " Old Momentum Direction: " << oldMomentum << G4endl;
211 }
212
213 // Definition of the parameters for the surface
214 G4ThreeVector theGlobalPoint = pPostStepPoint->GetPosition();
215
216 G4Navigator* theNavigator = G4TransportationManager::
218
219 G4bool valid;
220 theGlobalNormal = theNavigator->GetGlobalExitNormal(theGlobalPoint, &valid);
221
222 if (valid)
223 {
224 theGlobalNormal = -theGlobalNormal;
225 }
226 else
227 {
229 ed << " G4MicroElecSurface/PostStepDoIt(): "
230 << " The Navigator reports that it returned an invalid normal"
231 << G4endl;
232 G4Exception("G4MuElecSurf::PostStepDoIt", "OpBoun01",
234 "Invalid Surface Normal - Geometry must return valid surface normal");
235 }
236
237 // Exception: the particle is not in the right direction
238 if (oldMomentum * theGlobalNormal > 0.0)
239 {
240 theGlobalNormal = -theGlobalNormal;
241 }
242
243 if (aTrack.GetStepLength()<=kCarTolerance/2 * 0.0000000001)
244 {
245 if (flag_reflexion == true)
246 {
247 flag_reflexion = false;
248 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
249 }
250
251 theStatus = StepTooSmallSurf;
252
253 G4double energyThreshold_surface = 0.0*eV;
254
255 WorkFunctionTable::iterator postStepWF;
256 postStepWF = tableWF.find(pPostStepPoint->GetMaterial()->GetName());
257 WorkFunctionTable::iterator preStepWF;
258 preStepWF = tableWF.find(pPreStepPoint->GetMaterial()->GetName());
259
260 if (postStepWF == tableWF.end())
261 {
262 G4String str = "Material ";
263 str += pPostStepPoint->GetMaterial()->GetName() + " not found!";
264 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
265 return nullptr;
266 }
267 else if (preStepWF == tableWF.end())
268 {
269 G4String str = "Material ";
270 str += pPreStepPoint->GetMaterial()->GetName() + " not found!";
271 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
272 return nullptr;
273 }
274 else
275 {
276 G4double thresholdNew_surface = postStepWF->second;
277 G4double thresholdOld_surface = preStepWF->second;
278 energyThreshold_surface = thresholdNew_surface - thresholdOld_surface;
279 }
280
281 if (flag_franchissement_surface == true)
282 {
283 aParticleChange.ProposeEnergy(aStep.GetPostStepPoint()->GetKineticEnergy() + energyThreshold_surface);
284 flag_franchissement_surface = false;
285 }
286 if (flag_reflexion == true && flag_normal == true)
287 {
288 aParticleChange.ProposeMomentumDirection(-Reflexion(aStep.GetPostStepPoint()));
289 flag_reflexion = false;
290 flag_normal = false;
291 }
292 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
293 }
294
295 flag_normal = (theGlobalNormal == G4ThreeVector(0, 0, 1)
296 || theGlobalNormal == G4ThreeVector(0, 0, -1));
297
298 G4LogicalSurface* Surface = nullptr;
299
301 (pPreStepPoint ->GetPhysicalVolume(),
302 pPostStepPoint->GetPhysicalVolume());
303
304 if (Surface == nullptr)
305 {
306 G4bool enteredDaughter=(pPostStepPoint->GetPhysicalVolume()->GetMotherLogical()
307 == pPreStepPoint->GetPhysicalVolume()->GetLogicalVolume());
308 if(enteredDaughter)
309 {
311 (pPostStepPoint->GetPhysicalVolume()->GetLogicalVolume());
312
313 if(Surface == nullptr)
315 (pPreStepPoint->GetPhysicalVolume()->GetLogicalVolume());
316 }
317 else
318 {
320 (pPreStepPoint->GetPhysicalVolume()->GetLogicalVolume());
321
322 if(Surface == nullptr)
324 (pPostStepPoint->GetPhysicalVolume()->GetLogicalVolume());
325 }
326 }
327
328 // Definition of the parameters for the surface crossing
329 G4VPhysicalVolume* thePrePV = pPreStepPoint->GetPhysicalVolume();
330 G4VPhysicalVolume* thePostPV = pPostStepPoint->GetPhysicalVolume();
331
332 energyThreshold = 0.0*eV;
333 G4double energyDelta = 0;
334
335 if ((thePrePV)&&(thePostPV))
336 {
337 WorkFunctionTable::iterator postStepWF;
338 postStepWF = tableWF.find(thePostPV->GetLogicalVolume()->GetMaterial()->GetName());
339 WorkFunctionTable::iterator preStepWF;
340 preStepWF = tableWF.find(thePrePV->GetLogicalVolume()->GetMaterial()->GetName());
341
342 if (postStepWF == tableWF.end())
343 {
344 G4String str = "Material ";
345 str += thePostPV->GetLogicalVolume()->GetMaterial()->GetName() + " not found!";
346 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
347 return nullptr;
348 }
349 else if (preStepWF == tableWF.end())
350 {
351 G4String str = "Material ";
352 str += thePrePV->GetLogicalVolume()->GetMaterial()->GetName() + " not found!";
353 G4Exception("G4Surface::G4Surface", "em0002", FatalException, str);
354 return nullptr;
355 }
356 else
357 {
358 G4double thresholdNew = postStepWF->second;
359 G4double thresholdOld = preStepWF->second;
360
361 energyThreshold = thresholdNew - thresholdOld;
362 energyDelta = thresholdOld- thresholdNew;
363 }
364 }
365
366 ekint=aStep.GetPostStepPoint()->GetKineticEnergy();
367 thetat= GetIncidentAngle(); //angle d'incidence
368 G4double ekinNormalt=ekint*std::cos(thetat)*std::cos(thetat);
369 G4double sin_thetaft = std::sqrt(ekint/(ekint+energyThreshold))*std::sin(thetat);
370 G4double cos_thetaft = 0.0;
371 //Refraction angle
372
373 if (1.0-sin_thetaft*sin_thetaft>0) {
374 cos_thetaft = std::sqrt(1.0-sin_thetaft*sin_thetaft);
375 }
376 else {
377 cos_thetaft = 0.0;
378 }
379
380 G4double aleat=G4UniformRand();
381 G4double waveVectort=std::sqrt(2*9.1093826E-31*1.602176487E-19)/(6.6260755E-34/(2.0*pi));
382
383 // Parameter for an exponential barrier of potential (Thesis P68)
384 G4double at=0.5E-10;
385
386 crossingProbability=0;
387
388 G4double kft=waveVectort*std::sqrt(ekint+energyThreshold)*cos_thetaft;
389 G4double kit=waveVectort*std::sqrt(ekinNormalt);
390
391 crossingProbability=1-(std::pow(std::sinh(pi*at*(kit-kft)), 2.0)/std::pow(std::sinh(pi*at*(kit+kft)), 2.0));
392
393 // First case: the electron crosses the surface
394 if((aleat<=crossingProbability)&&(ekint> energyDelta))
395 {
396 if (aStep.GetPreStepPoint()->GetMaterial()->GetName()
397 != aStep.GetPostStepPoint()->GetMaterial()->GetName())
398 {
399 aParticleChange.ProposeEnergy(ekint - energyDelta);
400 flag_franchissement_surface = true;
401 }
402
403 // calculation of cos_thetaft for thetaft=std::abs(thetaft-thetat);
404 cos_thetaft = cos_thetaft*std::cos(thetat)+sin_thetaft*std::sin(thetat);
405
407 G4ThreeVector xVerst = zVerst.orthogonal();
408 G4ThreeVector yVerst = zVerst.cross(xVerst);
409
410 G4double xDirt = std::sqrt(1. - cos_thetaft*cos_thetaft);
411 G4double yDirt = xDirt;
412
413 G4ThreeVector zPrimeVerst=((xDirt*xVerst + yDirt*yVerst + cos_thetaft*zVerst));
414
415 aParticleChange.ProposeMomentumDirection(zPrimeVerst.unit());
416 }
417 else if ((aleat > crossingProbability) && (ekint> energyDelta))
418 {
419 flag_reflexion = true;
420 if (flag_normal)
421 {
422 aParticleChange.ProposeMomentumDirection(-oldMomentum.unit());
423 }
424 else
425 {
426 aParticleChange.ProposeMomentumDirection(Reflexion(aStep.GetPostStepPoint()));
427 }
428 }
429 else
430 {
431 if (flag_normal)
432 {
433 aParticleChange.ProposeMomentumDirection(-oldMomentum.unit());
434 }
435 else
436 {
437 aParticleChange.ProposeMomentumDirection(Reflexion(aStep.GetPostStepPoint()));
438 }
439 flag_reflexion = true;
440 }
441 return G4VDiscreteProcess::PostStepDoIt(aTrack, aStep);
442}
@ FatalException
@ EventMustBeAborted
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
std::ostringstream G4ExceptionDescription
@ StepTooSmallSurf
@ SameMaterialSurf
@ NotAtBoundarySurf
@ fGeomBoundary
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition G4Types.hh:83
bool G4bool
Definition G4Types.hh:86
#define G4UniformRand()
Definition Randomize.hh:52
Hep3Vector unit() const
Hep3Vector orthogonal() const
Hep3Vector cross(const Hep3Vector &) const
const G4ThreeVector & GetMomentumDirection() const
G4double GetTotalMomentum() const
static G4LogicalBorderSurface * GetSurface(const G4VPhysicalVolume *vol1, const G4VPhysicalVolume *vol2)
static G4LogicalSkinSurface * GetSurface(const G4LogicalVolume *vol)
G4Material * GetMaterial() const
virtual G4ThreeVector GetGlobalExitNormal(const G4ThreeVector &point, G4bool *valid)
G4StepStatus GetStepStatus() const
G4Material * GetMaterial() const
const G4ThreeVector & GetPosition() const
const G4ThreeVector & GetMomentumDirection() const
G4VPhysicalVolume * GetPhysicalVolume() const
G4double GetKineticEnergy() const
G4StepPoint * GetPreStepPoint() const
G4StepPoint * GetPostStepPoint() const
G4double GetVelocity() const
const G4DynamicParticle * GetDynamicParticle() const
G4double GetStepLength() const
static G4TransportationManager * GetTransportationManager()
G4Navigator * GetNavigatorForTracking() const
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
G4LogicalVolume * GetLogicalVolume() const
const G4String & GetName() const
G4ParticleChange aParticleChange
const G4double pi

◆ SetFlagFranchissement()

void G4MicroElecSurface::SetFlagFranchissement ( )

Definition at line 539 of file G4MicroElecSurface.cc.

540{
541 flag_franchissement_surface = false;
542}

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