Message: Re: Near thermal Neutron simulation problem Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Near thermal Neutron simulation problem 

Forum: Physics List
Re: None Near thermal Neutron simulation problem (Alan Shippen)
Date: 01 Mar, 2008
From: Vladimir IVANTCHENKO <vnivanch@mail.cern.ch>

On Thu, 28 Feb 2008, Alan wrote:

> *** Discussion title: Physics List
> Email replies to PublicHyperNews@slac.stanford.edu must include:
>   In-Reply-To: <"/phys-list/458"@geant4-hn.slac.stanford.edu>
>   Subject: ...change this to be about your reply.
> 
> I have been trying to simulate the creation of neutrons from an Al+H2O
> target hit by a 100MeV proton beam line. The goal of the simulations is
> to obtain a value for neutron flux at the end of a path which is offset
> from the target by about 20cm. The target is surrounded by a water bath
> 40cm in radius and 80cm in depth with the pipe measuring neutrons coming
> directly off this pipe. The measurement of neutron flux occurs exactly
> 1m down this beam line. Other that the aforementioned structures the
> rest of the setup is either shielding concrete or steal. During the
> simulation neutrons are produced as you should expect however when the
> simulation completes I get results for flux for energies above 2keV
> however any bin below 2keV is 0. This is in stark contrast to other
> simulation toolboxes I use (MCNPX and FLUKA) which populate the below
> 2keV bins. Iím wondering if Iím doing something wrong with the physics
> list, making the below 2keV neutrons not be simulated. With this in mind
> I have included the physics list I use. Otherwise I have played with the
> tracking cuts on many different simulations but to no success, though I
> could do doing something wrong there too. Does anybody have any
> suggestions?
> 
> TriumfPhysicsList::TriumfPhysicsList() : G4VUserPhysicsList()
> {
>   SetVerboseLevel(1);
> }
> 
> TriumfPhysicsList::~TriumfPhysicsList() {}
> 
> void TriumfPhysicsList::ConstructParticle()
> {
>   // In this method, static member functions should be called
>   // for all particles which you want to use.
>   // This ensures that objects of these particle types will be
>   // created in the program. 
> 
>   ConstructAllBosons();
>   ConstructAllLeptons();
>   ConstructAllMesons();
>   ConstructAllBaryons();
>   ConstructAllIons();
>   ConstructAllShortLiveds();
> }
> 
> void TriumfPhysicsList::ConstructAllBosons()
> {
>   // Construct all bosons
>   G4BosonConstructor pConstructor;
>   pConstructor.ConstructParticle();
> }
> 
> void TriumfPhysicsList::ConstructAllLeptons()
> {
>   // Construct all leptons
>   G4LeptonConstructor pConstructor;
>   pConstructor.ConstructParticle();
> }
> 
> void TriumfPhysicsList::ConstructAllMesons()
> {
>   //  Construct all mesons
>   G4MesonConstructor pConstructor;
>   pConstructor.ConstructParticle();
> }
> 
> void TriumfPhysicsList::ConstructAllBaryons()
> {
>   //  Construct all barions
>   G4BaryonConstructor pConstructor;
>   pConstructor.ConstructParticle();
> }
> 
> void TriumfPhysicsList::ConstructAllIons()
> {
>   //  Construct light ions
>   G4IonConstructor pConstructor;
>   pConstructor.ConstructParticle();  
> }
> 
> void TriumfPhysicsList::ConstructAllShortLiveds()
> {
>   //  Construct  resonaces and quarks
>   G4ShortLivedConstructor pConstructor;
>   pConstructor.ConstructParticle();  
> }
> 
> void TriumfPhysicsList::ConstructProcess()
> {
>   AddTransportation();
>   ConstructEM();
>   ConstructLeptHad();
>   ConstructHad();
>   ConstructGeneral();
> }
> 
> #include "G4ComptonScattering.hh"
> #include "G4GammaConversion.hh"
> #include "G4PhotoElectricEffect.hh"
> 
> #include "G4MultipleScattering.hh"
> 
> #include "G4eIonisation.hh"
> #include "G4eBremsstrahlung.hh"
> #include "G4eplusAnnihilation.hh"
> 
> #include "G4MuIonisation.hh"
> #include "G4MuBremsstrahlung.hh"
> #include "G4MuPairProduction.hh"
> 
> #include "G4hIonisation.hh"
> #include "G4ionIonisation.hh"
> 
> void TriumfPhysicsList::ConstructEM()
> {
>   theParticleIterator->reset();
>   while( (*theParticleIterator)() ){
>     G4ParticleDefinition* particle = theParticleIterator->value();
>     G4ProcessManager* pmanager = particle->GetProcessManager();
>     G4String particleName = particle->GetParticleName();
> 
>     if (particleName == "gamma") {
>     // gamma
>       // Construct processes for gamma
>       pmanager->AddDiscreteProcess(new G4GammaConversion());
>       pmanager->AddDiscreteProcess(new G4ComptonScattering());      
>       pmanager->AddDiscreteProcess(new G4PhotoElectricEffect());
> 
>     } else if (particleName == "e-") {
>     //electron
>       // Construct processes for electron
>       pmanager->AddProcess(new G4MultipleScattering(),-1,1,1);
>       pmanager->AddProcess(new G4eIonisation(),-1,2,2);
>       pmanager->AddProcess(new G4eBremsstrahlung(),-1,-1,3);
> 
>     } else if (particleName == "e+") {
>     //positron
>       // Construct processes for positron
>      pmanager->AddProcess(new G4MultipleScattering(),-1,1,1);
> 
>      pmanager->AddProcess(new G4eIonisation(),-1,2,2);
>      pmanager->AddProcess(new G4eBremsstrahlung(),-1,-1,3);      
>      pmanager->AddProcess(new G4eplusAnnihilation(),0,-1,4);
> 
>     } else if( particleName == "mu+" || 
>                particleName == "mu-"    ) {
>     //muon  
>      // Construct processes for muon+
>      pmanager->AddProcess(new G4MultipleScattering(),-1,1,1);
>      pmanager->AddProcess(new G4MuIonisation(),-1,2,2);
>      pmanager->AddProcess(new G4MuBremsstrahlung(),-1,-1,3);
>      pmanager->AddProcess(new G4MuPairProduction(),-1,-1,4);       
> 
>     } else if( particleName == "GenericIon" ) {
>       pmanager->AddProcess(new G4MultipleScattering(),-1,1,1);
>       pmanager->AddProcess(new G4ionIonisation(),-1,2,2); 
>     } else { 
>       if ((particle->GetPDGCharge() != 0.0) && 
>           (particle->GetParticleName() != "chargedgeantino")&&
>           (!particle->IsShortLived()) ) {
>      // all others charged particles except geantino
>        pmanager->AddProcess(new G4MultipleScattering(),-1,1,1);
>        pmanager->AddProcess(new G4hIonisation(),-1,2,2);       
>      }
>     }
>   }
> }
> 
> // Hadron Processes
> 
> #include "G4HadronElasticProcess.hh"
> #include "G4HadronFissionProcess.hh"
> #include "G4HadronCaptureProcess.hh"
> 
> #include "G4PionPlusInelasticProcess.hh"
> #include "G4PionMinusInelasticProcess.hh"
> #include "G4KaonPlusInelasticProcess.hh"
> #include "G4KaonZeroSInelasticProcess.hh"
> #include "G4KaonZeroLInelasticProcess.hh"
> #include "G4KaonMinusInelasticProcess.hh"
> #include "G4ProtonInelasticProcess.hh"
> #include "G4AntiProtonInelasticProcess.hh"
> #include "G4NeutronInelasticProcess.hh"
> #include "G4AntiNeutronInelasticProcess.hh"
> #include "G4LambdaInelasticProcess.hh"
> #include "G4AntiLambdaInelasticProcess.hh"
> #include "G4SigmaPlusInelasticProcess.hh"
> #include "G4SigmaMinusInelasticProcess.hh"
> #include "G4AntiSigmaPlusInelasticProcess.hh"
> #include "G4AntiSigmaMinusInelasticProcess.hh"
> #include "G4XiZeroInelasticProcess.hh"
> #include "G4XiMinusInelasticProcess.hh"
> #include "G4AntiXiZeroInelasticProcess.hh"
> #include "G4AntiXiMinusInelasticProcess.hh"
> #include "G4DeuteronInelasticProcess.hh"
> #include "G4TritonInelasticProcess.hh"
> #include "G4AlphaInelasticProcess.hh"
> #include "G4OmegaMinusInelasticProcess.hh"
> #include "G4AntiOmegaMinusInelasticProcess.hh"
> 
> // Low-energy Models
> 
> #include "G4LElastic.hh"
> #include "G4LFission.hh"
> #include "G4LCapture.hh"
> 
> #include "G4LEPionPlusInelastic.hh"
> #include "G4LEPionMinusInelastic.hh"
> #include "G4LEKaonPlusInelastic.hh"
> #include "G4LEKaonZeroSInelastic.hh"
> #include "G4LEKaonZeroLInelastic.hh"
> #include "G4LEKaonMinusInelastic.hh"
> #include "G4LEProtonInelastic.hh"
> #include "G4LEAntiProtonInelastic.hh"
> #include "G4LENeutronInelastic.hh"
> #include "G4LEAntiNeutronInelastic.hh"
> #include "G4LELambdaInelastic.hh"
> #include "G4LEAntiLambdaInelastic.hh"
> #include "G4LESigmaPlusInelastic.hh"
> #include "G4LESigmaMinusInelastic.hh"
> #include "G4LEAntiSigmaPlusInelastic.hh"
> #include "G4LEAntiSigmaMinusInelastic.hh"
> #include "G4LEXiZeroInelastic.hh"
> #include "G4LEXiMinusInelastic.hh"
> #include "G4LEAntiXiZeroInelastic.hh"
> #include "G4LEAntiXiMinusInelastic.hh"
> #include "G4LEDeuteronInelastic.hh"
> #include "G4LETritonInelastic.hh"
> #include "G4LEAlphaInelastic.hh"
> #include "G4LEOmegaMinusInelastic.hh"
> #include "G4LEAntiOmegaMinusInelastic.hh"
> 
> // -- generator models
> 
> #include "G4TheoFSGenerator.hh"
> #include "G4ExcitationHandler.hh"
> #include "G4Evaporation.hh"
> #include "G4CompetitiveFission.hh"
> #include "G4FermiBreakUp.hh"
> #include "G4StatMF.hh"
> #include "G4GeneratorPrecompoundInterface.hh"
> #include "G4Fancy3DNucleus.hh"
> #include "G4LEProtonInelastic.hh"
> #include "G4StringModel.hh"
> #include "G4PreCompoundModel.hh"
> #include "G4FTFModel.hh"
> #include "G4QGSMFragmentation.hh"
> #include "G4ExcitedStringDecay.hh"
> 
> // // ConstructHad() // // Makes discrete physics processes for the
> hadrons, at present limited // to those particles with GHEISHA
> interactions (INTRC > 0). // The processes are: Elastic scattering,
> Inelastic scattering, // Fission (for neutron only), and Capture
> (neutron). // // F.W.Jones 06-JUL-1998 //
> 
> void TriumfPhysicsList::ConstructHad()
> {
>     // this will be the model class for high energies
>     G4TheoFSGenerator * theTheoModel = new G4TheoFSGenerator;
> 
>     // all models for treatment of thermal nucleus 
>     G4Evaporation * theEvaporation = new G4Evaporation;
>     G4FermiBreakUp * theFermiBreakUp = new G4FermiBreakUp;
>     G4StatMF * theMF = new G4StatMF;
> 
>     // Evaporation logic
>     G4ExcitationHandler * theHandler = new G4ExcitationHandler;
>         theHandler->SetEvaporation(theEvaporation);
>         theHandler->SetFermiModel(theFermiBreakUp);
>         theHandler->SetMultiFragmentation(theMF);
>         theHandler->SetMaxAandZForFermiBreakUp(12, 6);
>         theHandler->SetMinEForMultiFrag(3*MeV);
> 
>     // Pre equilibrium stage 
>     G4PreCompoundModel * thePreEquilib = new G4PreCompoundModel(theHandler);
> 
>     // a no-cascade generator-precompound interaface
>     G4GeneratorPrecompoundInterface * theCascade = new G4GeneratorPrecompoundInterface;
>             theCascade->SetDeExcitation(thePreEquilib);  
> 
>     // here come the high energy parts
>     // the string model; still not quite according to design - Explicite use of the forseen interfaces 
>     // will be tested and documented in this program by beta-02 at latest.
>     G4VPartonStringModel * theStringModel;
>     theStringModel = new G4FTFModel;
>     theTheoModel->SetTransport(theCascade);
>     theTheoModel->SetHighEnergyGenerator(theStringModel);
>     theTheoModel->SetMinEnergy(19*GeV);
>     theTheoModel->SetMaxEnergy(100*TeV);
> 
>       G4VLongitudinalStringDecay * theFragmentation = new G4QGSMFragmentation;
>       G4ExcitedStringDecay * theStringDecay = new G4ExcitedStringDecay(theFragmentation);
>       theStringModel->SetFragmentationModel(theStringDecay);
> 
> // done with the generator model (most of the above is also available as default)
>    G4HadronElasticProcess* theElasticProcess = 
>                                     new G4HadronElasticProcess;
>    G4LElastic* theElasticModel = new G4LElastic;
>    theElasticProcess->RegisterMe(theElasticModel);
>    G4HadronElasticProcess* theElasticProcess1 = 
>                                     new G4HadronElasticProcess;
>    theParticleIterator->reset();
>    while ((*theParticleIterator)()) {
>       G4ParticleDefinition* particle = theParticleIterator->value();
>       G4ProcessManager* pmanager = particle->GetProcessManager();
>       G4String particleName = particle->GetParticleName();
> 
>       if (particleName == "pi+") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4PionPlusInelasticProcess* theInelasticProcess = 
>                                 new G4PionPlusInelasticProcess("inelastic");
>          G4LEPionPlusInelastic* theInelasticModel = 
>                                 new G4LEPionPlusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "pi-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4PionMinusInelasticProcess* theInelasticProcess = 
>                                 new G4PionMinusInelasticProcess("inelastic");
>          G4LEPionMinusInelastic* theInelasticModel = 
>                                 new G4LEPionMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "kaon+") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4KaonPlusInelasticProcess* theInelasticProcess = 
>                                   new G4KaonPlusInelasticProcess("inelastic");
>          G4LEKaonPlusInelastic* theInelasticModel = new G4LEKaonPlusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "kaon0S") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4KaonZeroSInelasticProcess* theInelasticProcess = 
>                              new G4KaonZeroSInelasticProcess("inelastic");
>          G4LEKaonZeroSInelastic* theInelasticModel = 
>                              new G4LEKaonZeroSInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "kaon0L") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4KaonZeroLInelasticProcess* theInelasticProcess = 
>                              new G4KaonZeroLInelasticProcess("inelastic");
>          G4LEKaonZeroLInelastic* theInelasticModel = 
>                              new G4LEKaonZeroLInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "kaon-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4KaonMinusInelasticProcess* theInelasticProcess = 
>                                  new G4KaonMinusInelasticProcess("inelastic");
>          G4LEKaonMinusInelastic* theInelasticModel = 
>                                  new G4LEKaonMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "proton") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4ProtonInelasticProcess* theInelasticProcess = 
>                                     new G4ProtonInelasticProcess("inelastic");
>          G4LEProtonInelastic* theInelasticModel = new G4LEProtonInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_proton") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiProtonInelasticProcess* theInelasticProcess = 
>                                 new G4AntiProtonInelasticProcess("inelastic");
>          G4LEAntiProtonInelastic* theInelasticModel = 
>                                 new G4LEAntiProtonInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "neutron") {
> 
>           // elastic scattering
>          G4LElastic* theElasticModel1 = new G4LElastic;
>          theElasticProcess1->RegisterMe(theElasticModel1);
>          pmanager->AddDiscreteProcess(theElasticProcess1);
>           // inelastic scattering
>          G4NeutronInelasticProcess* theInelasticProcess = 
>                                     new G4NeutronInelasticProcess("inelastic");
>          G4LENeutronInelastic* theInelasticModel = new G4LENeutronInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>           // fission
>          G4HadronFissionProcess* theFissionProcess =
>                                     new G4HadronFissionProcess;
>          G4LFission* theFissionModel = new G4LFission;
>          theFissionProcess->RegisterMe(theFissionModel);
>          pmanager->AddDiscreteProcess(theFissionProcess);
>          // capture
>          G4HadronCaptureProcess* theCaptureProcess =
>                                     new G4HadronCaptureProcess;
>          G4LCapture* theCaptureModel = new G4LCapture;
>          theCaptureProcess->RegisterMe(theCaptureModel);
>          pmanager->AddDiscreteProcess(theCaptureProcess);
>                  G4cout << "Neutron Created" << G4endl;
>       }  
>       else if (particleName == "anti_neutron") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiNeutronInelasticProcess* theInelasticProcess = 
>                                new G4AntiNeutronInelasticProcess("inelastic");
>          G4LEAntiNeutronInelastic* theInelasticModel = 
>                                new G4LEAntiNeutronInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "lambda") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4LambdaInelasticProcess* theInelasticProcess = 
>                                     new G4LambdaInelasticProcess("inelastic");
>          G4LELambdaInelastic* theInelasticModel = new G4LELambdaInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_lambda") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiLambdaInelasticProcess* theInelasticProcess = 
>                                 new G4AntiLambdaInelasticProcess("inelastic");
>          G4LEAntiLambdaInelastic* theInelasticModel = 
>                                 new G4LEAntiLambdaInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "sigma+") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4SigmaPlusInelasticProcess* theInelasticProcess = 
>                                  new G4SigmaPlusInelasticProcess("inelastic");
>          G4LESigmaPlusInelastic* theInelasticModel = 
>                                  new G4LESigmaPlusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "sigma-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4SigmaMinusInelasticProcess* theInelasticProcess = 
>                                  new G4SigmaMinusInelasticProcess("inelastic");
>          G4LESigmaMinusInelastic* theInelasticModel = 
>                                  new G4LESigmaMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_sigma+") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiSigmaPlusInelasticProcess* theInelasticProcess = 
>                              new G4AntiSigmaPlusInelasticProcess("inelastic");
>          G4LEAntiSigmaPlusInelastic* theInelasticModel = 
>                                  new G4LEAntiSigmaPlusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_sigma-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiSigmaMinusInelasticProcess* theInelasticProcess = 
>                             new G4AntiSigmaMinusInelasticProcess("inelastic");
>          G4LEAntiSigmaMinusInelastic* theInelasticModel = 
>                                  new G4LEAntiSigmaMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "xi0") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4XiZeroInelasticProcess* theInelasticProcess = 
>                             new G4XiZeroInelasticProcess("inelastic");
>          G4LEXiZeroInelastic* theInelasticModel = 
>                                  new G4LEXiZeroInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "xi-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4XiMinusInelasticProcess* theInelasticProcess = 
>                             new G4XiMinusInelasticProcess("inelastic");
>          G4LEXiMinusInelastic* theInelasticModel = 
>                                  new G4LEXiMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_xi0") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiXiZeroInelasticProcess* theInelasticProcess = 
>                             new G4AntiXiZeroInelasticProcess("inelastic");
>          G4LEAntiXiZeroInelastic* theInelasticModel = 
>                                  new G4LEAntiXiZeroInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_xi-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiXiMinusInelasticProcess* theInelasticProcess = 
>                             new G4AntiXiMinusInelasticProcess("inelastic");
>          G4LEAntiXiMinusInelastic* theInelasticModel = 
>                                  new G4LEAntiXiMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "deuteron") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4DeuteronInelasticProcess* theInelasticProcess = 
>                             new G4DeuteronInelasticProcess("inelastic");
>          G4LEDeuteronInelastic* theInelasticModel = 
>                                  new G4LEDeuteronInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "triton") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4TritonInelasticProcess* theInelasticProcess = 
>                             new G4TritonInelasticProcess("inelastic");
>          G4LETritonInelastic* theInelasticModel = 
>                                  new G4LETritonInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "alpha") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AlphaInelasticProcess* theInelasticProcess = 
>                             new G4AlphaInelasticProcess("inelastic");
>          G4LEAlphaInelastic* theInelasticModel = 
>                                  new G4LEAlphaInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "omega-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4OmegaMinusInelasticProcess* theInelasticProcess = 
>                             new G4OmegaMinusInelasticProcess("inelastic");
>          G4LEOmegaMinusInelastic* theInelasticModel = 
>                                  new G4LEOmegaMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>       else if (particleName == "anti_omega-") {
>          pmanager->AddDiscreteProcess(theElasticProcess);
>          G4AntiOmegaMinusInelasticProcess* theInelasticProcess = 
>                             new G4AntiOmegaMinusInelasticProcess("inelastic");
>          G4LEAntiOmegaMinusInelastic* theInelasticModel = 
>                                  new G4LEAntiOmegaMinusInelastic;
>          theInelasticProcess->RegisterMe(theInelasticModel);
>          theInelasticProcess->RegisterMe(theTheoModel);
>          pmanager->AddDiscreteProcess(theInelasticProcess);
>       }
>    }
> }
> 
> void TriumfPhysicsList::ConstructLeptHad() {;}
> 
> #include "G4Decay.hh"
> void TriumfPhysicsList::ConstructGeneral()
> {
>   G4Decay* theDecayProcess = new G4Decay();
>   theParticleIterator->reset();
>   while( (*theParticleIterator)() ){
>     G4ParticleDefinition* particle = theParticleIterator->value();
>     G4ProcessManager* pmanager = particle->GetProcessManager();
>     if (theDecayProcess->IsApplicable(*particle)) { 
>       pmanager ->AddProcess(theDecayProcess);
>       pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
>       pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
>     }
>   }
> }
> 
> void TriumfPhysicsList::SetCuts()
> {
>   if (verboseLevel >0)
>   {
>     G4cout << "TriumfPhysicsList::SetCuts:";
>     G4cout << "CutLength : " << defaultCutValue/mm << " (mm)" << G4endl;
>   }  
>   //   "G4VUserPhysicsList::SetCutsWithDefault" method sets 
>   //   the default cut value for all particle types 
>   SetCutValue(1*mm,"gamma");  
>   SetCutValue(1*mm,"e+");
>   SetCutValue(1*mm,"e-");
>   SetCutValue(1*mm,"proton");
>   SetCutValue(0.1*eV,"neutron");
> }
> 

Hello,

I personally would move out this structure of Physics List to use Physics 
Lists components from physics_list library, which is now a part of the 
toolkit. You can take one of the existing Physics Lists (QGSP_BRT, 
QGSP_BIC, QGSP_BIC_HP, QGSP_BERT_HP....) or use an approach shown in 
$G4INSTALL/examples/extended/hadronic/Hadr01 

VI

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

 Add Message Add Message
to: "Re: Near thermal Neutron simulation problem"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]