fof.h File Reference

#include "forcetree.h"
#include "utils/paramset.h"
#include "timestep.h"
#include "slotsmanager.h"

Include dependency graph for fof.h:

This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes
struct	BaseGroup
struct	Group
struct	FOFGroups

Typedefs
typedef struct FOFGroups	FOFGroups

Functions
void	set_fof_params (ParameterSet *ps)
void	fof_init (double DMMeanSeparation)
FOFGroups	fof_fof (DomainDecomp *ddecomp, const int StoreGrNr, MPI_Comm Comm)
void	fof_finish (FOFGroups *fof)
void	fof_seed (FOFGroups *fof, ActiveParticles *act, double atime, MPI_Comm Comm)
void	fof_save_groups (FOFGroups *fof, const char *OutputDir, const char *FOFFileBase, int num, Cosmology *CP, double atime, const double *MassTable, int MetalReturnOn, int BlackholeOn, MPI_Comm Comm)
void	fof_save_particles (FOFGroups *fof, const char *OutputDir, const char *FOFFileBase, int num, int SaveParticles, Cosmology *CP, double atime, const double *MassTable, int MetalReturnOn, int BlackholeOn, MPI_Comm Comm)

Typedef Documentation

FOFGroups

typedef struct FOFGroups FOFGroups

Function Documentation

fof_finish()

void fof_finish

(

FOFGroups *

fof

)

Definition at line 257 of file fof.c.

{
    myfree(fof->Group);
 
    message(0, "Finished computing FoF groups.  (presently allocated=%g MB)\n",
            mymalloc_usedbytes() / (1024.0 * 1024.0));
 
    walltime_measure("/FOF/MISC");
 
    MPI_Type_free(&MPI_TYPE_GROUP);
}

References FOFGroups::Group, message(), MPI_TYPE_GROUP, myfree, mymalloc_usedbytes, and walltime_measure.

Referenced by run(), runfof(), and test_fof().

Here is the call graph for this function:

Here is the caller graph for this function:

fof_fof()

FOFGroups fof_fof	(	DomainDecomp *	ddecomp,
		const int	StoreGrNr,
		MPI_Comm	Comm
	)

Definition at line 151 of file fof.c.

