Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Static Public Member Functions
G4AdjointCSManager Class Reference

#include <G4AdjointCSManager.hh>

Public Member Functions

 ~G4AdjointCSManager ()
 
G4int GetNbProcesses ()
 
size_t RegisterEmAdjointModel (G4VEmAdjointModel *)
 
void RegisterEmProcess (G4VEmProcess *aProcess, G4ParticleDefinition *aPartDef)
 
void RegisterEnergyLossProcess (G4VEnergyLossProcess *aProcess, G4ParticleDefinition *aPartDef)
 
void RegisterAdjointParticle (G4ParticleDefinition *aPartDef)
 
void BuildCrossSectionMatrices ()
 
void BuildTotalSigmaTables ()
 
G4double GetTotalAdjointCS (G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
 
G4double GetTotalForwardCS (G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
 
G4double GetAdjointSigma (G4double Ekin_nuc, size_t index_model, G4bool is_scat_proj_to_proj, const G4MaterialCutsCouple *aCouple)
 
void GetEminForTotalCS (G4ParticleDefinition *aPartDef, const G4MaterialCutsCouple *aCouple, G4double &emin_adj, G4double &emin_fwd)
 
void GetMaxFwdTotalCS (G4ParticleDefinition *aPartDef, const G4MaterialCutsCouple *aCouple, G4double &e_sigma_max, G4double &sigma_max)
 
void GetMaxAdjTotalCS (G4ParticleDefinition *aPartDef, const G4MaterialCutsCouple *aCouple, G4double &e_sigma_max, G4double &sigma_max)
 
G4double GetCrossSectionCorrection (G4ParticleDefinition *aPartDef, G4double PreStepEkin, const G4MaterialCutsCouple *aCouple, G4bool &fwd_is_used, G4double &fwd_TotCS)
 
void SetFwdCrossSectionMode (G4bool aBool)
 
G4double GetContinuousWeightCorrection (G4ParticleDefinition *aPartDef, G4double PreStepEkin, G4double AfterStepEkin, const G4MaterialCutsCouple *aCouple, G4double step_length)
 
G4double GetPostStepWeightCorrection ()
 
G4double ComputeAdjointCS (G4Material *aMaterial, G4VEmAdjointModel *aModel, G4double PrimEnergy, G4double Tcut, G4bool IsScatProjToProjCase, std::vector< G4double > &AdjointCS_for_each_element)
 
G4ElementSampleElementFromCSMatrices (G4Material *aMaterial, G4VEmAdjointModel *aModel, G4double PrimEnergy, G4double Tcut, G4bool IsScatProjToProjCase)
 
G4double ComputeTotalAdjointCS (const G4MaterialCutsCouple *aMatCutCouple, G4ParticleDefinition *aPart, G4double PrimEnergy)
 
G4ParticleDefinitionGetAdjointParticleEquivalent (G4ParticleDefinition *theFwdPartDef)
 
G4ParticleDefinitionGetForwardParticleEquivalent (G4ParticleDefinition *theAdjPartDef)
 
void SetTmin (G4double aVal)
 
void SetTmax (G4double aVal)
 
void SetNbins (G4int aInt)
 
void SetIon (G4ParticleDefinition *adjIon, G4ParticleDefinition *fwdIon)
 

Static Public Member Functions

static G4AdjointCSManagerGetAdjointCSManager ()
 

Detailed Description

Definition at line 69 of file G4AdjointCSManager.hh.

Constructor & Destructor Documentation

G4AdjointCSManager::~G4AdjointCSManager ( )

Definition at line 116 of file G4AdjointCSManager.cc.

117 {;
118 }

Member Function Documentation

void G4AdjointCSManager::BuildCrossSectionMatrices ( )

Definition at line 180 of file G4AdjointCSManager.cc.

References G4cout, G4endl, G4Element::GetA(), G4Element::GetElementTable(), G4Material::GetMaterialTable(), G4VEmAdjointModel::GetName(), G4VEmAdjointModel::GetUseMatrix(), G4VEmAdjointModel::GetUseMatrixPerElement(), G4VEmAdjointModel::GetUseOnlyOneMatrixForAllElements(), G4Element::GetZ(), int(), and G4VEmAdjointModel::SetCSMatrices().

Referenced by G4VAdjointReverseReaction::BuildPhysicsTable().

181 {
182  if (CrossSectionMatrixesAreBuilt) return;
183  //Tcut, Tmax
184  //The matrices will be computed probably just once
185  //When Tcut will change some PhysicsTable will be recomputed
186  // for each MaterialCutCouple but not all the matrices
187  //The Tcut defines a lower limit in the energy of the Projectile before the scattering
188  //In the Projectile to Scattered Projectile case we have
189  // E_ScatProj<E_Proj-Tcut
190  //Therefore in the adjoint case we have
191  // Eproj> E_ScatProj+Tcut
192  //This implies that when computing the adjoint CS we should integrate over Epro
193  // from E_ScatProj+Tcut to Emax
194  //In the Projectile to Secondary case Tcut plays a role only in the fact that
195  // Esecond should be greater than Tcut to have the possibility to have any adjoint
196  //process
197  //To avoid to recompute the matrices for all changes of MaterialCutCouple
198  //We propose to compute the matrices only once for the minimum possible Tcut and then
199  //to interpolate the probability for a new Tcut (implemented in G4VAdjointEmModel)
200 
201 
202  theAdjointCSMatricesForScatProjToProj.clear();
203  theAdjointCSMatricesForProdToProj.clear();
204  const G4ElementTable* theElementTable = G4Element::GetElementTable();
205  const G4MaterialTable* theMaterialTable = G4Material::GetMaterialTable();
206 
207  G4cout<<"========== Computation of cross section matrices for adjoint models =========="<<G4endl;
208  for (size_t i=0; i<listOfAdjointEMModel.size();i++){
209  G4VEmAdjointModel* aModel =listOfAdjointEMModel[i];
210  G4cout<<"Build adjoint cross section matrices for "<<aModel->GetName()<<G4endl;
211  if (aModel->GetUseMatrix()){
212  std::vector<G4AdjointCSMatrix*>* aListOfMat1 = new std::vector<G4AdjointCSMatrix*>();
213  std::vector<G4AdjointCSMatrix*>* aListOfMat2 = new std::vector<G4AdjointCSMatrix*>();
214  aListOfMat1->clear();
215  aListOfMat2->clear();
216  if (aModel->GetUseMatrixPerElement()){
217  if (aModel->GetUseOnlyOneMatrixForAllElements()){
218  std::vector<G4AdjointCSMatrix*>
219  two_matrices=BuildCrossSectionsMatricesForAGivenModelAndElement(aModel,1, 1, 80);
220  aListOfMat1->push_back(two_matrices[0]);
221  aListOfMat2->push_back(two_matrices[1]);
222  }
223  else {
224  for (size_t j=0; j<theElementTable->size();j++){
225  G4Element* anElement=(*theElementTable)[j];
226  G4int Z = int(anElement->GetZ());
227  G4int A = int(anElement->GetA());
228  std::vector<G4AdjointCSMatrix*>
229  two_matrices=BuildCrossSectionsMatricesForAGivenModelAndElement(aModel,Z, A, 40);
230  aListOfMat1->push_back(two_matrices[0]);
231  aListOfMat2->push_back(two_matrices[1]);
232  }
233  }
234  }
235  else { //Per material case
236  for (size_t j=0; j<theMaterialTable->size();j++){
237  G4Material* aMaterial=(*theMaterialTable)[j];
238  std::vector<G4AdjointCSMatrix*>
239  two_matrices=BuildCrossSectionsMatricesForAGivenModelAndMaterial(aModel,aMaterial, 40);
240  aListOfMat1->push_back(two_matrices[0]);
241  aListOfMat2->push_back(two_matrices[1]);
242  }
243 
244  }
245  theAdjointCSMatricesForProdToProj.push_back(*aListOfMat1);
246  theAdjointCSMatricesForScatProjToProj.push_back(*aListOfMat2);
247  aModel->SetCSMatrices(aListOfMat1, aListOfMat2);
248  }
249  else { G4cout<<"The model "<<aModel->GetName()<<" does not use cross section matrices"<<G4endl;
250  std::vector<G4AdjointCSMatrix*> two_empty_matrices;
251  theAdjointCSMatricesForProdToProj.push_back(two_empty_matrices);
252  theAdjointCSMatricesForScatProjToProj.push_back(two_empty_matrices);
253 
254  }
255  }
256  G4cout<<" All adjoint cross section matrices are computed!"<<G4endl;
257  G4cout<<"======================================================================"<<G4endl;
258 
259  CrossSectionMatrixesAreBuilt = true;
260 
261 
262 }
typedef int(XMLCALL *XML_NotStandaloneHandler)(void *userData)
G4double GetZ() const
Definition: G4Element.hh:131
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:564
std::vector< G4Material * > G4MaterialTable
G4bool GetUseOnlyOneMatrixForAllElements()
G4double GetA() const
Definition: G4Element.hh:138
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
void SetCSMatrices(std::vector< G4AdjointCSMatrix * > *Vec1CSMatrix, std::vector< G4AdjointCSMatrix * > *Vec2CSMatrix)
#define G4endl
Definition: G4ios.hh:61
std::vector< G4Element * > G4ElementTable
static G4ElementTable * GetElementTable()
Definition: G4Element.cc:395
void G4AdjointCSManager::BuildTotalSigmaTables ( )

Definition at line 267 of file G4AdjointCSManager.cc.

References ComputeTotalAdjointCS(), G4VEmModel::GetChargeSquareRatio(), G4MaterialCutsCouple::GetIndex(), G4PhysicsVector::GetLowEdgeEnergy(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4PhysicsVector::GetVectorLength(), and G4PhysicsVector::PutValue().

Referenced by G4VAdjointReverseReaction::BuildPhysicsTable().

268 { if (TotalSigmaTableAreBuilt) return;
269 
270 
272 
273 
274  //Prepare the Sigma table for all AdjointEMModel, will be filled later on
275  for (size_t i=0; i<listOfAdjointEMModel.size();i++){
276  listSigmaTableForAdjointModelScatProjToProj[i]->clearAndDestroy();
277  listSigmaTableForAdjointModelProdToProj[i]->clearAndDestroy();
278  for (size_t j=0;j<theCoupleTable->GetTableSize();j++){
279  listSigmaTableForAdjointModelScatProjToProj[i]->push_back(new G4PhysicsLogVector(Tmin, Tmax, nbins));
280  listSigmaTableForAdjointModelProdToProj[i]->push_back(new G4PhysicsLogVector(Tmin, Tmax, nbins));
281  }
282  }
283 
284 
285 
286  for (size_t i=0;i<theListOfAdjointParticlesInAction.size();i++){
287  G4ParticleDefinition* thePartDef = theListOfAdjointParticlesInAction[i];
288  DefineCurrentParticle(thePartDef);
289  theTotalForwardSigmaTableVector[i]->clearAndDestroy();
290  theTotalAdjointSigmaTableVector[i]->clearAndDestroy();
291  EminForFwdSigmaTables[i].clear();
292  EminForAdjSigmaTables[i].clear();
293  EkinofFwdSigmaMax[i].clear();
294  EkinofAdjSigmaMax[i].clear();
295  //G4cout<<thePartDef->GetParticleName();
296 
297  for (size_t j=0;j<theCoupleTable->GetTableSize();j++){
298  const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(j);
299 
300  /*
301  G4String file_name1=couple->GetMaterial()->GetName()+"_"+thePartDef->GetParticleName()+"_adj_totCS.txt";
302  G4String file_name2=couple->GetMaterial()->GetName()+"_"+thePartDef->GetParticleName()+"_fwd_totCS.txt";
303 
304  std::fstream FileOutputAdjCS(file_name1, std::ios::out);
305  std::fstream FileOutputFwdCS(file_name2, std::ios::out);
306 
307 
308 
309  FileOutputAdjCS<<std::setiosflags(std::ios::scientific);
310  FileOutputAdjCS<<std::setprecision(6);
311  FileOutputFwdCS<<std::setiosflags(std::ios::scientific);
312  FileOutputFwdCS<<std::setprecision(6);
313  */
314 
315 
316  //make first the total fwd CS table for FwdProcess
317  G4PhysicsVector* aVector = new G4PhysicsLogVector(Tmin, Tmax, nbins);
318  G4bool Emin_found=false;
319  G4double sigma_max =0.;
320  G4double e_sigma_max =0.;
321  for(size_t l=0; l<aVector->GetVectorLength(); l++) {
322  G4double totCS=0.;
323  G4double e=aVector->GetLowEdgeEnergy(l);
324  for (size_t k=0; k<listOfForwardEmProcess[i]->size(); k++){
325  totCS+=(*listOfForwardEmProcess[i])[k]->GetLambda(e, couple);
326  }
327  for (size_t k=0; k<listOfForwardEnergyLossProcess[i]->size(); k++){
328  if (thePartDef == theAdjIon) { // e is considered already as the scaled energy
329  size_t mat_index = couple->GetIndex();
330  G4VEmModel* currentModel = (*listOfForwardEnergyLossProcess[i])[k]->SelectModelForMaterial(e,mat_index);
331  G4double chargeSqRatio = currentModel->GetChargeSquareRatio(theFwdIon,couple->GetMaterial(),e/massRatio);
332  (*listOfForwardEnergyLossProcess[i])[k]->SetDynamicMassCharge(massRatio,chargeSqRatio);
333  }
334  G4double e1=e/massRatio;
335  totCS+=(*listOfForwardEnergyLossProcess[i])[k]->GetLambda(e1, couple);
336  }
337  aVector->PutValue(l,totCS);
338  if (totCS>sigma_max){
339  sigma_max=totCS;
340  e_sigma_max = e;
341 
342  }
343  //FileOutputFwdCS<<e<<'\t'<<totCS<<G4endl;
344 
345  if (totCS>0 && !Emin_found) {
346  EminForFwdSigmaTables[i].push_back(e);
347  Emin_found=true;
348  }
349 
350 
351  }
352  //FileOutputFwdCS.close();
353 
354  EkinofFwdSigmaMax[i].push_back(e_sigma_max);
355 
356 
357  if(!Emin_found) EminForFwdSigmaTables[i].push_back(Tmax);
358 
359  theTotalForwardSigmaTableVector[i]->push_back(aVector);
360 
361 
362  Emin_found=false;
363  sigma_max=0;
364  e_sigma_max =0.;
365  G4PhysicsVector* aVector1 = new G4PhysicsLogVector(Tmin, Tmax, nbins);
366  for(eindex=0; eindex<aVector->GetVectorLength(); eindex++) {
367  G4double e=aVector->GetLowEdgeEnergy(eindex);
368  G4double totCS =ComputeTotalAdjointCS(couple,thePartDef,e*0.9999999/massRatio); //massRatio needed for ions
369  aVector1->PutValue(eindex,totCS);
370  if (totCS>sigma_max){
371  sigma_max=totCS;
372  e_sigma_max = e;
373 
374  }
375  //FileOutputAdjCS<<e<<'\t'<<totCS<<G4endl;
376  if (totCS>0 && !Emin_found) {
377  EminForAdjSigmaTables[i].push_back(e);
378  Emin_found=true;
379  }
380 
381  }
382  //FileOutputAdjCS.close();
383  EkinofAdjSigmaMax[i].push_back(e_sigma_max);
384  if(!Emin_found) EminForAdjSigmaTables[i].push_back(Tmax);
385 
386  theTotalAdjointSigmaTableVector[i]->push_back(aVector1);
387 
388  }
389  }
390  TotalSigmaTableAreBuilt =true;
391 
392 }
G4double ComputeTotalAdjointCS(const G4MaterialCutsCouple *aMatCutCouple, G4ParticleDefinition *aPart, G4double PrimEnergy)
size_t GetVectorLength() const
G4double GetLowEdgeEnergy(size_t binNumber) const
bool G4bool
Definition: G4Types.hh:79
void PutValue(size_t index, G4double theValue)
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:322
static G4ProductionCutsTable * GetProductionCutsTable()
const G4MaterialCutsCouple * GetMaterialCutsCouple(G4int i) const
double G4double
Definition: G4Types.hh:76
const G4Material * GetMaterial() const
G4double G4AdjointCSManager::ComputeAdjointCS ( G4Material aMaterial,
G4VEmAdjointModel aModel,
G4double  PrimEnergy,
G4double  Tcut,
G4bool  IsScatProjToProjCase,
std::vector< G4double > &  AdjointCS_for_each_element 
)

Definition at line 532 of file G4AdjointCSManager.cc.

References G4VEmAdjointModel::AdjointCrossSection(), G4VEmAdjointModel::GetApplyCutInRange(), G4Material::GetElement(), G4Element::GetIndex(), G4Material::GetIndex(), G4Material::GetNumberOfElements(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMaxForScatProjToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForProdToProjCase(), G4VEmAdjointModel::GetSecondAdjEnergyMinForScatProjToProjCase(), G4VEmAdjointModel::GetUseMatrix(), G4VEmAdjointModel::GetUseMatrixPerElement(), G4VEmAdjointModel::GetUseOnlyOneMatrixForAllElements(), G4Material::GetVecNbOfAtomsPerVolume(), and G4Element::GetZ().

Referenced by G4VEmAdjointModel::AdjointCrossSection(), ComputeTotalAdjointCS(), and SampleElementFromCSMatrices().

538 {
539 
540  G4double EminSec=0;
541  G4double EmaxSec=0;
542 
543  if (IsScatProjToProjCase){
544  EminSec= aModel->GetSecondAdjEnergyMinForScatProjToProjCase(PrimEnergy,Tcut);
545  EmaxSec= aModel->GetSecondAdjEnergyMaxForScatProjToProjCase(PrimEnergy);
546  }
547  else if (PrimEnergy > Tcut || !aModel->GetApplyCutInRange()) {
548  EminSec= aModel->GetSecondAdjEnergyMinForProdToProjCase(PrimEnergy);
549  EmaxSec= aModel->GetSecondAdjEnergyMaxForProdToProjCase(PrimEnergy);
550  }
551  if (EminSec >= EmaxSec) return 0.;
552 
553 
554  G4bool need_to_compute=false;
555  if ( aMaterial!= lastMaterial || PrimEnergy != lastPrimaryEnergy || Tcut != lastTcut){
556  lastMaterial =aMaterial;
557  lastPrimaryEnergy = PrimEnergy;
558  lastTcut=Tcut;
559  listOfIndexOfAdjointEMModelInAction.clear();
560  listOfIsScatProjToProjCase.clear();
561  lastAdjointCSVsModelsAndElements.clear();
562  need_to_compute=true;
563 
564  }
565  size_t ind=0;
566  if (!need_to_compute){
567  need_to_compute=true;
568  for (size_t i=0;i<listOfIndexOfAdjointEMModelInAction.size();i++){
569  size_t ind1=listOfIndexOfAdjointEMModelInAction[i];
570  if (aModel == listOfAdjointEMModel[ind1] && IsScatProjToProjCase == listOfIsScatProjToProjCase[i]){
571  need_to_compute=false;
572  CS_Vs_Element = lastAdjointCSVsModelsAndElements[ind];
573  }
574  ind++;
575  }
576  }
577 
578  if (need_to_compute){
579  size_t ind_model=0;
580  for (size_t i=0;i<listOfAdjointEMModel.size();i++){
581  if (aModel == listOfAdjointEMModel[i]){
582  ind_model=i;
583  break;
584  }
585  }
586  G4double Tlow=Tcut;
587  if (!listOfAdjointEMModel[ind_model]->GetApplyCutInRange()) Tlow =listOfAdjointEMModel[ind_model]->GetLowEnergyLimit();
588  listOfIndexOfAdjointEMModelInAction.push_back(ind_model);
589  listOfIsScatProjToProjCase.push_back(IsScatProjToProjCase);
590  CS_Vs_Element.clear();
591  if (!aModel->GetUseMatrix()){
592  CS_Vs_Element.push_back(aModel->AdjointCrossSection(currentCouple,PrimEnergy,IsScatProjToProjCase));
593 
594 
595  }
596  else if (aModel->GetUseMatrixPerElement()){
597  size_t n_el = aMaterial->GetNumberOfElements();
598  if (aModel->GetUseOnlyOneMatrixForAllElements()){
599  G4AdjointCSMatrix* theCSMatrix;
600  if (IsScatProjToProjCase){
601  theCSMatrix=theAdjointCSMatricesForScatProjToProj[ind_model][0];
602  }
603  else theCSMatrix=theAdjointCSMatricesForProdToProj[ind_model][0];
604  G4double CS =0.;
605  if (PrimEnergy > Tlow)
606  CS = ComputeAdjointCS(PrimEnergy,theCSMatrix,Tlow);
607  G4double factor=0.;
608  for (size_t i=0;i<n_el;i++){ //this could be computed only once
609  //size_t ind_el = aMaterial->GetElement(i)->GetIndex();
610  factor+=aMaterial->GetElement(i)->GetZ()*aMaterial->GetVecNbOfAtomsPerVolume()[i];
611  }
612  CS *=factor;
613  CS_Vs_Element.push_back(CS);
614 
615  }
616  else {
617  for (size_t i=0;i<n_el;i++){
618  size_t ind_el = aMaterial->GetElement(i)->GetIndex();
619  //G4cout<<aMaterial->GetName()<<G4endl;
620  G4AdjointCSMatrix* theCSMatrix;
621  if (IsScatProjToProjCase){
622  theCSMatrix=theAdjointCSMatricesForScatProjToProj[ind_model][ind_el];
623  }
624  else theCSMatrix=theAdjointCSMatricesForProdToProj[ind_model][ind_el];
625  G4double CS =0.;
626  if (PrimEnergy > Tlow)
627  CS = ComputeAdjointCS(PrimEnergy,theCSMatrix,Tlow);
628  //G4cout<<CS<<G4endl;
629  CS_Vs_Element.push_back(CS*(aMaterial->GetVecNbOfAtomsPerVolume()[i]));
630  }
631  }
632 
633  }
634  else {
635  size_t ind_mat = aMaterial->GetIndex();
636  G4AdjointCSMatrix* theCSMatrix;
637  if (IsScatProjToProjCase){
638  theCSMatrix=theAdjointCSMatricesForScatProjToProj[ind_model][ind_mat];
639  }
640  else theCSMatrix=theAdjointCSMatricesForProdToProj[ind_model][ind_mat];
641  G4double CS =0.;
642  if (PrimEnergy > Tlow)
643  CS = ComputeAdjointCS(PrimEnergy,theCSMatrix,Tlow);
644  CS_Vs_Element.push_back(CS);
645 
646 
647  }
648  lastAdjointCSVsModelsAndElements.push_back(CS_Vs_Element);
649 
650  }
651 
652 
653  G4double CS=0;
654  for (size_t i=0;i<CS_Vs_Element.size();i++){
655  CS+=CS_Vs_Element[i]; //We could put the progressive sum of the CS instead of the CS of an element itself
656 
657  }
658  return CS;
659 }
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
size_t GetIndex() const
Definition: G4Material.hh:260
G4double GetZ() const
Definition: G4Element.hh:131
G4bool GetUseOnlyOneMatrixForAllElements()
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:200
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase(G4double PrimAdjEnergy)
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:204
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
bool G4bool
Definition: G4Types.hh:79
size_t GetIndex() const
Definition: G4Element.hh:181
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase(G4double PrimAdjEnergy, G4double Tcut=0)
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
double G4double
Definition: G4Types.hh:76
G4double ComputeAdjointCS(G4Material *aMaterial, G4VEmAdjointModel *aModel, G4double PrimEnergy, G4double Tcut, G4bool IsScatProjToProjCase, std::vector< G4double > &AdjointCS_for_each_element)
G4double G4AdjointCSManager::ComputeTotalAdjointCS ( const G4MaterialCutsCouple aMatCutCouple,
G4ParticleDefinition aPart,
G4double  PrimEnergy 
)

Definition at line 687 of file G4AdjointCSManager.cc.

References ComputeAdjointCS(), G4ProductionCutsTable::GetEnergyCutsVector(), GetForwardParticleEquivalent(), G4MaterialCutsCouple::GetIndex(), G4ParticleDefinition::GetParticleName(), and G4ProductionCutsTable::GetProductionCutsTable().

Referenced by BuildTotalSigmaTables().

690 {
691  G4double TotalCS=0.;
692 
693  DefineCurrentMaterial(aCouple);
694 
695 
696  std::vector<G4double> CS_Vs_Element;
697  G4double CS;
698  for (size_t i=0; i<listOfAdjointEMModel.size();i++){
699 
700  G4double Tlow=0;
701  if (!listOfAdjointEMModel[i]->GetApplyCutInRange()) Tlow =listOfAdjointEMModel[i]->GetLowEnergyLimit();
702  else {
703  G4ParticleDefinition* theDirSecondPartDef =
704  GetForwardParticleEquivalent(listOfAdjointEMModel[i]->GetAdjointEquivalentOfDirectSecondaryParticleDefinition());
705  size_t idx=56;
706  if (theDirSecondPartDef->GetParticleName() == "gamma") idx = 0;
707  else if (theDirSecondPartDef->GetParticleName() == "e-") idx = 1;
708  else if (theDirSecondPartDef->GetParticleName() == "e+") idx = 2;
709  if (idx <56) {
710  const std::vector<G4double>* aVec = G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(idx);
711  Tlow =(*aVec)[aCouple->GetIndex()];
712  }
713 
714 
715  }
716  if ( Ekin<=listOfAdjointEMModel[i]->GetHighEnergyLimit() && Ekin>=listOfAdjointEMModel[i]->GetLowEnergyLimit()){
717  if (aPartDef == listOfAdjointEMModel[i]->GetAdjointEquivalentOfDirectPrimaryParticleDefinition()){
718  CS=ComputeAdjointCS(currentMaterial,
719  listOfAdjointEMModel[i],
720  Ekin, Tlow,true,CS_Vs_Element);
721  TotalCS += CS;
722  (*listSigmaTableForAdjointModelScatProjToProj[i])[currentMatIndex]->PutValue(eindex,CS);
723  }
724  if (aPartDef == listOfAdjointEMModel[i]->GetAdjointEquivalentOfDirectSecondaryParticleDefinition()){
725  CS = ComputeAdjointCS(currentMaterial,
726  listOfAdjointEMModel[i],
727  Ekin, Tlow,false, CS_Vs_Element);
728  TotalCS += CS;
729  (*listSigmaTableForAdjointModelProdToProj[i])[currentMatIndex]->PutValue(eindex,CS);
730  }
731 
732  }
733  else {
734  (*listSigmaTableForAdjointModelScatProjToProj[i])[currentMatIndex]->PutValue(eindex,0.);
735  (*listSigmaTableForAdjointModelProdToProj[i])[currentMatIndex]->PutValue(eindex,0.);
736 
737  }
738  }
739  return TotalCS;
740 
741 
742 }
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
const G4String & GetParticleName() const
G4ParticleDefinition * GetForwardParticleEquivalent(G4ParticleDefinition *theAdjPartDef)
static G4ProductionCutsTable * GetProductionCutsTable()
double G4double
Definition: G4Types.hh:76
G4double ComputeAdjointCS(G4Material *aMaterial, G4VEmAdjointModel *aModel, G4double PrimEnergy, G4double Tcut, G4bool IsScatProjToProjCase, std::vector< G4double > &AdjointCS_for_each_element)
G4AdjointCSManager * G4AdjointCSManager::GetAdjointCSManager ( )
static
G4ParticleDefinition * G4AdjointCSManager::GetAdjointParticleEquivalent ( G4ParticleDefinition theFwdPartDef)

Definition at line 943 of file G4AdjointCSManager.cc.

References G4AdjointElectron::AdjointElectron(), G4AdjointGamma::AdjointGamma(), G4AdjointProton::AdjointProton(), and G4ParticleDefinition::GetParticleName().

Referenced by RegisterEmProcess(), and RegisterEnergyLossProcess().

944 {
945  if (theFwdPartDef->GetParticleName() == "e-") return G4AdjointElectron::AdjointElectron();
946  else if (theFwdPartDef->GetParticleName() == "gamma") return G4AdjointGamma::AdjointGamma();
947  else if (theFwdPartDef->GetParticleName() == "proton") return G4AdjointProton::AdjointProton();
948  else if (theFwdPartDef ==theFwdIon) return theAdjIon;
949 
950  return 0;
951 }
static G4AdjointGamma * AdjointGamma()
static G4AdjointElectron * AdjointElectron()
const G4String & GetParticleName() const
static G4AdjointProton * AdjointProton()
G4double G4AdjointCSManager::GetAdjointSigma ( G4double  Ekin_nuc,
size_t  index_model,
G4bool  is_scat_proj_to_proj,
const G4MaterialCutsCouple aCouple 
)

Definition at line 418 of file G4AdjointCSManager.cc.

References test::b.

420 { DefineCurrentMaterial(aCouple);
421  G4bool b;
422  if (is_scat_proj_to_proj) return (((*listSigmaTableForAdjointModelScatProjToProj[index_model])[currentMatIndex])->GetValue(Ekin_nuc, b));
423  else return (((*listSigmaTableForAdjointModelProdToProj[index_model])[currentMatIndex])->GetValue(Ekin_nuc, b));
424 }
bool G4bool
Definition: G4Types.hh:79
G4double G4AdjointCSManager::GetContinuousWeightCorrection ( G4ParticleDefinition aPartDef,
G4double  PreStepEkin,
G4double  AfterStepEkin,
const G4MaterialCutsCouple aCouple,
G4double  step_length 
)

Definition at line 502 of file G4AdjointCSManager.cc.

References GetTotalAdjointCS(), and GetTotalForwardCS().

Referenced by G4AdjointAlongStepWeightCorrection::AlongStepDoIt().

504 { G4double corr_fac = 1.;
505  //return corr_fac;
506  //G4double after_adjCS = GetTotalAdjointCS(aPartDef, AfterStepEkin,aCouple);
507  G4double after_fwdCS = GetTotalForwardCS(aPartDef, AfterStepEkin,aCouple);
508  G4double pre_adjCS = GetTotalAdjointCS(aPartDef, PreStepEkin,aCouple);
509  if (!forward_CS_is_used || pre_adjCS ==0. || after_fwdCS==0.) {
510  forward_CS_is_used=false;
511  G4double pre_fwdCS = GetTotalForwardCS(aPartDef, PreStepEkin,aCouple);
512  corr_fac *=std::exp((pre_adjCS-pre_fwdCS)*step_length);
513  LastCSCorrectionFactor = 1.;
514  }
515  else {
516  LastCSCorrectionFactor = after_fwdCS/pre_adjCS;
517  }
518 
519 
520 
521  return corr_fac;
522 }
G4double GetTotalAdjointCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
G4double GetTotalForwardCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
double G4double
Definition: G4Types.hh:76
G4double G4AdjointCSManager::GetCrossSectionCorrection ( G4ParticleDefinition aPartDef,
G4double  PreStepEkin,
const G4MaterialCutsCouple aCouple,
G4bool fwd_is_used,
G4double fwd_TotCS 
)

Definition at line 465 of file G4AdjointCSManager.cc.

References GetTotalAdjointCS(), and GetTotalForwardCS().

Referenced by G4VAdjointReverseReaction::GetMeanFreePath().

467 { G4double corr_fac = 1.;
468  if (forward_CS_mode && aPartDef ) {
469  fwd_TotCS=PrefwdCS;
470  if (LastEkinForCS != PreStepEkin || aPartDef != lastPartDefForCS || aCouple!=currentCouple) {
471  DefineCurrentMaterial(aCouple);
472  PreadjCS = GetTotalAdjointCS(aPartDef, PreStepEkin,aCouple);
473  PrefwdCS = GetTotalForwardCS(aPartDef, PreStepEkin,aCouple);
474  LastEkinForCS = PreStepEkin;
475  lastPartDefForCS = aPartDef;
476  if (PrefwdCS >0. && PreadjCS >0.) {
477  forward_CS_is_used = true;
478  LastCSCorrectionFactor = PrefwdCS/PreadjCS;
479  }
480  else {
481  forward_CS_is_used = false;
482  LastCSCorrectionFactor = 1.;
483 
484  }
485 
486  }
487  corr_fac =LastCSCorrectionFactor;
488 
489 
490 
491  }
492  else {
493  forward_CS_is_used = false;
494  LastCSCorrectionFactor = 1.;
495  }
496  fwd_TotCS=PrefwdCS;
497  fwd_is_used = forward_CS_is_used;
498  return corr_fac;
499 }
G4double GetTotalAdjointCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
G4double GetTotalForwardCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
double G4double
Definition: G4Types.hh:76
void G4AdjointCSManager::GetEminForTotalCS ( G4ParticleDefinition aPartDef,
const G4MaterialCutsCouple aCouple,
G4double emin_adj,
G4double emin_fwd 
)

Definition at line 427 of file G4AdjointCSManager.cc.

429 { DefineCurrentMaterial(aCouple);
430  DefineCurrentParticle(aPartDef);
431  emin_adj = EminForAdjSigmaTables[currentParticleIndex][currentMatIndex]/massRatio;
432  emin_fwd = EminForFwdSigmaTables[currentParticleIndex][currentMatIndex]/massRatio;
433 
434 
435 
436 }
G4ParticleDefinition * G4AdjointCSManager::GetForwardParticleEquivalent ( G4ParticleDefinition theAdjPartDef)

Definition at line 954 of file G4AdjointCSManager.cc.

References G4Electron::Electron(), G4Gamma::Gamma(), G4ParticleDefinition::GetParticleName(), and G4Proton::Proton().

Referenced by ComputeTotalAdjointCS().

955 {
956  if (theAdjPartDef->GetParticleName() == "adj_e-") return G4Electron::Electron();
957  else if (theAdjPartDef->GetParticleName() == "adj_gamma") return G4Gamma::Gamma();
958  else if (theAdjPartDef->GetParticleName() == "adj_proton") return G4Proton::Proton();
959  else if (theAdjPartDef == theAdjIon) return theFwdIon;
960  return 0;
961 }
const G4String & GetParticleName() const
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
static G4Electron * Electron()
Definition: G4Electron.cc:94
void G4AdjointCSManager::GetMaxAdjTotalCS ( G4ParticleDefinition aPartDef,
const G4MaterialCutsCouple aCouple,
G4double e_sigma_max,
G4double sigma_max 
)

Definition at line 452 of file G4AdjointCSManager.cc.

References test::b.

454 { DefineCurrentMaterial(aCouple);
455  DefineCurrentParticle(aPartDef);
456  e_sigma_max = EkinofAdjSigmaMax[currentParticleIndex][currentMatIndex];
457  G4bool b;
458  sigma_max =((*theTotalAdjointSigmaTableVector[currentParticleIndex])[currentMatIndex])->GetValue(e_sigma_max, b);
459  e_sigma_max/=massRatio;
460 
461 
462 }
bool G4bool
Definition: G4Types.hh:79
void G4AdjointCSManager::GetMaxFwdTotalCS ( G4ParticleDefinition aPartDef,
const G4MaterialCutsCouple aCouple,
G4double e_sigma_max,
G4double sigma_max 
)

Definition at line 439 of file G4AdjointCSManager.cc.

References test::b.

441 { DefineCurrentMaterial(aCouple);
442  DefineCurrentParticle(aPartDef);
443  e_sigma_max = EkinofFwdSigmaMax[currentParticleIndex][currentMatIndex];
444  G4bool b;
445  sigma_max =((*theTotalForwardSigmaTableVector[currentParticleIndex])[currentMatIndex])->GetValue(e_sigma_max, b);
446  e_sigma_max/=massRatio;
447 
448 
449 }
bool G4bool
Definition: G4Types.hh:79
G4int G4AdjointCSManager::GetNbProcesses ( )
G4double G4AdjointCSManager::GetPostStepWeightCorrection ( )
G4double G4AdjointCSManager::GetTotalAdjointCS ( G4ParticleDefinition aPartDef,
G4double  Ekin,
const G4MaterialCutsCouple aCouple 
)

Definition at line 395 of file G4AdjointCSManager.cc.

References test::b.

Referenced by GetContinuousWeightCorrection(), and GetCrossSectionCorrection().

397 { DefineCurrentMaterial(aCouple);
398  DefineCurrentParticle(aPartDef);
399  G4bool b;
400  return (((*theTotalAdjointSigmaTableVector[currentParticleIndex])[currentMatIndex])->GetValue(Ekin*massRatio, b));
401 
402 
403 
404 }
bool G4bool
Definition: G4Types.hh:79
G4double G4AdjointCSManager::GetTotalForwardCS ( G4ParticleDefinition aPartDef,
G4double  Ekin,
const G4MaterialCutsCouple aCouple 
)

Definition at line 407 of file G4AdjointCSManager.cc.

References test::b.

Referenced by GetContinuousWeightCorrection(), and GetCrossSectionCorrection().

409 { DefineCurrentMaterial(aCouple);
410  DefineCurrentParticle(aPartDef);
411  G4bool b;
412  return (((*theTotalForwardSigmaTableVector[currentParticleIndex])[currentMatIndex])->GetValue(Ekin*massRatio, b));
413 
414 
415 }
bool G4bool
Definition: G4Types.hh:79
void G4AdjointCSManager::RegisterAdjointParticle ( G4ParticleDefinition aPartDef)

Definition at line 159 of file G4AdjointCSManager.cc.

References G4ParticleDefinition::GetParticleName().

Referenced by G4AdjointPhysicsList::ConstructEM(), RegisterEmProcess(), and RegisterEnergyLossProcess().

160 { G4int index=-1;
161  for (size_t i=0;i<theListOfAdjointParticlesInAction.size();i++){
162  if (aPartDef->GetParticleName() == theListOfAdjointParticlesInAction[i]->GetParticleName()) index=i;
163  }
164 
165  if (index ==-1){
166  listOfForwardEnergyLossProcess.push_back(new std::vector<G4VEnergyLossProcess*>());
167  theTotalForwardSigmaTableVector.push_back(new G4PhysicsTable);
168  theTotalAdjointSigmaTableVector.push_back(new G4PhysicsTable);
169  listOfForwardEmProcess.push_back(new std::vector<G4VEmProcess*>());
170  theListOfAdjointParticlesInAction.push_back(aPartDef);
171  EminForFwdSigmaTables.push_back(std::vector<G4double> ());
172  EminForAdjSigmaTables.push_back(std::vector<G4double> ());
173  EkinofFwdSigmaMax.push_back(std::vector<G4double> ());
174  EkinofAdjSigmaMax.push_back(std::vector<G4double> ());
175 
176  }
177 }
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
size_t G4AdjointCSManager::RegisterEmAdjointModel ( G4VEmAdjointModel aModel)

Definition at line 121 of file G4AdjointCSManager.cc.

Referenced by G4VEmAdjointModel::G4VEmAdjointModel().

122 {listOfAdjointEMModel.push_back(aModel);
123  listSigmaTableForAdjointModelScatProjToProj.push_back(new G4PhysicsTable);
124  listSigmaTableForAdjointModelProdToProj.push_back(new G4PhysicsTable);
125  return listOfAdjointEMModel.size() -1;
126 
127 }
void G4AdjointCSManager::RegisterEmProcess ( G4VEmProcess aProcess,
G4ParticleDefinition aPartDef 
)

Definition at line 130 of file G4AdjointCSManager.cc.

References GetAdjointParticleEquivalent(), G4ParticleDefinition::GetParticleName(), and RegisterAdjointParticle().

131 {
132  G4ParticleDefinition* anAdjPartDef = GetAdjointParticleEquivalent(aFwdPartDef);
133  if (anAdjPartDef && aProcess){
134  RegisterAdjointParticle(anAdjPartDef);
135  G4int index=-1;
136 
137  for (size_t i=0;i<theListOfAdjointParticlesInAction.size();i++){
138  if (anAdjPartDef->GetParticleName() == theListOfAdjointParticlesInAction[i]->GetParticleName()) index=i;
139  }
140  listOfForwardEmProcess[index]->push_back(aProcess);
141  }
142 }
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void RegisterAdjointParticle(G4ParticleDefinition *aPartDef)
G4ParticleDefinition * GetAdjointParticleEquivalent(G4ParticleDefinition *theFwdPartDef)
void G4AdjointCSManager::RegisterEnergyLossProcess ( G4VEnergyLossProcess aProcess,
G4ParticleDefinition aPartDef 
)

Definition at line 145 of file G4AdjointCSManager.cc.

References GetAdjointParticleEquivalent(), G4ParticleDefinition::GetParticleName(), and RegisterAdjointParticle().

146 {
147  G4ParticleDefinition* anAdjPartDef = GetAdjointParticleEquivalent(aFwdPartDef);
148  if (anAdjPartDef && aProcess){
149  RegisterAdjointParticle(anAdjPartDef);
150  G4int index=-1;
151  for (size_t i=0;i<theListOfAdjointParticlesInAction.size();i++){
152  if (anAdjPartDef->GetParticleName() == theListOfAdjointParticlesInAction[i]->GetParticleName()) index=i;
153  }
154  listOfForwardEnergyLossProcess[index]->push_back(aProcess);
155  }
156 }
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void RegisterAdjointParticle(G4ParticleDefinition *aPartDef)
G4ParticleDefinition * GetAdjointParticleEquivalent(G4ParticleDefinition *theFwdPartDef)
G4Element * G4AdjointCSManager::SampleElementFromCSMatrices ( G4Material aMaterial,
G4VEmAdjointModel aModel,
G4double  PrimEnergy,
G4double  Tcut,
G4bool  IsScatProjToProjCase 
)

Definition at line 662 of file G4AdjointCSManager.cc.

References ComputeAdjointCS(), G4UniformRand, and G4Material::GetElement().

667 { std::vector<G4double> CS_Vs_Element;
668  G4double CS = ComputeAdjointCS(aMaterial,aModel,PrimEnergy,Tcut,IsScatProjToProjCase,CS_Vs_Element);
669  G4double rand_var= G4UniformRand();
670  G4double SumCS=0.;
671  size_t ind=0;
672  for (size_t i=0;i<CS_Vs_Element.size();i++){
673  SumCS+=CS_Vs_Element[i];
674  if (rand_var<=SumCS/CS){
675  ind=i;
676  break;
677  }
678  }
679 
680  return const_cast<G4Element*>(aMaterial->GetElement(ind));
681 
682 
683 
684 }
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:200
#define G4UniformRand()
Definition: Randomize.hh:87
double G4double
Definition: G4Types.hh:76
G4double ComputeAdjointCS(G4Material *aMaterial, G4VEmAdjointModel *aModel, G4double PrimEnergy, G4double Tcut, G4bool IsScatProjToProjCase, std::vector< G4double > &AdjointCS_for_each_element)
void G4AdjointCSManager::SetFwdCrossSectionMode ( G4bool  aBool)
inline

Definition at line 122 of file G4AdjointCSManager.hh.

122 {forward_CS_mode=aBool;}
void G4AdjointCSManager::SetIon ( G4ParticleDefinition adjIon,
G4ParticleDefinition fwdIon 
)
inline

Definition at line 169 of file G4AdjointCSManager.hh.

170  {theAdjIon=adjIon; theFwdIon =fwdIon;}
void G4AdjointCSManager::SetNbins ( G4int  aInt)
inline

Definition at line 168 of file G4AdjointCSManager.hh.

168 {nbins=aInt;}
void G4AdjointCSManager::SetTmax ( G4double  aVal)
inline

Definition at line 167 of file G4AdjointCSManager.hh.

167 {Tmax=aVal;}
void G4AdjointCSManager::SetTmin ( G4double  aVal)
inline

Definition at line 166 of file G4AdjointCSManager.hh.

166 {Tmin=aVal;}

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