#include <G4Fancy3DNucleus.hh>

Inheritance diagram for G4Fancy3DNucleus:

Public Member Functions
G4double	AddExcitationEnergy (G4double)

void	CenterNucleons ()

std::pair< G4double, G4double >	ChooseImpactXandY (G4double maxImpact)

G4double	CoulombBarrier ()

void	DoLorentzBoost (const G4LorentzVector &theBoost)

void	DoLorentzBoost (const G4ThreeVector &theBeta)

void	DoLorentzContraction (const G4LorentzVector &theBoost)

void	DoLorentzContraction (const G4ThreeVector &theBeta)

void	DoTranslation (const G4ThreeVector &theShift)

	G4Fancy3DNucleus ()

G4int	GetCharge ()

G4double	GetExcitationEnergy ()

G4double	GetMass ()

G4int	GetMassNumber ()

G4Nucleon *	GetNextNucleon ()

const G4VNuclearDensity *	GetNuclearDensity () const

G4double	GetNuclearRadius ()

G4double	GetNuclearRadius (const G4double maxRelativeDensity)

const std::vector< G4Nucleon > &	GetNucleons ()

G4int	GetNumberOfLambdas ()

G4double	GetOuterRadius ()

void	Init (G4int theA, G4int theZ, G4int numberOfLambdas=0)

std::pair< G4double, G4double >	RefetchImpactXandY ()

void	SortNucleonsDecZ ()

void	SortNucleonsIncZ ()

G4bool	StartLoop ()

	~G4Fancy3DNucleus ()

Private Member Functions
G4double	BindingEnergy ()

void	ChooseFermiMomenta ()

void	ChooseNucleons ()

void	ChoosePositions ()

	G4Fancy3DNucleus (const G4Fancy3DNucleus &right)

G4bool	operator!= (const G4Fancy3DNucleus &right) const

const G4Fancy3DNucleus &	operator= (const G4Fancy3DNucleus &right)

G4bool	operator== (const G4Fancy3DNucleus &right) const

G4bool	ReduceSum ()

Private Attributes
G4int	currentNucleon

G4double	excitationEnergy

std::vector< G4double >	fermiM

std::vector< G4ThreeVector >	momentum

G4int	myA

G4int	myL

G4int	myZ

G4double	nucleondistance

std::vector< G4ThreeVector >	places

std::vector< G4Fancy3DNucleusHelper >	testSums

G4VNuclearDensity *	theDensity

G4FermiMomentum	theFermi

std::pair< G4double, G4double >	theImpactParameter

std::vector< G4Nucleon >	theNucleons

Detailed Description

Definition at line 54 of file G4Fancy3DNucleus.hh.

Constructor & Destructor Documentation

◆ G4Fancy3DNucleus() [1/2]

G4Fancy3DNucleus::G4Fancy3DNucleus ( )

Definition at line 60 of file G4Fancy3DNucleus.cc.

  : myA(0), myZ(0), myL(0), theNucleons(250), currentNucleon(-1), theDensity(0), 
    nucleondistance(0.8*fermi),excitationEnergy(0.),
    places(250), momentum(250), fermiM(250), testSums(250)
{
}

◆ ~G4Fancy3DNucleus()

G4Fancy3DNucleus::~G4Fancy3DNucleus ( )

Definition at line 67 of file G4Fancy3DNucleus.cc.

{
  if(theDensity) delete theDensity;
}

References theDensity.

◆ G4Fancy3DNucleus() [2/2]

G4Fancy3DNucleus::G4Fancy3DNucleus ( const G4Fancy3DNucleus & right )

private

Member Function Documentation

◆ AddExcitationEnergy()

G4double G4Fancy3DNucleus::AddExcitationEnergy ( G4double anE )

inline

Definition at line 138 of file G4Fancy3DNucleus.hh.

{
   excitationEnergy +=anE;
   return excitationEnergy;
}

References excitationEnergy.

◆ BindingEnergy()

G4double G4Fancy3DNucleus::BindingEnergy ( )

private

Definition at line 173 of file G4Fancy3DNucleus.cc.

{
  return G4NucleiProperties::GetBindingEnergy(myA,myZ);
}

References G4NucleiProperties::GetBindingEnergy(), myA, and myZ.

Referenced by ChooseFermiMomenta(), GetMass(), and Init().

◆ CenterNucleons()

void G4Fancy3DNucleus::CenterNucleons ( )

Definition at line 251 of file G4Fancy3DNucleus.cc.

{
    G4ThreeVector   center;
    
    for (G4int i=0; i<myA; i++ )
    {
       center+=theNucleons[i].GetPosition();
    }   
    center /= -myA;
    DoTranslation(center);
}

References DoTranslation(), myA, and theNucleons.

Referenced by Init(), and G4BinaryLightIonReaction::Interact().

◆ ChooseFermiMomenta()

void G4Fancy3DNucleus::ChooseFermiMomenta ( )

private

Definition at line 457 of file G4Fancy3DNucleus.cc.

