Skip to content

Commit

Permalink
Merge pull request #622 from galacticusorg/featMergedSubhaloProperties
Browse files Browse the repository at this point in the history 
Add functionality to record (and output) the orbital properties of merged subhalos
  • Loading branch information
@abensonca
abensonca authored Jun 8, 2024
2 parents 6e05898 + d86c915 commit fc98fac
Show file tree
Hide file tree
Showing 4 changed files with 403 additions and 29 deletions.
71 changes: 70 additions & 1 deletion perl/Galacticus/Build/SourceTree/Process/Enumeration.pm
View file
Original file line number Diff line number Diff line change
Expand Up @@ -305,14 +305,83 @@ sub Process_Enumerations {
&Galacticus::Build::SourceTree::SetVisibility($node->{'parent'},$decodeFunctionName,$visibility);
}
# Create description function.
if ( exists($node->{'directive'}->{'decodeFunction'}) && $node->{'directive'}->{'decodeFunction'} eq "yes" ) {
my $functionName = "enumeration".ucfirst($node->{'directive'}->{'name'})."Description";
my $interface;
$interface .= " interface ".$functionName."\n";
$interface .= " module procedure ".$functionName."Enumerator\n";
$interface .= " module procedure ".$functionName."ID\n";
$interface .= " end interface ".$functionName."\n\n";
my $descriptorFunctionWrapper;
$descriptorFunctionWrapper .= "\n";
$descriptorFunctionWrapper .= " ! Auto-generated enumeration function\n";
$descriptorFunctionWrapper .= " function ".$functionName."Enumerator(enumerationValue)\n";
$descriptorFunctionWrapper .= " !!{\n";
$descriptorFunctionWrapper .= " Return a description of a {\\normalfont \\ttfamily ".$node->{'directive'}->{'name'}."} enumeration member.\n";
$descriptorFunctionWrapper .= " !!}\n";
$descriptorFunctionWrapper .= " use ISO_Varying_String\n";
$descriptorFunctionWrapper .= " implicit none\n\n";
$descriptorFunctionWrapper .= " type (varying_string) :: ".$functionName."Enumerator\n";
$descriptorFunctionWrapper .= " type (enumeration".$node->{'directive'}->{'name'}."Type), intent(in ) :: enumerationValue\n\n";
$descriptorFunctionWrapper .= " ".$functionName."Enumerator=".$functionName."(enumerationValue%ID)\n";
$descriptorFunctionWrapper .= " return\n";
$descriptorFunctionWrapper .= " end function ".$functionName."Enumerator\n";
my $descriptorFunction;
$descriptorFunction .= " function ".$functionName."ID(enumerationValue) result(description)\n";
$descriptorFunction .= " !!{\n";
$descriptorFunction .= " Return a description of a {\\normalfont \\ttfamily ".$node->{'directive'}->{'name'}."} enumeration value.\n";
$descriptorFunction .= " !!}\n";
$descriptorFunction .= " use :: ISO_Varying_String, only : varying_string, assignment(=)\n";
$descriptorFunction .= " implicit none\n";
$descriptorFunction .= " type(varying_string) :: description\n";
$descriptorFunction .= " integer, intent(in) :: enumerationValue\n\n";
my $description = " select case (enumerationValue)\n";
my @entries = &List::ExtraUtils::as_array($node->{'directive'}->{'entry'});
my $i = $indexing-1;
foreach my $entry ( @entries ) {
$description .= " case (".++$i.")\n";
$description .= " description='".(exists($entry->{'description'}) ? $entry->{'description'} : "")."'\n";
}
$description .= " end select\n";
$description .= " return\n";
$descriptorFunction .= " end function ".$functionName."ID\n";
$descriptorFunction .= " ! End auto-generated enumeration function\n";
# Insert into the module.
my $descriptorWrapperTree = &Galacticus::Build::SourceTree::ParseCode($descriptorFunctionWrapper,"Galacticus::Build::SourceTree::Process::Enumeration()", instrument => 0);
my $descriptorTree = &Galacticus::Build::SourceTree::ParseCode($descriptorFunction ,"Galacticus::Build::SourceTree::Process::Enumeration()", instrument => 0);
my @descriptorNodes = &Galacticus::Build::SourceTree::Children ($descriptorTree );
my @descriptorWrapperNodes = &Galacticus::Build::SourceTree::Children ($descriptorWrapperTree );
my $newNode = $descriptorNodes[0];
while ( $newNode->{'type'} ne "function" ) {
$newNode = $newNode->{'sibling'};
}
$newNode = $newNode->{'firstChild'};
while ( defined($newNode->{'sibling'}) ) {
$newNode = $newNode->{'sibling'};
}
my $describeNode =
{
type => "code",
content => $description
};
&Galacticus::Build::SourceTree::InsertAfterNode($newNode,[$describeNode]);
&Galacticus::Build::SourceTree::InsertPostContains($node->{'parent'},\@descriptorNodes );
&Galacticus::Build::SourceTree::InsertPostContains($node->{'parent'},\@descriptorWrapperNodes);
my $interfaceTree = &Galacticus::Build::SourceTree::ParseCode($interface,"Galacticus::Build::SourceTree::Process::Enumeration()");
my @interfaceNodes = &Galacticus::Build::SourceTree::Children($interfaceTree);
&Galacticus::Build::SourceTree::InsertPreContains($node->{'parent'},\@interfaceNodes);
# Set the visibility.
&Galacticus::Build::SourceTree::SetVisibility($node->{'parent'},$functionName,$visibility);
}
# Create describe function.
{
my $functionName = "enumeration".ucfirst($node->{'directive'}->{'name'})."Describe";
my $descriptorFunction;
$descriptorFunction .= "\n";
$descriptorFunction .= " ! Auto-generated enumeration function\n";
$descriptorFunction .= " function ".$functionName."() result(description)\n";
$descriptorFunction .= " !!{\n";
$descriptorFunction .= " Return a description of a {\\normalfont \\ttfamily ".$node->{'directive'}->{'name'}."} enumeration value.\n";
$descriptorFunction .= " Return a description of the {\\normalfont \\ttfamily ".$node->{'directive'}->{'name'}."} enumeration.\n";
$descriptorFunction .= " !!}\n";
$descriptorFunction .= " use :: ISO_Varying_String, only : varying_string, var_str, operator(//)\n";
$descriptorFunction .= " implicit none\n";
Expand Down
92 changes: 81 additions & 11 deletions source/nodes.operators.physics.satellite_merging.radius_trigger.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@

use :: Dark_Matter_Halo_Scales, only : darkMatterHaloScaleClass
use :: Galactic_Structure , only : galacticStructureClass
use :: Kepler_Orbits , only : keplerOrbitCount

!![
<nodeOperator name="nodeOperatorSatelliteMergingRadiusTrigger">
Expand All @@ -34,9 +35,11 @@
A node operator class that triggers merging of satellites based on their orbital radius.
!!}
private
class (darkMatterHaloScaleClass), pointer :: darkMatterHaloScale_ => null()
class (galacticStructureClass ), pointer :: galacticStructure_ => null()
class (darkMatterHaloScaleClass), pointer :: darkMatterHaloScale_ => null()
class (galacticStructureClass ), pointer :: galacticStructure_ => null()
double precision :: radiusVirialFraction
logical :: recordMergedSubhaloProperties
integer :: mergedSubhaloIDs (keplerOrbitCount)
contains
!![
<methods>
Expand All @@ -55,6 +58,10 @@
module procedure satelliteMergingRadiusTriggerConstructorParameters
module procedure satelliteMergingRadiusTriggerConstructorInternal
end interface nodeOperatorSatelliteMergingRadiusTrigger

! Sub-module-scope pointer to self used in callback function.
class(nodeOperatorSatelliteMergingRadiusTrigger), pointer :: self_
!$omp threadprivate(self)

contains

Expand All @@ -69,6 +76,7 @@ function satelliteMergingRadiusTriggerConstructorParameters(parameters) result(s
class (darkMatterHaloScaleClass ), pointer :: darkMatterHaloScale_
class (galacticStructureClass ), pointer :: galacticStructure_
double precision :: radiusVirialFraction
logical :: recordMergedSubhaloProperties

!![
<inputParameter>
Expand All @@ -77,10 +85,16 @@ function satelliteMergingRadiusTriggerConstructorParameters(parameters) result(s
<description>The fraction of the virial radius below which satellites are merged.</description>
<source>parameters</source>
</inputParameter>
<inputParameter>
<name>recordMergedSubhaloProperties</name>
<defaultValue>.false.</defaultValue>
<description>If true, record the orbital properties of subhalo that merge.</description>
<source>parameters</source>
</inputParameter>
<objectBuilder class="darkMatterHaloScale" name="darkMatterHaloScale_" source="parameters"/>
<objectBuilder class="galacticStructure" name="galacticStructure_" source="parameters"/>
!!]
self=nodeOperatorSatelliteMergingRadiusTrigger(radiusVirialFraction,darkMatterHaloScale_,galacticStructure_)
self=nodeOperatorSatelliteMergingRadiusTrigger(radiusVirialFraction,recordMergedSubhaloProperties,darkMatterHaloScale_,galacticStructure_)
!![
<inputParametersValidate source="parameters"/>
<objectDestructor name="darkMatterHaloScale_"/>
Expand All @@ -89,19 +103,31 @@ function satelliteMergingRadiusTriggerConstructorParameters(parameters) result(s
return
end function satelliteMergingRadiusTriggerConstructorParameters

function satelliteMergingRadiusTriggerConstructorInternal(radiusVirialFraction,darkMatterHaloScale_,galacticStructure_) result(self)
function satelliteMergingRadiusTriggerConstructorInternal(radiusVirialFraction,recordMergedSubhaloProperties,darkMatterHaloScale_,galacticStructure_) result(self)
!!{
Internal constructor for the {\normalfont \ttfamily satelliteMergingRadiusTrigger} node operator class.
!!}
use :: Kepler_Orbits, only : keplerOrbitTimeInitial , keplerOrbitMassSatellite, keplerOrbitMassHost, keplerOrbitRadius, &
& keplerOrbitRadiusPericenter
implicit none
type (nodeOperatorSatelliteMergingRadiusTrigger) :: self
double precision , intent(in ) :: radiusVirialFraction
logical , intent(in ) :: recordMergedSubhaloProperties
class (darkMatterHaloScaleClass ), intent(in ), target :: darkMatterHaloScale_
class (galacticStructureClass ), intent(in ), target :: galacticStructure_
!![
<constructorAssign variables="radiusVirialFraction, *darkMatterHaloScale_, *galacticStructure_"/>
<constructorAssign variables="radiusVirialFraction, recordMergedSubhaloProperties, *darkMatterHaloScale_, *galacticStructure_"/>
!!]

if (recordMergedSubhaloProperties) then
!![
<addMetaProperty component="basic" name="mergedSubhaloTimeInitial" id="self%mergedSubhaloIDs(keplerOrbitTimeInitial %ID)" rank="1" isCreator="yes"/>
<addMetaProperty component="basic" name="mergedSubhaloMassSatellite" id="self%mergedSubhaloIDs(keplerOrbitMassSatellite %ID)" rank="1" isCreator="yes"/>
<addMetaProperty component="basic" name="mergedSubhaloMassHost" id="self%mergedSubhaloIDs(keplerOrbitMassHost %ID)" rank="1" isCreator="yes"/>
<addMetaProperty component="basic" name="mergedSubhaloRadius" id="self%mergedSubhaloIDs(keplerOrbitRadius %ID)" rank="1" isCreator="yes"/>
<addMetaProperty component="basic" name="mergedSubhaloRadiusPericenter" id="self%mergedSubhaloIDs(keplerOrbitRadiusPericenter%ID)" rank="1" isCreator="yes"/>
!!]
end if
return
end function satelliteMergingRadiusTriggerConstructorInternal

Expand All @@ -111,7 +137,7 @@ subroutine satelliteMergingRadiusTriggerDestructor(self)
!!}
implicit none
type(nodeOperatorSatelliteMergingRadiusTrigger), intent(inout) :: self

!![
<objectDestructor name="self%darkMatterHaloScale_"/>
<objectDestructor name="self%galacticStructure_" />
Expand Down Expand Up @@ -149,6 +175,7 @@ subroutine satelliteMergingRadiusTriggerDifferentialEvolution(self,node,interrup
! Merging criterion met - trigger an interrupt.
interrupt = .true.
functionInterrupt => mergerTrigger
self_ => self
end if
return
end subroutine satelliteMergingRadiusTriggerDifferentialEvolution
Expand All @@ -157,17 +184,60 @@ subroutine mergerTrigger(node,timeEnd)
!!{
Trigger a merger of the satellite by setting the time until merging to zero.
!!}
use :: Galacticus_Nodes, only : nodeComponentSatellite, nodeComponentBasic, treeNode
use :: Galacticus_Nodes, only : nodeComponentSatellite , nodeComponentBasic , treeNode
use :: Kepler_Orbits , only : keplerOrbit , keplerOrbitTimeInitial, keplerOrbitMassSatellite, keplerOrbitMassHost, &
& keplerOrbitRadiusPericenter, keplerOrbitRadius
implicit none
type (treeNode ), intent(inout), target :: node
double precision , intent(in ), optional :: timeEnd
class (nodeComponentBasic ) , pointer :: basic
class (nodeComponentSatellite) , pointer :: satellite
type (treeNode ), intent(inout), target :: node
double precision , intent(in ), optional :: timeEnd
type (treeNode ) , pointer :: nodeHost
class (nodeComponentBasic ) , pointer :: basic , basicHost
class (nodeComponentSatellite) , pointer :: satellite
double precision , dimension(:) , allocatable :: propertyCurrent, propertyNew
double precision :: property
type (keplerOrbit ) :: orbit
integer :: i , ID
!$GLC attributes unused :: timeEnd

! Set the time of merging to the current time.
basic => node%basic ()
satellite => node%satellite()
call satellite%timeOfMergingSet(basic%time())
! Record properties of the merging subhalo if necessary.
if (self_%recordMergedSubhaloProperties) then
! Find the node to merge with.
nodeHost => node %mergesWith()
basicHost => nodeHost%basic ()
! Get the virial orbit of the halo about to merge.
orbit=satellite%virialOrbit()
! Append the orbit data.
do i=1,5
select case (i)
case (1)
ID =keplerOrbitTimeInitial %ID
property=basic%timeLastIsolated()
case (2)
ID =keplerOrbitMassSatellite %ID
property=orbit%massSatellite ()
case (3)
ID =keplerOrbitMassHost %ID
property=orbit%massHost ()
case (4)
ID =keplerOrbitRadius %ID
property=orbit%radius ()
case (5)
ID =keplerOrbitRadiusPericenter%ID
property=orbit%radiusPericenter()
end select
propertyCurrent=basicHost%floatRank1MetaPropertyGet(self_%mergedSubhaloIDs(ID))
allocate(propertyNew(size(propertyCurrent)+1_c_size_t))
propertyNew(1_c_size_t:size(propertyCurrent))=propertyCurrent(:)
propertyNew(size(propertyNew))=property
call basicHost%floatRank1MetaPropertySet(self_%mergedSubhaloIDs(ID),propertyNew)
deallocate(propertyCurrent)
deallocate(propertyNew )
end do
end if
return
end subroutine mergerTrigger

Expand Down
Loading

0 comments on commit fc98fac

Please sign in to comment.