Skip to content

Commit

Permalink
Merge pull request #656 from galacticusorg/fixRamPressureSegFault
Browse files Browse the repository at this point in the history
Avoid segfaults when computing ram pressure forces
  • Loading branch information
abensonca authored Jul 27, 2024
2 parents b90646e + ea070cd commit d38fdea
Show file tree
Hide file tree
Showing 5 changed files with 326 additions and 10 deletions.
16 changes: 16 additions & 0 deletions parameters/reference/cosmologyMillennium.xml
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
<?xml version="1.0" encoding="UTF-8"?>
<!-- Millennium simulation cosmological model -->
<parameters>
<formatVersion>2</formatVersion>

<!-- Cosmological parameters and options -->
<cosmologyFunctions value="matterLambda"/>
<cosmologyParameters value="simple" >
<HubbleConstant value="73.00000"/>
<OmegaMatter value=" 0.25000"/>
<OmegaDarkEnergy value=" 0.75000"/>
<OmegaBaryon value=" 0.04500"/>
<temperatureCMB value=" 2.72548"/>
</cosmologyParameters>

</parameters>
28 changes: 28 additions & 0 deletions parameters/reference/darkMatterHalosStructureNBody.xml
Original file line number Diff line number Diff line change
@@ -0,0 +1,28 @@
<?xml version="1.0" encoding="UTF-8"?>
<!-- Reference dark matter halo structure model -->
<parameters>
<formatVersion>2</formatVersion>