{
    int i;
 
    message(0, "Begin to compute FoF group catalogues. (allocated: %g MB)\n",
            mymalloc_usedbytes() / (1024.0 * 1024.0));
 
    walltime_measure("/Misc");
 
    message(0, "Comoving linking length: %g\n", fof_params.FOFHaloComovingLinkingLength);
 
    struct fof_particle_list * HaloLabel = (struct fof_particle_list *) mymalloc("HaloLabel", PartManager->NumPart * sizeof(struct fof_particle_list));
 
    /* HaloLabel stores the MinID and MinIDTask of particles, this pair serves as a halo label. */
    #pragma omp parallel for
    for(i = 0; i < PartManager->NumPart; i++) {
        HaloLabel[i].Pindex = i;
    }
 
    /* We only need a tree containing primary linking particles only. No moments*/
    ForceTree dmtree = {0};
    force_tree_rebuild_mask(&dmtree, ddecomp, fof_params.FOFPrimaryLinkTypes, 0, NULL);
 
    /* Fill FOFP_List of primary */
    fof_label_primary(HaloLabel, &dmtree, Comm);
 
    MPIU_Barrier(Comm);
    message(0, "Group finding done.\n");
    walltime_measure("/FOF/Primary");
 
    /* Fill FOFP_List of secondary */
    fof_label_secondary(HaloLabel, &dmtree);
    force_tree_free(&dmtree);
 
    MPIU_Barrier(Comm);
    message(0, "Attached gas and star particles to nearest dm particles.\n");
 
    walltime_measure("/FOF/Secondary");
 
    /* sort HaloLabel according to MinID, because we need that for compiling catalogues */
    qsort_openmp(HaloLabel, PartManager->NumPart, sizeof(struct fof_particle_list), fof_compare_HaloLabel_MinID);
 
    int NgroupsExt = 0;
 
    for(i = 0; i < PartManager->NumPart; i ++) {
        if(i == 0 || HaloLabel[i].MinID != HaloLabel[i - 1].MinID) NgroupsExt ++;
    }
 
    /* The first round is to eliminate groups that are too short. */
    /* We create the smaller 'BaseGroup' data set for this. */
    struct BaseGroup * base = (struct BaseGroup *) mymalloc("BaseGroup", sizeof(struct BaseGroup) * NgroupsExt);
 
    NgroupsExt = fof_compile_base(base, NgroupsExt, HaloLabel, Comm);
 
    MPIU_Barrier(Comm);
    message(0, "Compiled local group data and catalogue.\n");
 
    walltime_measure("/FOF/Compile");
 
    fof_assign_grnr(base, NgroupsExt, Comm);
 
    /*Store the group number in the particle struct*/
    if(StoreGrNr) {
        #pragma omp parallel for
        for(i = 0; i < PartManager->NumPart; i++)
            P[i].GrNr = -1; /* will mark particles that are not in any group */
 
        int64_t start = 0;
        for(i = 0; i < NgroupsExt; i++)
        {
            for(;start < PartManager->NumPart; start++) {
                if (HaloLabel[start].MinID >= base[i].MinID)
                    break;
            }
 
            for(;start < PartManager->NumPart; start++) {
                if (HaloLabel[start].MinID != base[i].MinID)
                    break;
                P[HaloLabel[start].Pindex].GrNr = base[i].GrNr;
            }
        }
    }
 
    /*Initialise the Group object from the BaseGroup*/
    FOFGroups fof;
    MPI_Type_contiguous(sizeof(fof.Group[0]), MPI_BYTE, &MPI_TYPE_GROUP);
    MPI_Type_commit(&MPI_TYPE_GROUP);
 
    fof.Group = fof_alloc_group(base, NgroupsExt);
 
    myfree(base);
 
    fof_compile_catalogue(&fof, NgroupsExt, HaloLabel, Comm);
 
    MPIU_Barrier(Comm);
    message(0, "Finished FoF. Group properties are now allocated.. (presently allocated=%g MB)\n",
            mymalloc_usedbytes() / (1024.0 * 1024.0));
 
    walltime_measure("/FOF/Prop");
 
    myfree(HaloLabel);
 
    return fof;
}

References fof_alloc_group(), fof_assign_grnr(), fof_compare_HaloLabel_MinID(), fof_compile_base(), fof_compile_catalogue(), fof_label_primary(), fof_label_secondary(), fof_params, FOFParams::FOFHaloComovingLinkingLength, FOFParams::FOFPrimaryLinkTypes, force_tree_free(), force_tree_rebuild_mask(), BaseGroup::GrNr, FOFGroups::Group, message(), fof_particle_list::MinID, BaseGroup::MinID, MPI_TYPE_GROUP, MPIU_Barrier, myfree, mymalloc, mymalloc_usedbytes, part_manager_type::NumPart, P, PartManager, fof_particle_list::Pindex, qsort_openmp, and walltime_measure.

Referenced by run(), runfof(), and test_fof().

Here is the call graph for this function:

Here is the caller graph for this function:

fof_init()

void fof_init

(

double

DMMeanSeparation

)

Definition at line 66 of file fof.c.

{
    fof_params.FOFHaloComovingLinkingLength = fof_params.FOFHaloLinkingLength * DMMeanSeparation;
}

References fof_params, FOFParams::FOFHaloComovingLinkingLength, and FOFParams::FOFHaloLinkingLength.

Referenced by init().

Here is the caller graph for this function:

fof_save_groups()

void fof_save_groups	(	FOFGroups *	fof,
		const char *	OutputDir,
		const char *	FOFFileBase,
		int	num,
		Cosmology *	CP,
		double	atime,
		const double *	MassTable,
		int	MetalReturnOn,
		int	BlackholeOn,
		MPI_Comm	Comm
	)

Definition at line 1113 of file fof.c.

{
    fof_save_particles(fof, OutputDir, FOFFileBase, num, fof_params.FOFSaveParticles, CP, atime, MassTable, MetalReturnOn, BlackholeOn, Comm);
}