{
    G4int i;
    G4double density;
 
    // Pre-allocate buffers for filling by index
    momentum.resize(myA, G4ThreeVector(0.,0.,0.));
    fermiM.resize(myA, 0.*GeV);
 
    for (G4int ntry=0; ntry<1 ; ntry ++ )
    {
    for (i=0; i < myA; i++ )    // momenta for all, including last, in case we swap nucleons
    {
       density = theDensity->GetDensity(theNucleons[i].GetPosition());
       fermiM[i] = theFermi.GetFermiMomentum(density);
       G4ThreeVector mom=theFermi.GetMomentum(density);
       if (theNucleons[i].GetDefinition() == G4Proton::Proton())
       {
          G4double eMax = std::sqrt(sqr(fermiM[i]) +sqr(theNucleons[i].GetDefinition()->GetPDGMass()) )
                          - CoulombBarrier();
          if ( eMax > theNucleons[i].GetDefinition()->GetPDGMass() )
          {
              G4double pmax2= sqr(eMax) - sqr(theNucleons[i].GetDefinition()->GetPDGMass());
          fermiM[i] = std::sqrt(pmax2);
          while ( mom.mag2() > pmax2 )  /* Loop checking, 30-Oct-2015, G.Folger */
          {
              mom=theFermi.GetMomentum(density, fermiM[i]);
          }
          }  else
          {
                  //AR-21Dec2017 : emit a "JustWarning" exception instead of writing on the error stream.
              //G4cerr << "G4Fancy3DNucleus: difficulty finding proton momentum" << G4endl;
                  G4ExceptionDescription ed;
                  ed << "Nucleus Z A " << myZ << " " << myA << G4endl;
                  ed << "proton with eMax=" << eMax << G4endl;
                  G4Exception( "G4Fancy3DNucleus::ChooseFermiMomenta(): difficulty finding proton momentum, set it to (0,0,0)",
                               "HAD_FANCY3DNUCLEUS_001", JustWarning, ed );
          mom=G4ThreeVector(0,0,0);
          }
 
       }
       momentum[i]= mom;
    }
 
    if ( ReduceSum() ) break;
//       G4cout <<" G4FancyNucleus: iterating to find momenta: "<< ntry<< G4endl;
    }
 
//     G4ThreeVector sum;
//     for (G4int index=0; index<myA;sum+=momentum[index++])
//     ;
//     G4cout << "final sum / mag() " << sum << " / " << sum.mag() << G4endl;
 
    G4double energy;
    for ( i=0; i< myA ; i++ )
    {
       energy = theNucleons[i].GetParticleType()->GetPDGMass()
            - BindingEnergy()/myA;
       G4LorentzVector tempV(momentum[i],energy);
       theNucleons[i].SetMomentum(tempV);
       // GF 11-05-2011: set BindingEnergy to be T of Nucleon with p , ~ p**2/2m
       //theNucleons[i].SetBindingEnergy(
       //     0.5*sqr(fermiM[i])/theNucleons[i].GetParticleType()->GetPDGMass());
    }
}

References BindingEnergy(), CoulombBarrier(), G4INCL::KinematicsUtils::energy(), fermiM, G4endl, G4Exception(), G4VNuclearDensity::GetDensity(), G4FermiMomentum::GetFermiMomentum(), G4FermiMomentum::GetMomentum(), GeV, JustWarning, CLHEP::Hep3Vector::mag2(), momentum, myA, myZ, G4Proton::Proton(), ReduceSum(), sqr(), theDensity, theFermi, and theNucleons.

Referenced by Init().

◆ ChooseImpactXandY()

std::pair< G4double, G4double > G4V3DNucleus::ChooseImpactXandY ( G4double maxImpact )

inlineinherited

Definition at line 86 of file G4V3DNucleus.hh.

{
  G4double x,y;
  do
  {
    x = 2*G4UniformRand() - 1;
    y = 2*G4UniformRand() - 1;
  }
  while(x*x + y*y > 1);  /* Loop checking, 30-Oct-2015, G.Folger */
 
  G4double impactX = x*(maxImpact); 
  G4double impactY = y*(maxImpact);
  theImpactParameter.first = impactX;
  theImpactParameter.second = impactY;
  return theImpactParameter;
}

References G4UniformRand, and G4V3DNucleus::theImpactParameter.

Referenced by G4QGSParticipants::GetList(), and G4FTFParticipants::GetList().

◆ ChooseNucleons()

void G4Fancy3DNucleus::ChooseNucleons ( )

private

Definition at line 280 of file G4Fancy3DNucleus.cc.

{
  G4int protons=0, nucleons=0, lambdas=0;
  G4double probProton = ( G4double(myZ) )/( G4double(myA) );
  G4double probLambda = myL > 0 ? ( G4double(myL) )/( G4double(myA) ) : 0.0;
  while ( nucleons < myA ) {  /* Loop checking, 30-Oct-2015, G.Folger */
    G4double rnd = G4UniformRand();
    if ( rnd < probProton ) {
      if ( protons < myZ ) {
    protons++;
    theNucleons[nucleons++].SetParticleType(G4Proton::Proton());
      }
    } else if ( rnd < probProton + probLambda ) {
      if ( lambdas < myL ) {
    lambdas++;
    theNucleons[nucleons++].SetParticleType(G4Lambda::Lambda());
      }
    } else {
      if ( (nucleons - protons - lambdas) < (myA - myZ - myL) ) {
    theNucleons[nucleons++].SetParticleType(G4Neutron::Neutron());
      }
    }
  }
  return;
}

