Ignore:
Timestamp:
Apr 20, 2009, 5:54:05 PM (15 years ago)
Author:
garnier
Message:

update to geant4.9.2

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/source/processes/hadronic/models/pre_equilibrium/exciton_model/src/G4PreCompoundDeuteron.cc

    r968 r1007  
    2424// ********************************************************************
    2525//
    26 // $Id: G4PreCompoundDeuteron.cc,v 1.5 2009/02/13 18:57:32 vnivanch Exp $
    27 // GEANT4 tag $Name: geant4-09-02-ref-02 $
    28 //
    29 // -------------------------------------------------------------------
    30 //
    31 // GEANT4 Class file
    32 //
    33 //
    34 // File name:     G4PreCompoundDeuteron
    35 //
    36 // Author:         V.Lara
    37 //
    38 // Modified: 
    39 // 21.08.2008 J. M. Quesada add choice of options 
    40 // 10.02.2009 J. M. Quesada set default opt1 
    41 //
     26//J.M.Quesada (August 08). New source file
     27//
     28// Modif (21 August 2008) by J. M. Quesada for external choice of inverse
     29// cross section option
    4230
    4331#include "G4PreCompoundDeuteron.hh"
    4432
    45 G4ReactionProduct * G4PreCompoundDeuteron::GetReactionProduct() const
    46 {
    47   G4ReactionProduct * theReactionProduct =
    48     new G4ReactionProduct(G4Deuteron::DeuteronDefinition());
    49   theReactionProduct->SetMomentum(GetMomentum().vect());
    50   theReactionProduct->SetTotalEnergy(GetMomentum().e());
     33
     34
     35  G4ReactionProduct * G4PreCompoundDeuteron::GetReactionProduct() const
     36  {
     37    G4ReactionProduct * theReactionProduct =
     38      new G4ReactionProduct(G4Deuteron::DeuteronDefinition());
     39    theReactionProduct->SetMomentum(GetMomentum().vect());
     40    theReactionProduct->SetTotalEnergy(GetMomentum().e());
    5141#ifdef PRECOMPOUND_TEST
    52   theReactionProduct->SetCreatorModel("G4PrecompoundModel");
     42    theReactionProduct->SetCreatorModel("G4PrecompoundModel");
    5343#endif
    54   return theReactionProduct;
    55 }   
    56 
    57  
    58 G4double G4PreCompoundDeuteron::FactorialFactor(const G4double N, const G4double P)
    59 {
    60   return (N-1.0)*(N-2.0)*(P-1.0)*P/2.0;
    61 }
     44    return theReactionProduct;
     45  }   
     46
     47 
     48   G4double G4PreCompoundDeuteron::FactorialFactor(const G4double N, const G4double P)
     49  {
     50      return (N-1.0)*(N-2.0)*(P-1.0)*P/2.0;
     51  }
    6252 
    63 G4double G4PreCompoundDeuteron::CoalescenceFactor(const G4double A)
    64 {
    65   return 16.0/A;
    66 }   
    67 
    68 G4double G4PreCompoundDeuteron::GetRj(const G4int NumberParticles, const G4int NumberCharged)
    69 {
    70   G4double rj = 0.0;
    71   G4double denominator = NumberParticles*(NumberParticles-1);
    72   if(NumberCharged >=1 && (NumberParticles-NumberCharged) >=1) {
    73     rj = 2.0*static_cast<G4double>(NumberCharged*(NumberParticles-NumberCharged))
    74       / static_cast<G4double>(denominator);
    75   }
    76   return rj;
    77 }
     53   G4double G4PreCompoundDeuteron::CoalescenceFactor(const G4double A)
     54  {
     55    return 16.0/A;
     56  }   
     57
     58  G4double G4PreCompoundDeuteron::GetRj(const G4int NumberParticles, const G4int NumberCharged)
     59  {
     60    G4double rj = 0.0;
     61    G4double denominator = NumberParticles*(NumberParticles-1);
     62   if(NumberCharged >=1 && (NumberParticles-NumberCharged) >=1) rj = 2.0*static_cast<G4double>(NumberCharged*(NumberParticles-NumberCharged))/static_cast<G4double>(denominator);
     63
     64    return rj;
     65  }
    7866
    7967////////////////////////////////////////////////////////////////////////////////////
    … …  
    8371//OPT=3,4 Kalbach's parameterization
    8472//
    85 G4double G4PreCompoundDeuteron::CrossSection(const  G4double K)
     73 G4double G4PreCompoundDeuteron::CrossSection(const  G4double K)
    8674{
    8775
    … …  
    118106//---------
    119107//
    120 G4double G4PreCompoundDeuteron::GetAlpha()
    121 {
    122   G4double C = 0.0;
    123   G4double aZ = GetZ() + GetRestZ();
    124   if (aZ >= 70)
    125     {
    126       C = 0.10;
    127     }
    128   else
    129     {
    130       C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375;
    131     }
    132   return 1.0 + C/2.0;
    133 }
     108 G4double G4PreCompoundDeuteron::GetAlpha()
     109  {
     110G4double C = 0.0;
     111    G4double aZ = GetZ() + GetRestZ();
     112    if (aZ >= 70)
     113      {
     114        C = 0.10;
     115      }
     116    else
     117      {
     118        C = ((((0.15417e-06*aZ) - 0.29875e-04)*aZ + 0.21071e-02)*aZ - 0.66612e-01)*aZ + 0.98375;
     119      }
     120    return 1.0 + C/2.0;
     121  }
    134122//
    135123//---------
    136124//
    137 G4double G4PreCompoundDeuteron::GetBeta()
    138 {
    139   return -GetCoulombBarrier();
    140 }
     125  G4double G4PreCompoundDeuteron::GetBeta()
     126  {
     127      return -GetCoulombBarrier();
     128  }
    141129//
    142130//********************* OPT=1,2 : Chatterjee's cross section ************************
    … …  
    152140
    153141  G4double landa ,mu ,nu ,p , Ec,q,r,ji,xs;
    154   //G4double Eo(0),epsilon1(0),epsilon2(0),discri(0);
     142//G4double Eo(0),epsilon1(0),epsilon2(0),discri(0);
    155143
    156144 
    … …  
    186174}
    187175
     176
     177
     178
    188179// *********** OPT=3,4 : Kalbach's cross sections (from PRECO code)*************
    189180G4double G4PreCompoundDeuteron::GetOpt34(const  G4double K)
    … …  
    214205  G4double     ra=0.80;
    215206       
    216   //JMQ 13/02/09 increase of reduced radius to lower the barrier
    217   // ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
    218   ec = 1.44 * theZ * ResidualZ / (1.7*ResidualAthrd+ra);
     207  ec = 1.44 * theZ * ResidualZ / (1.5*ResidualAthrd+ra);
    219208  ecsq = ec * ec;
    220209  p = p0 + p1/ec + p2/ecsq;
    … …  
    237226  elab = K * FragmentA / ResidualA;
    238227  sig = 0.;
    239 
     228 
    240229  if (elab <= ec) { //start for E<Ec
    241     if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;
     230    if (elab > ecut2)  sig = (p*elab*elab+a*elab+b) * signor;   
    242231  }           //end for E<Ec
    243232  else {           //start for E>Ec
Note: See TracChangeset for help on using the changeset viewer.