References CP, fof_params, fof_save_particles(), and FOFParams::FOFSaveParticles.

Referenced by run(), and runfof().

Here is the call graph for this function:

Here is the caller graph for this function:

fof_save_particles()

void fof_save_particles	(	FOFGroups *	fof,
		const char *	OutputDir,
		const char *	FOFFileBase,
		int	num,
		int	SaveParticles,
		Cosmology *	CP,
		double	atime,
		const double *	MassTable,
		int	MetalReturnOn,
		int	BlackholeOn,
		MPI_Comm	Comm
	)

Definition at line 32 of file fofpetaio.c.

                                                                                                                                                                                                                                  {
    int i;
    struct IOTable FOFIOTable = {0};
    char * fname = fastpm_strdup_printf("%s/%s_%03d", OutputDir, FOFFileBase, num);
    message(0, "Saving particle groups into %s\n", fname);
 
    fof_register_io_blocks(MetalReturnOn, BlackholeOn, &FOFIOTable);
    /* sort the groups according to group-number */
    mpsort_mpi(fof->Group, fof->Ngroups, sizeof(struct Group),
            fof_radix_Group_GrNr, 8, NULL, Comm);
 
    BigFile bf = {0};
    if(0 != big_file_mpi_create(&bf, fname, Comm)) {
        endrun(0, "Failed to open file at %s\n", fname);
    }
    myfree(fname);
    struct conversions conv = {0};
    conv.atime = atime;
    conv.hubble = hubble_function(CP, atime);
 
    MPIU_Barrier(Comm);
    fof_write_header(&bf, fof->TotNgroups, atime, MassTable, CP, Comm);
 
    for(i = 0; i < FOFIOTable.used; i ++) {
        /* only process the particle blocks */
        char blockname[128];
        int ptype = FOFIOTable.ent[i].ptype;
        BigArray array = {0};
        if(ptype == PTYPE_FOF_GROUP) {
            sprintf(blockname, "FOFGroups/%s", FOFIOTable.ent[i].name);
            build_buffer_fof(fof, &array, &FOFIOTable.ent[i], &conv);
            message(0, "Writing Block %s\n", blockname);
 
            petaio_save_block(&bf, blockname, &array, 1);
            petaio_destroy_buffer(&array);
        }
    }
    destroy_io_blocks(&FOFIOTable);
    walltime_measure("/FOF/IO/WriteFOF");
 
    if(SaveParticles) {
        struct IOTable IOTable = {0};
        register_io_blocks(&IOTable, 1, MetalReturnOn);
        struct part_manager_type halo_pman = {0};
        struct slots_manager_type halo_sman = {0};
        if(fof_distribute_particles(&halo_pman, &halo_sman, Comm)) {
            myfree(halo_sman.Base);
            myfree(halo_pman.Base);
            destroy_io_blocks(&IOTable);
            return;
        }
 
        int * selection = (int *) mymalloc("Selection", sizeof(int) * halo_pman.NumPart);
 
        int ptype_offset[6]={0};
        int ptype_count[6]={0};
        petaio_build_selection(selection, ptype_offset, ptype_count, halo_pman.Base, halo_pman.NumPart, NULL);
 
        walltime_measure("/FOF/IO/argind");
 
        for(i = 0; i < IOTable.used; i ++) {
            /* only process the particle blocks */
            char blockname[128];
            int ptype = IOTable.ent[i].ptype;
            BigArray array = {0};
            if(ptype < 6 && ptype >= 0) {
                sprintf(blockname, "%d/%s", ptype, IOTable.ent[i].name);
                petaio_build_buffer(&array, &IOTable.ent[i], selection + ptype_offset[ptype], ptype_count[ptype], halo_pman.Base, &halo_sman, &conv);
 
                message(0, "Writing Block %s\n", blockname);
 
                petaio_save_block(&bf, blockname, &array, 1);
                petaio_destroy_buffer(&array);
            }
        }
        myfree(selection);
        myfree(halo_sman.Base);
        myfree(halo_pman.Base);
        walltime_measure("/FOF/IO/WriteParticles");
        destroy_io_blocks(&IOTable);
    }
 
    big_file_mpi_close(&bf, Comm);
 
    message(0, "Group catalogues saved.\n");
}

