Skip to content

Commit

Permalink
fix: Pass "node" object to critical overdensity functions to allow en…
Browse files Browse the repository at this point in the history 
…vironment to be accessed

The "environmental" critical overdensity class must be able to access the halo environment, which requires it to be passed a treeNode object, "node", which has a pointer to the host tree (which stores the environment). This fix ensures that happens for several calling functions, and adds an error message from the "environmental" critical overdensity class is "node" is not present.
  • Loading branch information
@abensonca
abensonca committed Aug 2, 2021
1 parent 2a20104 commit 1218d96
Show file tree
Hide file tree
Showing 5 changed files with 41 additions and 32 deletions.
11 changes: 6 additions & 5 deletions source/dark_matter_halos.mass_accretion_history.Wechsler2002.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -171,7 +171,7 @@ double precision function wechsler2002Time(self,node,mass)
select case (self%formationRedshiftCompute)
case (.true.)
! Compute the expansion factor at formation.
mergerTreeFormationExpansionFactor=self%expansionFactorAtFormation(baseBasic%mass())
mergerTreeFormationExpansionFactor=self%expansionFactorAtFormation(baseBasic%mass(),node)
case (.false.)
! Use the specified formation redshift.
mergerTreeFormationExpansionFactor=self%cosmologyFunctions_%expansionFactorFromRedshift(self%formationRedshift)
Expand Down Expand Up @@ -203,7 +203,7 @@ double precision function wechsler2002MassAccretionRate(self,node,time)
select case (self%formationRedshiftCompute)
case (.true.)
! Compute the expansion factor at formation.
mergerTreeFormationExpansionFactor=self%expansionFactorAtFormation(baseBasic%mass())
mergerTreeFormationExpansionFactor=self%expansionFactorAtFormation(baseBasic%mass(),node)
case (.false.)
! Use the specified formation redshift.
mergerTreeFormationExpansionFactor=self%cosmologyFunctions_%expansionFactorFromRedshift(self%formationRedshift)
Expand Down Expand Up @@ -231,23 +231,24 @@ double precision function wechsler2002MassAccretionRate(self,node,time)
return
end function wechsler2002MassAccretionRate

double precision function wechsler2002ExpansionFactorAtFormation(self,haloMass)
double precision function wechsler2002ExpansionFactorAtFormation(self,haloMass,node)
!!{
Computes the expansion factor at formation using the simple model of \cite{bullock_profiles_2001}.
!!}
implicit none
class (darkMatterHaloMassAccretionHistoryWechsler2002), intent(inout) :: self
type (treeNode ), intent(inout) :: node
double precision , intent(in ) :: haloMass
double precision , parameter :: haloMassFraction =0.015d0 ! Wechsler et al. (2002; Astrophysical Journal, 568:52-70).
double precision :: formationTime , haloMassCharacteristic, &
& sigmaCharacteristic

! Compute the characteristic mass at formation time.
haloMassCharacteristic=haloMassFraction*haloMass
! Compute the corresponding rms fluctuation in the density field (i.e. sigma(M)).
! Compute the corresponding rms fluctuation in the density field (i.e. σ(M)).
sigmaCharacteristic=self%cosmologicalMassVariance_%rootVariance(haloMassCharacteristic,self%cosmologyFunctions_%cosmicTime(1.0d0))
! Get the time at which this equals the critical overdensity for collapse.
formationTime=self%criticalOverdensity_%timeOfCollapse(criticalOverdensity=sigmaCharacteristic,mass=haloMass)
formationTime=self%criticalOverdensity_%timeOfCollapse(criticalOverdensity=sigmaCharacteristic,mass=haloMass,node=node)
! Get the corresponding expansion factor.
wechsler2002ExpansionFactorAtFormation=self%cosmologyFunctions_%expansionFactor(formationTime)
return
Expand Down
38 changes: 19 additions & 19 deletions source/dark_matter_profiles.structure.concentration.Diemer-Joyce2019.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -311,31 +311,31 @@ double precision function diemerJoyce2019ConcentrationMean(self,node)
& .or. &
& basic%time() /= self%timePrevious &
& ) then
peakHeight =+self%criticalOverdensity_ %value (time=basic%time(),mass=basic%mass()) &
& /self%cosmologicalMassVariance_%rootVariance (time=basic%time(),mass=basic%mass())
massHalo =+self%kappa**3 &
peakHeight =+self%criticalOverdensity_ %value (time=basic%time(),mass=basic%mass(),node=node) &
& /self%cosmologicalMassVariance_%rootVariance (time=basic%time(),mass=basic%mass() )
massHalo =+self%kappa**3 &
& *basic%mass()
diemerJoyce2019PowerSpectrumSlope=-6.0d0 &
& *self%cosmologicalMassVariance_%rootVarianceLogarithmicGradient(time=basic%time(),mass=massHalo ) &
diemerJoyce2019PowerSpectrumSlope=-6.0d0 &
& *self%cosmologicalMassVariance_%rootVarianceLogarithmicGradient(time=basic%time(),mass=massHalo ) &
& -3.0d0
alphaEffective =+self%linearGrowth_%logarithmicDerivativeExpansionFactor(time=basic%time())
A =+self%a0 &
& *( &
& +1.0d0 &
& +self%a1*(diemerJoyce2019PowerSpectrumSlope+3.0d0) &
A =+self%a0 &
& *( &
& +1.0d0 &
& +self%a1*(diemerJoyce2019PowerSpectrumSlope+3.0d0) &
& )
B =+self%b0 &
& *( &
& +1.0d0 &
& +self%b1*(diemerJoyce2019PowerSpectrumSlope+3.0d0) &
B =+self%b0 &
& *( &
& +1.0d0 &
& +self%b1*(diemerJoyce2019PowerSpectrumSlope+3.0d0) &
& )
C =+1.0d0 &
C =+1.0d0 &
& -self%cAlpha*(1.0d0-alphaEffective)
diemerJoyce2019GRoot =+A &
& /peakHeight &
& *( &
& +1.0d0 &
& +peakHeight**2/B &
diemerJoyce2019GRoot =+A &
& /peakHeight &
& *( &
& +1.0d0 &
& +peakHeight**2/B &
& )
! Initial guess of the concentration.
if (self%concentrationMeanPrevious < 0.0d0) self%concentrationMeanPrevious=5.0d0
Expand Down
11 changes: 6 additions & 5 deletions source/dark_matter_profiles.structure.scale.concentration.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -238,11 +238,12 @@ double precision function concentrationRadius(self,node)
& ) &
& ) then
! Create a node and set the mass and time.
concentrationState_(concentrationStateCount)%self => self
concentrationState_(concentrationStateCount)%node => node
concentrationState_(concentrationStateCount)%nodeWork => treeNode ( )
concentrationState_(concentrationStateCount)%basic => concentrationState_(concentrationStateCount)%nodeWork%basic (autoCreate=.true.)
concentrationState_(concentrationStateCount)%darkMatterProfile => concentrationState_(concentrationStateCount)%nodeWork%darkMatterProfile(autoCreate=.true.)
concentrationState_(concentrationStateCount)%self => self
concentrationState_(concentrationStateCount)%node => node
concentrationState_(concentrationStateCount)%nodeWork => treeNode ( )
concentrationState_(concentrationStateCount)%nodeWork %hostTree => node%hostTree
concentrationState_(concentrationStateCount)%basic => concentrationState_(concentrationStateCount)%nodeWork%basic (autoCreate=.true.)
concentrationState_(concentrationStateCount)%darkMatterProfile => concentrationState_(concentrationStateCount)%nodeWork%darkMatterProfile(autoCreate=.true.)
call concentrationState_(concentrationStateCount)%basic%timeSet (basic%time())
call concentrationState_(concentrationStateCount)%basic%timeLastIsolatedSet(basic%time())
! The finder is initialized each time as it is allocated on the stack - this allows this function to be called recursively.
Expand Down
7 changes: 4 additions & 3 deletions source/nodes.property_extractor.peak_height.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -134,9 +134,10 @@ function peakHeightExtract(self,node,time,instance)
!$GLC attributes unused :: time, instance

basic => node%basic()
criticalOverdensityLastIsolated=self%criticalOverdensity_ %value (mass=basic%mass(),time=basic%timeLastIsolated())
densityFieldRootVariance =self%cosmologicalMassVariance_%rootVariance(mass=basic%mass(),time=basic%timeLastIsolated())
peakHeightNu =criticalOverdensityLastIsolated/densityFieldRootVariance
criticalOverdensityLastIsolated= self%criticalOverdensity_ %value (mass=basic%mass(),time=basic%timeLastIsolated(),node=node)
densityFieldRootVariance = self%cosmologicalMassVariance_%rootVariance(mass=basic%mass(),time=basic%timeLastIsolated() )
peakHeightNu =+criticalOverdensityLastIsolated &
& /densityFieldRootVariance
allocate(peakHeightExtract(3))
peakHeightExtract=[ &
& criticalOverdensityLastIsolated, &
Expand Down
6 changes: 6 additions & 0 deletions source/structure_formation.critical_overdensity.environmental.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -135,6 +135,7 @@ double precision function environmentalValue(self,time,expansionFactor,collapsin
!!{
Return the critical overdensity for collapse at the given time and mass.
!!}
use :: Galacticus_Error, only : Galacticus_Error_Report
use :: Galacticus_Nodes, only : nodeComponentBasic, treeNode
implicit none
class (criticalOverdensityEnvironmental), intent(inout) :: self
Expand All @@ -144,6 +145,7 @@ double precision function environmentalValue(self,time,expansionFactor,collapsin
type (treeNode ), intent(inout), optional :: node
class (nodeComponentBasic ), pointer :: basic

if (.not.present(node)) call Galacticus_Error_Report('"node" must be provided to give access to environment'//{introspection:location})
! Get the critical overdensity at zero environmental overdensity and scale by some power of the linear growth factor.
environmentalValue = + self%criticalOverdensity_ %value (time,expansionFactor,collapsing,mass,node )
basic => node%hostTree %baseNode%basic ( )
Expand All @@ -163,6 +165,7 @@ double precision function environmentalGradientTime(self,time,expansionFactor,co
Return the gradient with respect to time of critical overdensity at the given time and mass.
!!}
use :: Galacticus_Nodes, only : nodeComponentBasic, treeNode
use :: Galacticus_Error, only : Galacticus_Error_Report
implicit none
class (criticalOverdensityEnvironmental), intent(inout) :: self
double precision , intent(in ), optional :: time , expansionFactor
Expand All @@ -174,6 +177,7 @@ double precision function environmentalGradientTime(self,time,expansionFactor,co
<optionalArgument name="expansionFactor" defaultsTo="self%cosmologyFunctions_%expansionFactor(time)" />
!!]

if (.not.present(node)) call Galacticus_Error_Report('"node" must be provided to give access to environment'//{introspection:location})
basic => node%hostTree%baseNode%basic()
if (basic%mass() < self%massEnvironment) then
environmentalGradientTime=+ self%criticalOverdensity_%gradientTime (time,expansionFactor ,collapsing,mass,node ) &
Expand Down Expand Up @@ -203,6 +207,7 @@ double precision function environmentalGradientMass(self,time,expansionFactor,co
Return the gradient with respect to mass of critical overdensity at the given time and mass.
!!}
use :: Galacticus_Nodes, only : nodeComponentBasic, treeNode
use :: Galacticus_Error, only : Galacticus_Error_Report
implicit none
class (criticalOverdensityEnvironmental), intent(inout) :: self
double precision , intent(in ), optional :: time , expansionFactor
Expand All @@ -211,6 +216,7 @@ double precision function environmentalGradientMass(self,time,expansionFactor,co
type (treeNode ), intent(inout), optional :: node
class (nodeComponentBasic ), pointer :: basic

if (.not.present(node)) call Galacticus_Error_Report('"node" must be provided to give access to environment'//{introspection:location})
basic => node%hostTree%baseNode%basic()
if (basic%mass() < self%massEnvironment) then
environmentalGradientMass= +self%criticalOverdensity_%gradientMass (time,expansionFactor,collapsing,mass,node ) &
Expand Down

0 comments on commit 1218d96

Please sign in to comment.