References G4UniformRand, G4Lambda::Lambda(), myA, myL, myZ, G4Neutron::Neutron(), G4Proton::Proton(), and theNucleons.

Referenced by Init().

◆ ChoosePositions()

void G4Fancy3DNucleus::ChoosePositions ( )

private

Definition at line 306 of file G4Fancy3DNucleus.cc.

{
  if( myA != 12) {
 
    G4int i=0;
    G4ThreeVector   aPos, delta;
    G4bool      freeplace;
    const G4double nd2=sqr(nucleondistance);
    G4double maxR=GetNuclearRadius(0.001);   //  there are no nucleons at a
                                         //  relative Density of 0.01
    G4int jr=0;
    G4int jx,jy;
    G4double arand[600];
    G4double *prand=arand;
    places.clear();             // Reset data buffer
    G4int interationsLeft=1000*myA;
    while ( (i < myA) && (--interationsLeft>0))  /* Loop checking, 30-Oct-2015, G.Folger */
    {
      do
      {     
        if ( jr < 3 ) 
    {
          jr=std::min(600,9*(myA - i));
          G4RandFlat::shootArray(jr,prand);
      //CLHEP::RandFlat::shootArray(jr, prand );
    }
    jx=--jr;
    jy=--jr;
    aPos.set((2*arand[jx]-1.), (2*arand[jy]-1.), (2*arand[--jr]-1.));
      } while (aPos.mag2() > 1. );  /* Loop checking, 30-Oct-2015, G.Folger */
      aPos *=maxR;
      G4double density=theDensity->GetRelativeDensity(aPos);
      if (G4UniformRand() < density)
      {
        freeplace= true;
    std::vector<G4ThreeVector>::iterator iplace;
    for( iplace=places.begin(); iplace!=places.end() && freeplace;++iplace)
    {
      delta = *iplace - aPos;
      freeplace= delta.mag2() > nd2;
    }      
    if ( freeplace ) {
          G4double pFermi=theFermi.GetFermiMomentum(theDensity->GetDensity(aPos));
          // protons must at least have binding energy of CoulombBarrier, so
          //  assuming the Fermi energy corresponds to a potential, we must place these such
          //  that the Fermi Energy > CoulombBarrier
      if (theNucleons[i].GetDefinition() == G4Proton::Proton())
          {
            G4double nucMass = theNucleons[i].GetDefinition()->GetPDGMass();
            G4double eFermi= std::sqrt( sqr(pFermi) + sqr(nucMass) ) - nucMass;
        if (eFermi <= CoulombBarrier() ) freeplace=false;
      }
    } 
    if ( freeplace ) {
          theNucleons[i].SetPosition(aPos);
      places.push_back(aPos);
      ++i;
    }
      }
    }
    if (interationsLeft<=0) {
      G4Exception("model/util/G4Fancy3DNucleus.cc", "mod_util001", FatalException,
                  "Problem to place nucleons");
    }
 
  } else {
    // Start insertion
    // Alpha cluster structure of carbon nuclei, C-12, is implemented according to
    //       P. Bozek, W. Broniowski, E.R. Arriola and M. Rybczynski
    //                Phys. Rev. C90, 064902 (2014) 
    const G4double Lbase=3.05*fermi;
    const G4double Disp=0.552;        // 0.91^2*2/3 fermi^2
    const G4double nd2=sqr(nucleondistance);
    const G4ThreeVector Corner1=G4ThreeVector( Lbase/2.,         0., 0.);
    const G4ThreeVector Corner2=G4ThreeVector(-Lbase/2.,         0., 0.);
    const G4ThreeVector Corner3=G4ThreeVector(       0.,Lbase*0.866, 0.);  // 0.866=sqrt(3)/2
    G4ThreeVector R1;
    R1=G4ThreeVector(G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp))*fermi + Corner1;
    theNucleons[0].SetPosition(R1); // First nucleon of the first He-4
    G4int loopCounterLeft = 10000;
    for(G4int ii=1; ii<4; ii++)     // 2 - 4 nucleons of the first He-4
    {
      G4bool Continue;
      do
      {
        R1=G4ThreeVector(G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp))*fermi + Corner1;
        theNucleons[ii].SetPosition(R1);
        Continue=false;
        for(G4int jj=0; jj < ii; jj++)
        {
          if( (theNucleons[ii].GetPosition() - theNucleons[jj].GetPosition()).mag2() <= nd2 ) {Continue = true; break;}
        }
      } while( Continue && --loopCounterLeft > 0 );  /* Loop checking, 12-Dec-2017, A.Ribon */
    }
    if ( loopCounterLeft <= 0 ) {
      G4Exception("model/util/G4Fancy3DNucleus.cc", "mod_util002", FatalException,
                  "Unable to find a good position for the first alpha cluster");
    }
    loopCounterLeft = 10000;
    for(G4int ii=4; ii<8; ii++)     // 5 - 8 nucleons of the second He-4
    {
      G4bool Continue;
      do
      {
        R1=G4ThreeVector(G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp))*fermi + Corner2;
        theNucleons[ii].SetPosition(R1);
        Continue=false;
        for(G4int jj=0; jj < ii; jj++)
        {
          if( (theNucleons[ii].GetPosition() - theNucleons[jj].GetPosition()).mag2() <= nd2 ) {Continue = true; break;}
        }
      } while( Continue && --loopCounterLeft > 0 );  /* Loop checking, 12-Dec-2017, A.Ribon */
    }
    if ( loopCounterLeft <= 0 ) {
      G4Exception("model/util/G4Fancy3DNucleus.cc", "mod_util003", FatalException,
                  "Unable to find a good position for the second alpha cluster");
    }
    loopCounterLeft = 10000;
    for(G4int ii=8; ii<12; ii++)    // 9 - 12 nucleons of the third He-4
    {
      G4bool Continue;
      do
      {
        R1=G4ThreeVector(G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp), G4RandGauss::shoot(0.,Disp))*fermi + Corner3;
        theNucleons[ii].SetPosition(R1);
        Continue=false;
        for(G4int jj=0; jj < ii; jj++)
        {
          if( (theNucleons[ii].GetPosition() - theNucleons[jj].GetPosition()).mag2() <= nd2 ) {Continue = true; break;}
        }
      } while( Continue && --loopCounterLeft > 0 );  /* Loop checking, 12-Dec-2017, A.Ribon */
    }
    if ( loopCounterLeft <= 0 ) {
      G4Exception("model/util/G4Fancy3DNucleus.cc", "mod_util004", FatalException,
                  "Unable to find a good position for the third alpha cluster");
    }
    G4LorentzRotation RandomRotation;
    RandomRotation.rotateZ(2.*pi*G4UniformRand());
    RandomRotation.rotateY(std::acos(2.*G4UniformRand()-1.));
    // Randomly rotation of the created nucleus
    G4LorentzVector Pos;
    for(G4int ii=0; ii<myA; ii++ )
    {
      Pos=G4LorentzVector(theNucleons[ii].GetPosition(),0.); Pos *=RandomRotation;
      G4ThreeVector NewPos = Pos.vect();
      theNucleons[ii].SetPosition(NewPos);
    }
 
  }
}