References conversions::atime, part_manager_type::Base, slots_manager_type::Base, build_buffer_fof(), CP, destroy_io_blocks(), endrun(), IOTable::ent, fastpm_strdup_printf(), fof_distribute_particles(), fof_radix_Group_GrNr(), fof_register_io_blocks(), fof_write_header(), FOFGroups::Group, conversions::hubble, hubble_function(), message(), MPIU_Barrier, mpsort_mpi, myfree, mymalloc, IOTableEntry::name, FOFGroups::Ngroups, part_manager_type::NumPart, petaio_build_buffer(), petaio_build_selection(), petaio_destroy_buffer(), petaio_save_block(), IOTableEntry::ptype, ptype, PTYPE_FOF_GROUP, register_io_blocks(), FOFGroups::TotNgroups, IOTable::used, and walltime_measure.

Referenced by fof_save_groups().

Here is the call graph for this function:

Here is the caller graph for this function:

fof_seed()

void fof_seed	(	FOFGroups *	fof,
		ActiveParticles *	act,
		double	atime,
		MPI_Comm	Comm
	)

Definition at line 1290 of file fof.c.

{
    int i, j, n, ntot;
 
    int NTask;
    MPI_Comm_size(Comm, &NTask);
 
    char * Marked = (char *) mymalloc2("SeedMark", fof->Ngroups);
 
    int Nexport = 0;
    #pragma omp parallel for reduction(+:Nexport)
    for(i = 0; i < fof->Ngroups; i++)
    {
        Marked[i] =
            (fof->Group[i].Mass >= fof_params.MinFoFMassForNewSeed)
        &&  (fof->Group[i].MassType[4] >= fof_params.MinMStarForNewSeed)
        &&  (fof->Group[i].LenType[5] == 0)
        &&  (fof->Group[i].seed_index >= 0);
 
        if(Marked[i]) Nexport ++;
    }
    struct Group * ExportGroups = (struct Group *) mymalloc("Export", sizeof(fof->Group[0]) * Nexport);
    j = 0;
    for(i = 0; i < fof->Ngroups; i ++) {
        if(Marked[i]) {
            ExportGroups[j] = fof->Group[i];
            j++;
        }
    }
    myfree(Marked);
 
    qsort_openmp(ExportGroups, Nexport, sizeof(ExportGroups[0]), cmp_seed_task);
 
    int * Send_count = ta_malloc("Send_count", int, NTask);
    int * Recv_count = ta_malloc("Recv_count", int, NTask);
 
    memset(Send_count, 0, NTask * sizeof(int));
    for(i = 0; i < Nexport; i++) {
        Send_count[ExportGroups[i].seed_task]++;
    }
 
    MPI_Alltoall(Send_count, 1, MPI_INT, Recv_count, 1, MPI_INT, Comm);
 
    int Nimport = 0;
 
    for(j = 0;  j < NTask; j++)
    {
        Nimport += Recv_count[j];
    }
 
    struct Group * ImportGroups = (struct Group *)
            mymalloc2("ImportGroups", Nimport * sizeof(struct Group));
 
    MPI_Alltoallv_smart(ExportGroups, Send_count, NULL, MPI_TYPE_GROUP,
                        ImportGroups, Recv_count, NULL, MPI_TYPE_GROUP,
                        Comm);
 
    myfree(ExportGroups);
    ta_free(Recv_count);
    ta_free(Send_count);
 
    MPI_Allreduce(&Nimport, &ntot, 1, MPI_INT, MPI_SUM, Comm);
 
    message(0, "Making %d new black hole particles.\n", ntot);
 
    /* Do we have enough black hole slots to create this many black holes?
     * If not, allocate more slots. */
    if(Nimport + SlotsManager->info[5].size > SlotsManager->info[5].maxsize)
    {
        int *ActiveParticle_tmp=NULL;
        /* This is only called on a PM step, so the condition should never be true*/
        if(act->ActiveParticle) {
            ActiveParticle_tmp = (int *) mymalloc2("ActiveParticle_tmp", act->NumActiveParticle * sizeof(int));
            memmove(ActiveParticle_tmp, act->ActiveParticle, act->NumActiveParticle * sizeof(int));
            myfree(act->ActiveParticle);
        }
 
        /*Now we can extend the slots! */
        int64_t atleast[6];
        int64_t i;
        for(i = 0; i < 6; i++)
            atleast[i] = SlotsManager->info[i].maxsize;
        atleast[5] += ntot*1.1;
        slots_reserve(1, atleast, SlotsManager);
 
        /*And now we need our memory back in the right place*/
        if(ActiveParticle_tmp) {
            act->ActiveParticle = (int *) mymalloc("ActiveParticle", sizeof(int)*(act->NumActiveParticle + PartManager->MaxPart - PartManager->NumPart));
            memmove(act->ActiveParticle, ActiveParticle_tmp, act->NumActiveParticle * sizeof(int));
            myfree(ActiveParticle_tmp);
        }
    }
 
    int ThisTask;
    MPI_Comm_rank(Comm, &ThisTask);
 
    for(n = 0; n < Nimport; n++)
    {
        fof_seed_make_one(&ImportGroups[n], ThisTask, atime);
    }
 
    myfree(ImportGroups);
 
    walltime_measure("/FOF/Seeding");
}