<!-- Dark matter profile scale radii model -->
<darkMatterProfileScaleRadius value="concentrationLimiter">
<!-- Limit scale radii to keep concentrations within a reasonable range. -->
<concentrationMinimum value=" 3.0"/>
<concentrationMaximum value="100.0"/>
<!-- Use the Diemer & Joyce (2019; ApJ; 871; 168; http://adsabs.harvard.edu/abs/2019ApJ...871..168D) model for concentrations. -->
<darkMatterProfileScaleRadius value="concentration">
<correctForConcentrationDefinition value="true"/>
<darkMatterProfileConcentration value="diemerJoyce2019">
</darkMatterProfileConcentration>
</darkMatterProfileScaleRadius>
</darkMatterProfileScaleRadius>

<!-- Dark matter halo spin -->
<haloSpinDistribution value="bett2007">
<!-- For leaf nodes in the tree we fall back to drawing spins from the distribution function given by -->
<!-- Benson (2017; MNRAS; 471; 2871; http://adsabs.harvard.edu/abs/2017MNRAS.471.2871B). -->
<!-- Best fit paramter values are taken from that paper. -->
<alpha value="1.7091800"/>
<lambda0 value="0.0420190"/>
</haloSpinDistribution>

</parameters>
241 changes: 241 additions & 0 deletions parameters/reference/evolutionGalaxyFormationNBody.xml
Original file line number Diff line number Diff line change
@@ -0,0 +1,241 @@
<?xml version="1.0" encoding="UTF-8"?>
<!-- Reference subhalo orbits model -->
<parameters>
<version>0.9.4</version>
<formatVersion>2</formatVersion>

<!-- Component selection -->
<componentBasic value="standard"/>
<componentBlackHole value="standard" >
<massSeed value="100.0000"/>
<heatsHotHalo value="true" />
<efficiencyWind value=" 0.0024"/>
<efficiencyWindScalesWithEfficiencyRadiative value="true" />
<bondiHoyleAccretionEnhancementHotHalo value=" 6.0000"/>
<bondiHoyleAccretionEnhancementSpheroid value=" 5.0000"/>
<bondiHoyleAccretionTemperatureSpheroid value="100.0000"/>
<bondiHoyleAccretionHotModeOnly value="true" />
</componentBlackHole>
<componentDarkMatterProfile value="scale" />
<componentDisk value="standard" >
<massDistributionDisk value="exponentialDisk">
<dimensionless value="true"/>
</massDistributionDisk>
</componentDisk>
<componentHotHalo value="standard" >
<fractionLossAngularMomentum value="0.3" />
<starveSatellites value="false"/>
<efficiencyStrippingOutflow value="0.1" />
<trackStrippedGas value="true" />
</componentHotHalo>
<componentPosition value="cartesian"/>
<componentSatellite value="preset" />
<componentSpheroid value="standard" >
<ratioAngularMomentumScaleRadius value="0.5" />
<efficiencyEnergeticOutflow value="1.0e-2"/>
<massDistributionSpheroid value="hernquist">
<dimensionless value="true"/>
</massDistributionSpheroid>
</componentSpheroid>
<componentSpin value="vector" />

<!-- Intergalactic background radiation -->
<radiationFieldIntergalacticBackground value="summation">
<radiationField value="cosmicMicrowaveBackground"/>
<radiationField value="intergalacticBackgroundFile">
<fileName value="%DATASTATICPATH%/radiation/Cosmic_Background_Radiation_FG20.hdf5"/>
</radiationField>
</radiationFieldIntergalacticBackground>

<!-- Halo accretion options -->
<accretionHalo value="naozBarkana2007"/>
<intergalacticMediumFilteringMass value="gnedin2000"/>
<intergalacticMediumState value="metallicityPolynomial">
<coefficients value="-1.3 -1.9"/>
<intergalacticMediumState value="instantReionization">
<electronScatteringOpticalDepth value="0.0633" />
<reionizationTemperature value="2.0e4" />
<presentDayTemperature value="1.0e3" />
<intergalacticMediumState value="recFast"/>
</intergalacticMediumState>
</intergalacticMediumState>

<!-- Hot halo gas cooling model options -->
<hotHaloMassDistribution value="betaProfile" />
<hotHaloTemperatureProfile value="virial" />
<hotHaloMassDistributionCoreRadius value="virialFraction" >
<coreRadiusOverVirialRadius value="0.3"/>
</hotHaloMassDistributionCoreRadius>
<coolingSpecificAngularMomentum value="constantRotation">
<sourceAngularMomentumSpecificMean value="hotGas"/>
<sourceNormalizationRotation value="hotGas"/>
</coolingSpecificAngularMomentum>
<hotHaloOutflowReincorporation value="haloDynamicalTime" >
<multiplier value="5.0"/>
</hotHaloOutflowReincorporation>

<coolingFunction value="atomicCIECloudy"/>
<coolingRadius value="simple"/>
<coolingRate value="multiplier">
<multiplier value="0.5" />
<coolingRate value="whiteFrenk1991" >
<velocityCutOff value="10000"/>
</coolingRate>
</coolingRate>
<coolingTime value="simple">
<degreesOfFreedom value="3.0"/>
</coolingTime>
<coolingTimeAvailable value="whiteFrenk1991">
<ageFactor value="0.0"/>
</coolingTimeAvailable>

<!-- Hot halo ram pressure stripping options -->
<hotHaloRamPressureStripping value="font2008">
<solverFailureIsFatal value="false"/>
</hotHaloRamPressureStripping>
<hotHaloRamPressureForce value="relativePosition" />
<hotHaloRamPressureTimescale value="ramPressureAcceleration"/>

<!-- Galactic structure solver options -->
<galacticStructureSolver value="equilibrium"/>
<darkMatterProfile value="adiabaticGnedin2004">
<A value="0.73"/>
<omega value="0.70"/>
</darkMatterProfile>

<!-- Star formation rate options -->
<starFormationRateDisks value="intgrtdSurfaceDensity"/>
<starFormationRateSurfaceDensityDisks value="blitz2006" >
<assumeMonotonicSurfaceDensity value="true"/>
<assumeExponentialDisk value="true"/>
<useTabulation value="true"/>
</starFormationRateSurfaceDensityDisks>
<starFormationRateSpheroids value="timescale">
<starFormationTimescale value="dynamicalTime">
<efficiency value="0.040"/>
<exponentVelocity value="2.000"/>
<timescaleMinimum value="0.001"/>
</starFormationTimescale>
</starFormationRateSpheroids>

<!-- Stellar populations options -->
<stellarPopulationProperties value="instantaneous"/>
<stellarPopulationSpectra value="FSPS" />
<stellarPopulationSelector value="fixed" />

<initialMassFunction value="chabrier2001"/>
<stellarPopulation value="standard">
<recycledFraction value="0.460"/>
<metalYield value="0.035"/>
</stellarPopulation>

<!-- AGN feedback options -->
<hotHaloExcessHeatDrivesOutflow value="true"/>

<!-- Accretion disk properties -->
<accretionDisks value="switched">
<accretionRateThinDiskMaximum value="0.30" />
<accretionRateThinDiskMinimum value="0.01" />
<scaleADAFRadiativeEfficiency value="true" />
<accretionDisksShakuraSunyaev value="shakuraSunyaev"/>
<accretionDisksADAF value="ADAF" >
<efficiencyRadiationType value="thinDisk"/>
<adiabaticIndex value="1.444" />
<energyOption value="pureADAF"/>
<efficiencyRadiation value="0.01" />
<viscosityOption value="fit" />
</accretionDisksADAF>
</accretionDisks>

<!-- Black hole options -->
<blackHoleBinaryMerger value="rezzolla2008"/>

<!-- Galaxy merger options -->
<mergerMassMovements value="simple">
<destinationGasMinorMerger value="spheroid"/>
<massRatioMajorMerger value="0.25" />
</mergerMassMovements>
<mergerRemnantSize value="cole2000">
<energyOrbital value="1.0"/>
</mergerRemnantSize>

<!-- Node evolution and physics -->
<nodeOperator value="multi">
<!-- Cosmological epoch -->
<nodeOperator value="cosmicTime"/>
<!-- DMO evolution -->
<nodeOperator value="DMOInterpolate"/>
<!-- Halo concentrations -->
<nodeOperator value="darkMatterProfileScaleInterpolate"/>
<!-- Interpolate halo spins -->
<nodeOperator value="haloAngularMomentumRandom">
<factorReset value="2.0"/>
</nodeOperator>
<nodeOperator value="haloAngularMomentumInterpolate"/>
<!-- Star formation -->
<nodeOperator value="starFormationDisks" />
<nodeOperator value="starFormationSpheroids"/>
<!--Stellar feedback outflows-->
<nodeOperator value="stellarFeedbackDisks">
<stellarFeedbackOutflows value="rateLimit">
<timescaleOutflowFractionalMinimum value="0.001"/>
<stellarFeedbackOutflows value="powerLaw">
<velocityCharacteristic value="250.0"/>
<exponent value=" 2.0"/>
</stellarFeedbackOutflows>
</stellarFeedbackOutflows>
</nodeOperator>
<nodeOperator value="stellarFeedbackSpheroids">
<stellarFeedbackOutflows value="rateLimit">
<timescaleOutflowFractionalMinimum value="0.001"/>
<stellarFeedbackOutflows value="powerLaw">
<velocityCharacteristic value="100.0"/>
<exponent value=" 2.0"/>
</stellarFeedbackOutflows>
</stellarFeedbackOutflows>
</nodeOperator>
<!-- Bar instability in galactic disks -->
<nodeOperator value="barInstability">
<galacticDynamicsBarInstability value="efstathiou1982">
<stabilityThresholdGaseous value="0.7"/>
<stabilityThresholdStellar value="1.1"/>
<timescaleMinimum value="1.0e-3"/>
</galacticDynamicsBarInstability>
</nodeOperator>
<nodeOperator value="nodeFormationTimeMassFraction">
<fractionMassFormation value="0.5"/>
</nodeOperator>
</nodeOperator>

<!-- Merger tree evolution -->
<mergerTreeEvolver value="standard">
<!-- Standard merger tree evolver with parameters chosen to (somewhat) optimize the evolution. -->
<timestepHostAbsolute value="1.00" />
<timestepHostRelative value="0.10" />
<fractionTimestepSatelliteMinimum value="0.75" />
<backtrackToSatellites value="true" />
<allTreesExistAtFinalTime value="false"/>
</mergerTreeEvolver>
<mergerTreeNodeEvolver value="standard">
<!-- Standard node evolve with parameters chosen to (somewhat) optimize the evolution. -->
<odeToleranceAbsolute value="0.01" />
<odeToleranceRelative value="0.01" />
<reuseODEStepSize value="false"/>
</mergerTreeNodeEvolver>
<mergerTreeEvolveTimestep value="multi">
<!-- Standard time-stepping rule -->
<mergerTreeEvolveTimestep value="simple">
<timeStepAbsolute value="1.000"/>
<timeStepRelative value="0.100"/>
</mergerTreeEvolveTimestep>
<!-- Limit timesteps based on star formation history times - ensures accurate recording of star formation history.
(Will have no effect if star formation history is not being recorded.) -->
<mergerTreeEvolveTimestep value="starFormationHistory"/>
<!-- Limit timesteps based on satellite evolution -->
<mergerTreeEvolveTimestep value="satellite">
<timeOffsetMaximumAbsolute value="0.010"/>
<timeOffsetMaximumRelative value="0.001"/>
</mergerTreeEvolveTimestep>
</mergerTreeEvolveTimestep>

</parameters>
27 changes: 27 additions & 0 deletions parameters/reference/powerSpectrumMillennium.xml
Original file line number Diff line number Diff line change
@@ -0,0 +1,27 @@
<?xml version="1.0" encoding="UTF-8"?>
<!-- Millennium simulation power spectrum model -->
<parameters>
<formatVersion>2</formatVersion>

<!-- Power spectrum options -->
<cosmologicalMassVariance value="filteredPower">
<sigma_8 value="0.9" />
<tolerance value="3.0e-4"/>
<toleranceTopHat value="3.0e-4"/>
<nonMonotonicIsFatal value="false" />
<monotonicInterpolation value="true" />
<powerSpectrumWindowFunction value="topHat"/>
</cosmologicalMassVariance>
<transferFunction value="CAMB" >
<redshift value="100.0"/>
<!-- Include cosmological model explicitly here so that baryons are included in the transfer function even for dark matter-only models. -->
<xi:include href="cosmologyMillennium.xml" xpointer="xpointer(parameters/*)" xmlns:xi="http://www.w3.org/2001/XInclude"/>
</transferFunction>
<powerSpectrumPrimordial value="powerLaw" >
<index value="1.0"/>
<wavenumberReference value="1.0"/>
<running value="0.0"/>
</powerSpectrumPrimordial>
<powerSpectrumPrimordialTransferred value="simple" />

</parameters>
24 changes: 14 additions & 10 deletions source/hot_halo.ram_pressure_force.relative_position.F90
Original file line number Diff line number Diff line change
Expand Up @@ -151,21 +151,25 @@ double precision function relativePositionForce(self,node) result(force)

! Find the host node. Seek the descendant of the node closest in time to our satellite node. This is necessary as satellites
! can evolve ahead of their hosts.
basic => node%basic ()
nodeHostPrevious => node%parent
nodeHostCurrent => node%parent
basicPrevious => nodeHostPrevious%basic()
basic => node %basic()
basicCurrent => nodeHostCurrent %basic()
basicPrevious => nodeHostPrevious%basic()
do while (associated(nodeHostCurrent))
basicCurrent => nodeHostCurrent%basic()
if (basicCurrent%time() > basic%time()) exit
basicPrevious => basicCurrent
nodeHostPrevious => nodeHostCurrent
nodeHostCurrent => nodeHostCurrent%parent
nodeHostPrevious => nodeHostCurrent
nodeHostCurrent => nodeHostCurrent%parent
basicPrevious => basicCurrent
if (associated(nodeHostCurrent)) &
& basicCurrent => nodeHostCurrent%basic ()
end do
if ( &
& abs(basicPrevious%time()-basic%time()) &
& < &
& abs(basicCurrent %time()-basic%time()) &
if ( &
& abs(basicPrevious%time()-basic%time()) &
& < &
& abs(basicCurrent %time()-basic%time()) &
& .or. &
& .not.associated(nodeHostCurrent) &
& ) then
nodeHost => nodeHostPrevious
basicHost => basicPrevious
Expand Down

0 comments on commit d38fdea

Please sign in to comment.