References CoulombBarrier(), FatalException, fermi, G4Exception(), G4UniformRand, G4VNuclearDensity::GetDensity(), G4FermiMomentum::GetFermiMomentum(), GetNuclearRadius(), G4VNuclearDensity::GetRelativeDensity(), CLHEP::Hep3Vector::mag2(), G4INCL::Math::min(), myA, nucleondistance, pi, places, G4Proton::Proton(), CLHEP::HepLorentzRotation::rotateY(), CLHEP::HepLorentzRotation::rotateZ(), CLHEP::Hep3Vector::set(), G4INCL::DeJongSpin::shoot(), sqr(), theDensity, theFermi, and theNucleons.

Referenced by Init().

◆ CoulombBarrier()

G4double G4Fancy3DNucleus::CoulombBarrier ( )

virtual

Implements G4V3DNucleus.

Definition at line 607 of file G4Fancy3DNucleus.cc.

{
  static const G4double cfactor = (1.44/1.14) * MeV;
  return cfactor*myZ/(1.0 + G4Pow::GetInstance()->Z13(myA));
}

References G4Pow::GetInstance(), MeV, myA, myZ, and G4Pow::Z13().

Referenced by ChooseFermiMomenta(), and ChoosePositions().

◆ DoLorentzBoost() [1/2]

void G4Fancy3DNucleus::DoLorentzBoost ( const G4LorentzVector & theBoost )

virtual

Implements G4V3DNucleus.

Definition at line 213 of file G4Fancy3DNucleus.cc.

{
    for (G4int i=0; i<myA; i++){
        theNucleons[i].Boost(theBoost);
    }
}

References myA, and theNucleons.

◆ DoLorentzBoost() [2/2]

void G4Fancy3DNucleus::DoLorentzBoost ( const G4ThreeVector & theBeta )

virtual

Implements G4V3DNucleus.

Definition at line 220 of file G4Fancy3DNucleus.cc.

{
    for (G4int i=0; i<myA; i++){
        theNucleons[i].Boost(theBeta);
    }
}

References myA, and theNucleons.

◆ DoLorentzContraction() [1/2]

void G4Fancy3DNucleus::DoLorentzContraction ( const G4LorentzVector & theBoost )

virtual

Implements G4V3DNucleus.

Definition at line 241 of file G4Fancy3DNucleus.cc.

{
    if (theBoost.e() !=0 ) {
       G4ThreeVector beta = theBoost.vect()/theBoost.e();
       DoLorentzContraction(beta);
    }
}

References anonymous_namespace{G4PionRadiativeDecayChannel.cc}::beta, DoLorentzContraction(), CLHEP::HepLorentzVector::e(), and CLHEP::HepLorentzVector::vect().

Referenced by DoLorentzContraction().

◆ DoLorentzContraction() [2/2]

void G4Fancy3DNucleus::DoLorentzContraction ( const G4ThreeVector & theBeta )

virtual

Implements G4V3DNucleus.

Definition at line 227 of file G4Fancy3DNucleus.cc.

{
    G4double beta2=theBeta.mag2();
    if (beta2 > 0) {
       G4double factor=(1-std::sqrt(1-beta2))/beta2; // (gamma-1)/gamma/beta**2
       G4ThreeVector rprime;
       for (G4int i=0; i< myA; i++) {
         rprime = theNucleons[i].GetPosition() - 
           factor * (theBeta*theNucleons[i].GetPosition()) * theBeta;  
         theNucleons[i].SetPosition(rprime);
       }
    }    
}

References CLHEP::Hep3Vector::mag2(), myA, and theNucleons.

◆ DoTranslation()

void G4Fancy3DNucleus::DoTranslation ( const G4ThreeVector & theShift )

virtual

Implements G4V3DNucleus.

Definition at line 263 of file G4Fancy3DNucleus.cc.

{
  G4ThreeVector tempV;
  for (G4int i=0; i<myA; i++ )
    {
      tempV = theNucleons[i].GetPosition() + theShift;
      theNucleons[i].SetPosition(tempV);
    }   
}

References myA, and theNucleons.

Referenced by CenterNucleons().

◆ GetCharge()

G4int G4Fancy3DNucleus::GetCharge ( )

inlinevirtual

Implements G4V3DNucleus.

Definition at line 123 of file G4Fancy3DNucleus.hh.

{
    return myZ;
}

References myZ.

Referenced by G4BinaryLightIonReaction::Interact().

◆ GetExcitationEnergy()

G4double G4Fancy3DNucleus::GetExcitationEnergy ( void )

inline

Definition at line 144 of file G4Fancy3DNucleus.hh.

{
   return excitationEnergy;
}

References excitationEnergy.

◆ GetMass()

G4double G4Fancy3DNucleus::GetMass ( )

virtual

Implements G4V3DNucleus.

Definition at line 203 of file G4Fancy3DNucleus.cc.

{
  if ( myL <= 0 ) return myZ*G4Proton::Proton()->GetPDGMass() + 
                         (myA-myZ)*G4Neutron::Neutron()->GetPDGMass() -
                         BindingEnergy();
  else            return G4HyperNucleiProperties::GetNuclearMass(myA, myZ, myL); 
}

References BindingEnergy(), G4HyperNucleiProperties::GetNuclearMass(), G4ParticleDefinition::GetPDGMass(), myA, myL, myZ, G4Neutron::Neutron(), and G4Proton::Proton().

Referenced by G4BinaryLightIonReaction::Interact().

◆ GetMassNumber()

G4int G4Fancy3DNucleus::GetMassNumber ( )

inlinevirtual

Implements G4V3DNucleus.

Definition at line 128 of file G4Fancy3DNucleus.hh.

{
    return myA;
}

References myA.

Referenced by G4BinaryLightIonReaction::Interact().

◆ GetNextNucleon()

G4Nucleon * G4Fancy3DNucleus::GetNextNucleon ( )

virtual

Implements G4V3DNucleus.

Definition at line 138 of file G4Fancy3DNucleus.cc.

{
  return ( (currentNucleon>=0 && currentNucleon<myA) ? 
       &theNucleons[currentNucleon++] : 0 );
}

References currentNucleon, myA, and theNucleons.

Referenced by G4LowEIonFragmentation::ApplyYourself(), G4BinaryLightIonReaction::GetProjectileExcitation(), and G4BinaryLightIonReaction::Interact().

◆ GetNuclearDensity()

const G4VNuclearDensity * G4Fancy3DNucleus::GetNuclearDensity ( ) const

virtual

Implements G4V3DNucleus.

Definition at line 273 of file G4Fancy3DNucleus.cc.

{
    return theDensity;
}

References theDensity.

Referenced by G4BinaryLightIonReaction::GetProjectileExcitation(), and G4BinaryLightIonReaction::Interact().

◆ GetNuclearRadius() [1/2]

G4double G4Fancy3DNucleus::GetNuclearRadius ( )

virtual

Implements G4V3DNucleus.

Definition at line 179 of file G4Fancy3DNucleus.cc.

{
    return GetNuclearRadius(0.5);
}

References GetNuclearRadius().

Referenced by ChoosePositions(), and GetNuclearRadius().

◆ GetNuclearRadius() [2/2]

G4double G4Fancy3DNucleus::GetNuclearRadius ( const G4double maxRelativeDensity )

virtual

Implements G4V3DNucleus.

Definition at line 184 of file G4Fancy3DNucleus.cc.

{   
    return theDensity->GetRadius(maxRelativeDensity);
}