References ActiveParticles::ActiveParticle, cmp_seed_task(), fof_params, fof_seed_make_one(), FOFGroups::Group, slots_manager_type::info, Group::LenType, Group::Mass, Group::MassType, part_manager_type::MaxPart, slot_info::maxsize, message(), FOFParams::MinFoFMassForNewSeed, FOFParams::MinMStarForNewSeed, MPI_Alltoallv_smart(), MPI_TYPE_GROUP, myfree, mymalloc, mymalloc2, FOFGroups::Ngroups, NTask, ActiveParticles::NumActiveParticle, part_manager_type::NumPart, PartManager, qsort_openmp, Group::seed_index, Group::seed_task, slot_info::size, slots_reserve(), SlotsManager, ta_free, ta_malloc, ThisTask, and walltime_measure.

Referenced by run().

Here is the call graph for this function:

Here is the caller graph for this function:

set_fof_params()

void set_fof_params

(

ParameterSet *

ps

)

Definition at line 50 of file fof.c.

{
    int ThisTask;
    MPI_Comm_rank(MPI_COMM_WORLD, &ThisTask);
    if(ThisTask == 0) {
        fof_params.FOFSaveParticles = param_get_int(ps, "FOFSaveParticles");
        fof_params.FOFHaloLinkingLength = param_get_double(ps, "FOFHaloLinkingLength");
        fof_params.FOFHaloMinLength = param_get_int(ps, "FOFHaloMinLength");
        fof_params.MinFoFMassForNewSeed = param_get_double(ps, "MinFoFMassForNewSeed");
        fof_params.MinMStarForNewSeed = param_get_double(ps, "MinMStarForNewSeed");
        fof_params.FOFPrimaryLinkTypes = param_get_int(ps, "FOFPrimaryLinkTypes");
        fof_params.FOFSecondaryLinkTypes = param_get_int(ps, "FOFSecondaryLinkTypes");
    }
    MPI_Bcast(&fof_params, sizeof(struct FOFParams), MPI_BYTE, 0, MPI_COMM_WORLD);
}

References fof_params, FOFParams::FOFHaloLinkingLength, FOFParams::FOFHaloMinLength, FOFParams::FOFPrimaryLinkTypes, FOFParams::FOFSaveParticles, FOFParams::FOFSecondaryLinkTypes, FOFParams::MinFoFMassForNewSeed, FOFParams::MinMStarForNewSeed, param_get_double(), param_get_int(), and ThisTask.

Referenced by read_parameter_file().

Here is the call graph for this function:

Here is the caller graph for this function: