Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improvements to the Radial Range that Subhalo Density Profiles are Tabulated when Fitting TNFW Profile #483

Merged
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 23 additions & 11 deletions source/nodes.property_extractor.tidally_truncated_NFW_fit.F90
Original file line number Diff line number Diff line change
Expand Up @@ -165,18 +165,20 @@ function tidallyTruncatedNFWFitExtract(self,node,time,instance)
type (multiCounter ), intent(inout), optional :: instance
class (nodeComponentDarkMatterProfile ), pointer :: darkMatterProfile
class (nodeComponentSatellite ), pointer :: satellite
double precision , allocatable , dimension(:) :: radii , fractionDensity
double precision , allocatable , dimension(:) :: radii , fractionDensity
type (multiDMinimizer ), allocatable :: minimizer_
double precision , dimension(1) :: locationMinimum
double precision , parameter :: fractionRadiusOuter=0.95d0, fractionRadiusScale=0.1d0, &
& fractionMaximum =0.10d0, fractionStep =0.1d0
integer , parameter :: countRadiiPerDecade=10
integer :: countRadii , i , &
& iteration
double precision , parameter :: fractionRadiusScale =0.1d0, fractionMaximum =0.1d0, &
& radiusMaximumFractionDensityVirialMinimum=0.1d0, fractionStep =0.1d0, &
& radiusMaximumScaleVirialMaximum =1.0d1
integer , parameter :: radiusMaximumCountRadiiPerDecade =10 , countRadiiPerDecade=10
integer :: countRadii , i , &
& iteration , radiusMaximumCountRadii
logical :: converged
double precision :: radiusOuter , massTotal , &
& radiusMinimum , radiusMaximum , &
& radiusScale , radiusVirial
double precision :: radiusOuter , massTotal , &
& radiusMinimum , radiusMaximum , &
& radiusScale , radiusVirial , &
& radiusMaximumFractionDensityVirial , factorStepRadius
!$GLC attributes unused :: instance

allocate(tidallyTruncatedNFWFitExtract(3))
Expand All @@ -190,8 +192,18 @@ function tidallyTruncatedNFWFitExtract(self,node,time,instance)
radiusScale = darkMatterProfile %scale ( )
radiusVirial = self %darkMatterHaloScale_%radiusVirial (node )
! Choose radii for fitting.
radiusMaximum= fractionRadiusOuter*radiusOuter
radiusMinimum=min(fractionRadiusScale*radiusScale,fractionMaximum*radiusMaximum)
radiusMaximum=radiusVirial
if (radiusOuter > radiusVirial) then
radiusMaximumCountRadii=int(log10(radiusMaximumScaleVirialMaximum)*dble(radiusMaximumCountRadiiPerDecade)+1.0d0)
factorStepRadius = log10(radiusMaximumScaleVirialMaximum)/dble(radiusMaximumCountRadii )
do i=1,radiusMaximumCountRadii
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This was previously:

         do n=1,radiusMaxiumCountRadii
             radiusMaximum                     =10**(log10(radiusVirial)+step*dble(i))

So the loop variable was n but on the following line it used i to increment radiusMaximum. Most likely i=0 initially, so this would not have actually increased radiusMaximum, leaving it always equal to radiusVirial. I've corrected this (just using i as the loop variable), but it would be good if you can check if this still gives good results for the fits.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I recompiled and re-ran the parameter with the bug fixed. The fits look good.

radiusMaximum =10.0d0**(log10(radiusVirial)+factorStepRadius*dble(i))
radiusMaximumFractionDensityVirial=+self%darkMatterProfileDMO_%density(node,radiusMaximum) &
& /self%darkMatterProfileDMO_%density(node,radiusVirial )
if (radiusMaximumFractionDensityVirial < radiusMaximumFractionDensityVirialMinimum) exit
end do
end if
radiusMinimum=min(fractionRadiusScale*radiusScale,fractionMaximum*radiusMaximum,fractionMaximum*radiusOuter)
countRadii =int(log10(radiusMaximum/radiusMinimum)*dble(countRadiiPerDecade)+1.0d0)
radii =Make_Range(radiusMinimum,radiusMaximum,countRadii,rangeTypeLogarithmic)
! Tabulate the density ratio relative to an NFW profile.
Expand Down
Loading