References G4VNuclearDensity::GetRadius(), and theDensity.

◆ GetNucleons()

const std::vector< G4Nucleon > & G4Fancy3DNucleus::GetNucleons ( )

virtual

Implements G4V3DNucleus.

Definition at line 144 of file G4Fancy3DNucleus.cc.

{
  return theNucleons;
}

References theNucleons.

Referenced by G4MuMinusCapturePrecompound::ApplyYourself().

◆ GetNumberOfLambdas()

G4int G4Fancy3DNucleus::GetNumberOfLambdas ( )

inlinevirtual

Implements G4V3DNucleus.

Definition at line 133 of file G4Fancy3DNucleus.hh.

{
        return myL;
} 

References myL.

◆ GetOuterRadius()

G4double G4Fancy3DNucleus::GetOuterRadius ( )

virtual

Implements G4V3DNucleus.

Definition at line 189 of file G4Fancy3DNucleus.cc.

{
    G4double maxradius2=0;
    
    for (int i=0; i<myA; i++)
    {
       if ( theNucleons[i].GetPosition().mag2() > maxradius2 )
       {
        maxradius2=theNucleons[i].GetPosition().mag2();
       }
    }
    return std::sqrt(maxradius2)+nucleondistance;
}

References myA, nucleondistance, and theNucleons.

Referenced by G4LowEIonFragmentation::ApplyYourself(), and G4BinaryLightIonReaction::Interact().

◆ Init()

void G4Fancy3DNucleus::Init	(	G4int	theA,
		G4int	theZ,
		G4int	numberOfLambdas = `0`
	)

virtual

Implements G4V3DNucleus.

Definition at line 83 of file G4Fancy3DNucleus.cc.

{
  currentNucleon=-1;
  theNucleons.clear();
  nucleondistance = 0.8*fermi;
  places.clear();
  momentum.clear();
  fermiM.clear();
  testSums.clear();
 
  myZ = theZ;
  myA = theA;
  myL = std::max(numberOfLambdas, 0);  // Cannot be negative
  excitationEnergy=0;
 
  theNucleons.resize(myA);  // Pre-loads vector with empty elements
 
  // For simplicity, we neglect eventual Lambdas in the nucleus as far as the
  // density of nucler levels and the Fermi level are concerned.
  
  if(theDensity) delete theDensity;
  if ( myA < 17 ) {
     theDensity = new G4NuclearShellModelDensity(myA, myZ);
     if( myA == 12 ) nucleondistance=0.9*fermi; 
  } else {
     theDensity = new G4NuclearFermiDensity(myA, myZ);
  }
 
  theFermi.Init(myA, myZ);
  
  ChooseNucleons();
  
  ChoosePositions();
  
  if( myA == 12 ) CenterNucleons();   // This would introduce a bias
 
  ChooseFermiMomenta();
  
  G4double Ebinding= BindingEnergy()/myA;
  
  for (G4int aNucleon=0; aNucleon < myA; aNucleon++)
  {
    theNucleons[aNucleon].SetBindingEnergy(Ebinding);
  }
  
  return;
}

References BindingEnergy(), CenterNucleons(), ChooseFermiMomenta(), ChooseNucleons(), ChoosePositions(), currentNucleon, excitationEnergy, fermi, fermiM, G4FermiMomentum::Init(), G4INCL::Math::max(), momentum, myA, myL, myZ, nucleondistance, places, testSums, theDensity, theFermi, and theNucleons.

Referenced by G4MuMinusCapturePrecompound::ApplyYourself(), G4LowEIonFragmentation::ApplyYourself(), and G4BinaryLightIonReaction::Interact().

◆ operator!=()

G4bool G4Fancy3DNucleus::operator!= ( const G4Fancy3DNucleus & right ) const

private

◆ operator=()

const G4Fancy3DNucleus & G4Fancy3DNucleus::operator= ( const G4Fancy3DNucleus & right )

private

◆ operator==()

G4bool G4Fancy3DNucleus::operator== ( const G4Fancy3DNucleus & right ) const

private

◆ ReduceSum()

G4bool G4Fancy3DNucleus::ReduceSum ( )

private

Definition at line 524 of file G4Fancy3DNucleus.cc.

{
    G4ThreeVector sum;
    G4double PFermi=fermiM[myA-1];
 
    for (G4int i=0; i < myA-1 ; i++ )
         { sum+=momentum[i]; }
 
// check if have to do anything at all..
    if ( sum.mag() <= PFermi )
    {
        momentum[myA-1]=-sum;
        return true;
    }
 
// find all possible changes in momentum, changing only the component parallel to sum
    G4ThreeVector testDir=sum.unit();
    testSums.clear();
    testSums.resize(myA-1);     // Allocate block for filling below
 
    G4ThreeVector delta;
    for (G4int aNucleon=0; aNucleon < myA-1; aNucleon++) {
      delta = 2.*((momentum[aNucleon]*testDir)*testDir);
 
      testSums[aNucleon].Fill(delta, delta.mag(), aNucleon);
    }
 
    std::sort(testSums.begin(), testSums.end());
 
//    reduce Momentum Sum until the next would be allowed.
    G4int index=testSums.size();
    while ( (sum-testSums[--index].Vector).mag()>PFermi && index>0)  /* Loop checking, 30-Oct-2015, G.Folger */
    {
      // Only take one which improve, ie. don't change sign and overshoot...
      if ( sum.mag() > (sum-testSums[index].Vector).mag() ) {
        momentum[testSums[index].Index]-=testSums[index].Vector;
        sum-=testSums[index].Vector;
      }
    }
 
    if ( (sum-testSums[index].Vector).mag() <= PFermi )
    {
        G4int best=-1;
        G4double pBest=2*PFermi; // anything larger than PFermi
        for ( G4int aNucleon=0; aNucleon<=index; aNucleon++)
        {
            // find the momentum closest to choosen momentum for last Nucleon.
            G4double pTry=(testSums[aNucleon].Vector-sum).mag();
            if ( pTry < PFermi
             &&  std::abs(momentum[myA-1].mag() - pTry ) < pBest )
                {
               pBest=std::abs(momentum[myA-1].mag() - pTry );
               best=aNucleon;
            }
        }
        if ( best < 0 )  
        {
          G4String text = "G4Fancy3DNucleus.cc: Logic error in ReduceSum()";
              throw G4HadronicException(__FILE__, __LINE__, text);
        }
        momentum[testSums[best].Index]-=testSums[best].Vector;
        momentum[myA-1]=testSums[best].Vector-sum;
 
        return true;
 
    }
 
    // try to compensate momentum using another Nucleon....
    G4int swapit=-1;
    while (swapit<myA-1)  /* Loop checking, 30-Oct-2015, G.Folger */
    {
      if ( fermiM[++swapit] > PFermi ) break;
    }
    if (swapit == myA-1 ) return false;
 
    // Now we have a nucleon with a bigger Fermi Momentum.
    // Exchange with last nucleon.. and iterate.
    std::swap(theNucleons[swapit], theNucleons[myA-1]);
    std::swap(momentum[swapit], momentum[myA-1]);
    std::swap(fermiM[swapit], fermiM[myA-1]);
    return ReduceSum();
}

References fermiM, CLHEP::Hep3Vector::mag(), momentum, myA, ReduceSum(), testSums, theNucleons, and CLHEP::Hep3Vector::unit().

Referenced by ChooseFermiMomenta(), and ReduceSum().

◆ RefetchImpactXandY()

std::pair< G4double, G4double > G4V3DNucleus::RefetchImpactXandY ( )

inlineinherited

Definition at line 77 of file G4V3DNucleus.hh.

77{return theImpactParameter;}

References G4V3DNucleus::theImpactParameter.

◆ SortNucleonsDecZ()

void G4Fancy3DNucleus::SortNucleonsDecZ ( )

virtual

Implements G4V3DNucleus.

Definition at line 164 of file G4Fancy3DNucleus.cc.

{
  if (theNucleons.size() < 2 ) return;      // Avoid unnecessary work
  SortNucleonsIncZ();
 
  std::reverse(theNucleons.begin(), theNucleons.end());
}

References SortNucleonsIncZ(), and theNucleons.

◆ SortNucleonsIncZ()

void G4Fancy3DNucleus::SortNucleonsIncZ ( )

virtual

Implements G4V3DNucleus.

Definition at line 156 of file G4Fancy3DNucleus.cc.

{
  if (theNucleons.size() < 2 ) return;      // Avoid unnecesary work
 
  std::sort(theNucleons.begin(), theNucleons.end(),
        G4Fancy3DNucleusHelperForSortInZ); 
}

References G4Fancy3DNucleusHelperForSortInZ(), and theNucleons.

Referenced by SortNucleonsDecZ().

◆ StartLoop()

G4bool G4Fancy3DNucleus::StartLoop ( )

virtual

Implements G4V3DNucleus.

Definition at line 131 of file G4Fancy3DNucleus.cc.

{
    currentNucleon=0;
    return (theNucleons.size()>0);
}   

References currentNucleon, and theNucleons.

Referenced by G4LowEIonFragmentation::ApplyYourself(), G4BinaryLightIonReaction::GetProjectileExcitation(), and G4BinaryLightIonReaction::Interact().

Field Documentation

◆ currentNucleon

G4int G4Fancy3DNucleus::currentNucleon

private

Definition at line 110 of file G4Fancy3DNucleus.hh.

Referenced by GetNextNucleon(), Init(), and StartLoop().

◆ excitationEnergy

G4double G4Fancy3DNucleus::excitationEnergy

private

Definition at line 114 of file G4Fancy3DNucleus.hh.

Referenced by AddExcitationEnergy(), GetExcitationEnergy(), and Init().

◆ fermiM

std::vector<G4double> G4Fancy3DNucleus::fermiM

private

Definition at line 118 of file G4Fancy3DNucleus.hh.

Referenced by ChooseFermiMomenta(), Init(), and ReduceSum().

◆ momentum

std::vector<G4ThreeVector> G4Fancy3DNucleus::momentum

private

Definition at line 117 of file G4Fancy3DNucleus.hh.

Referenced by ChooseFermiMomenta(), Init(), and ReduceSum().

◆ myA

G4int G4Fancy3DNucleus::myA

private

Definition at line 105 of file G4Fancy3DNucleus.hh.

Referenced by BindingEnergy(), CenterNucleons(), ChooseFermiMomenta(), ChooseNucleons(), ChoosePositions(), CoulombBarrier(), DoLorentzBoost(), DoLorentzContraction(), DoTranslation(), GetMass(), GetMassNumber(), GetNextNucleon(), GetOuterRadius(), Init(), and ReduceSum().

◆ myL

G4int G4Fancy3DNucleus::myL

private

Definition at line 107 of file G4Fancy3DNucleus.hh.

Referenced by ChooseNucleons(), GetMass(), GetNumberOfLambdas(), and Init().

◆ myZ

G4int G4Fancy3DNucleus::myZ

private

Definition at line 106 of file G4Fancy3DNucleus.hh.

Referenced by BindingEnergy(), ChooseFermiMomenta(), ChooseNucleons(), CoulombBarrier(), GetCharge(), GetMass(), and Init().

◆ nucleondistance

G4double G4Fancy3DNucleus::nucleondistance

private

Definition at line 113 of file G4Fancy3DNucleus.hh.

Referenced by ChoosePositions(), GetOuterRadius(), and Init().

◆ places

std::vector<G4ThreeVector> G4Fancy3DNucleus::places

private

Definition at line 116 of file G4Fancy3DNucleus.hh.

Referenced by ChoosePositions(), and Init().

◆ testSums

std::vector<G4Fancy3DNucleusHelper> G4Fancy3DNucleus::testSums

private

Definition at line 119 of file G4Fancy3DNucleus.hh.

Referenced by Init(), and ReduceSum().

◆ theDensity

G4VNuclearDensity* G4Fancy3DNucleus::theDensity

private

Definition at line 111 of file G4Fancy3DNucleus.hh.

Referenced by ChooseFermiMomenta(), ChoosePositions(), GetNuclearDensity(), GetNuclearRadius(), Init(), and ~G4Fancy3DNucleus().

◆ theFermi

G4FermiMomentum G4Fancy3DNucleus::theFermi

private

Definition at line 112 of file G4Fancy3DNucleus.hh.

Referenced by ChooseFermiMomenta(), ChoosePositions(), and Init().

◆ theImpactParameter

std::pair<G4double, G4double> G4V3DNucleus::theImpactParameter

privateinherited

Definition at line 81 of file G4V3DNucleus.hh.

Referenced by G4V3DNucleus::ChooseImpactXandY(), and G4V3DNucleus::RefetchImpactXandY().

◆ theNucleons

std::vector<G4Nucleon> G4Fancy3DNucleus::theNucleons

private

Definition at line 108 of file G4Fancy3DNucleus.hh.

Referenced by CenterNucleons(), ChooseFermiMomenta(), ChooseNucleons(), ChoosePositions(), DoLorentzBoost(), DoLorentzContraction(), DoTranslation(), GetNextNucleon(), GetNucleons(), GetOuterRadius(), Init(), ReduceSum(), SortNucleonsDecZ(), SortNucleonsIncZ(), and StartLoop().

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

source/processes/hadronic/util/include/G4Fancy3DNucleus.hh
source/processes/hadronic/util/src/G4Fancy3DNucleus.cc

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ G4Fancy3DNucleus() [1/2]

◆ ~G4Fancy3DNucleus()

◆ G4Fancy3DNucleus() [2/2]

Member Function Documentation

◆ AddExcitationEnergy()

◆ BindingEnergy()

◆ CenterNucleons()

◆ ChooseFermiMomenta()

◆ ChooseImpactXandY()

◆ ChooseNucleons()

◆ ChoosePositions()

◆ CoulombBarrier()

◆ DoLorentzBoost() [1/2]

◆ DoLorentzBoost() [2/2]

◆ DoLorentzContraction() [1/2]

◆ DoLorentzContraction() [2/2]

◆ DoTranslation()

◆ GetCharge()

◆ GetExcitationEnergy()

◆ GetMass()

◆ GetMassNumber()

◆ GetNextNucleon()

◆ GetNuclearDensity()

◆ GetNuclearRadius() [1/2]

◆ GetNuclearRadius() [2/2]

◆ GetNucleons()

◆ GetNumberOfLambdas()

◆ GetOuterRadius()

◆ Init()

◆ operator!=()

◆ operator=()

◆ operator==()

◆ ReduceSum()

◆ RefetchImpactXandY()

◆ SortNucleonsDecZ()

◆ SortNucleonsIncZ()

◆ StartLoop()

Field Documentation

◆ currentNucleon

◆ excitationEnergy

◆ fermiM

◆ momentum

◆ myA

◆ myL

◆ myZ

◆ nucleondistance

◆ places

◆ testSums

◆ theDensity

◆ theFermi

◆ theImpactParameter

◆ theNucleons