hydra.h File Reference

#include "forcetree.h"
#include "types.h"
#include "timestep.h"
#include "density.h"
#include "utils/paramset.h"

Include dependency graph for hydra.h:

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

Go to the source code of this file.

Functions
void	hydro_force (const ActiveParticles *act, const double atime, struct sph_pred_data *SPH_predicted, double MinEgySpec, const DriftKickTimes times, Cosmology *CP, const ForceTree *const tree)
void	set_hydro_params (ParameterSet *ps)
int	DensityIndependentSphOn (void)

Function Documentation

DensityIndependentSphOn()

int DensityIndependentSphOn

(

void

)

Definition at line 49 of file hydra.c.

{
    return HydroParams.DensityIndependentSphOn;
}

References hydro_params::DensityIndependentSphOn, and HydroParams.

Referenced by register_io_blocks(), run(), setup_density_indep_entropy(), and setup_smoothinglengths().

Here is the caller graph for this function:

hydro_force()

void hydro_force	(	const ActiveParticles *	act,
		const double	atime,
		struct sph_pred_data *	SPH_predicted,
		double	MinEgySpec,
		const DriftKickTimes	times,
		Cosmology *	CP,
		const ForceTree *const	tree
	)

This function is the driver routine for the calculation of hydrodynamical force and rate of change of entropy due to shock heating for all active particles .

Definition at line 147 of file hydra.c.

{
    int i;
    TreeWalk tw[1] = {{0}};
 
    struct HydraPriv priv[1];
 
    tw->ev_label = "HYDRO";
    tw->visit = (TreeWalkVisitFunction) treewalk_visit_ngbiter;
    tw->ngbiter = (TreeWalkNgbIterFunction) hydro_ngbiter;
    tw->ngbiter_type_elsize = sizeof(TreeWalkNgbIterHydro);
    tw->haswork = hydro_haswork;
    tw->fill = (TreeWalkFillQueryFunction) hydro_copy;
    tw->reduce = (TreeWalkReduceResultFunction) hydro_reduce;
    tw->postprocess = (TreeWalkProcessFunction) hydro_postprocess;
    tw->query_type_elsize = sizeof(TreeWalkQueryHydro);
    tw->result_type_elsize = sizeof(TreeWalkResultHydro);
    tw->tree = tree;
    tw->priv = priv;
 
    if(!tree->hmax_computed_flag)
        endrun(5, "Hydro called before hmax computed\n");
    /* Cache the pressure for speed*/
    HYDRA_GET_PRIV(tw)->SPH_predicted = SPH_predicted;
    HYDRA_GET_PRIV(tw)->PressurePred = (double *) mymalloc("PressurePred", SlotsManager->info[0].size * sizeof(double));
    memset(HYDRA_GET_PRIV(tw)->PressurePred, 0, SlotsManager->info[0].size * sizeof(double));
 
    /* Compute pressure for active particles*/
    if(act->ActiveParticle) {
        #pragma omp parallel for
        for(i = 0; i < act->NumActiveParticle; i++) {
            int p_i = act->ActiveParticle[i];
            if(P[p_i].Type != 0)
                continue;
            int pi = P[p_i].PI;
            HYDRA_GET_PRIV(tw)->PressurePred[pi] = PressurePred(SPH_EOMDensity(&SphP[pi]), SPH_predicted->EntVarPred[pi]);
        }
    }
    else{
        /* Do it in slot order for memory locality*/
        #pragma omp parallel for
        for(i = 0; i < SlotsManager->info[0].size; i++)
            HYDRA_GET_PRIV(tw)->PressurePred[i] = PressurePred(SPH_EOMDensity(&SphP[i]), SPH_predicted->EntVarPred[i]);
    }
 
 
    double timeall = 0, timenetwork = 0;
    double timecomp, timecomm, timewait;
 
    walltime_measure("/SPH/Hydro/Init");
 
    /* Initialize some time factors*/
    const double hubble = hubble_function(CP, atime);
    HYDRA_GET_PRIV(tw)->fac_mu = pow(atime, 3 * (GAMMA - 1) / 2) / atime;
    HYDRA_GET_PRIV(tw)->fac_vsic_fix = hubble * pow(atime, 3 * GAMMA_MINUS1);
    HYDRA_GET_PRIV(tw)->atime = atime;
    HYDRA_GET_PRIV(tw)->hubble_a2 = hubble * atime * atime;
    priv->MinEgySpec = MinEgySpec;
    priv->times = &times;
    priv->FgravkickB = get_exact_gravkick_factor(CP, times.PM_kick, times.Ti_Current);
    memset(priv->gravkicks, 0, sizeof(priv->gravkicks[0])*(TIMEBINS+1));
    memset(priv->hydrokicks, 0, sizeof(priv->hydrokicks[0])*(TIMEBINS+1));
    memset(priv->drifts, 0, sizeof(priv->drifts[0])*(TIMEBINS+1));
    /* Compute the factors to move a current kick times velocity to the drift time velocity.
     * We need to do the computation for all timebins up to the maximum because even inactive
     * particles may have interactions. */
    #pragma omp parallel for
    for(i = times.mintimebin; i <= TIMEBINS; i++)
    {
        priv->gravkicks[i] = get_exact_gravkick_factor(CP, times.Ti_kick[i], times.Ti_Current);
        priv->hydrokicks[i] = get_exact_hydrokick_factor(CP, times.Ti_kick[i], times.Ti_Current);
        /* For density: last active drift time is Ti_kick - 1/2 timestep as the kick time is half a timestep ahead.
         * For active particles no density drift is needed.*/
        if(!is_timebin_active(i, times.Ti_Current))
            priv->drifts[i] = get_exact_drift_factor(CP, times.Ti_lastactivedrift[i], times.Ti_Current);
    }
 
    treewalk_run(tw, act->ActiveParticle, act->NumActiveParticle);
 
    myfree(HYDRA_GET_PRIV(tw)->PressurePred);
    /* collect some timing information */
 
    timeall += walltime_measure(WALLTIME_IGNORE);
 
    timecomp = tw->timecomp1 + tw->timecomp2 + tw->timecomp3;
    timewait = tw->timewait1 + tw->timewait2;
    timecomm = tw->timecommsumm1 + tw->timecommsumm2;
 
    walltime_add("/SPH/Hydro/Compute", timecomp);
    walltime_add("/SPH/Hydro/Wait", timewait);
    walltime_add("/SPH/Hydro/Comm", timecomm);
    walltime_add("/SPH/Hydro/Misc", timeall - (timecomp + timewait + timecomm + timenetwork));
}

References ActiveParticles::ActiveParticle, HydraPriv::atime, CP, HydraPriv::drifts, endrun(), sph_pred_data::EntVarPred, TreeWalk::ev_label, HydraPriv::FgravkickB, TreeWalk::fill, GAMMA, GAMMA_MINUS1, get_exact_drift_factor(), get_exact_gravkick_factor(), get_exact_hydrokick_factor(), HydraPriv::gravkicks, TreeWalk::haswork, ForceTree::hmax_computed_flag, hubble_function(), HYDRA_GET_PRIV, hydro_copy(), hydro_haswork(), hydro_ngbiter(), hydro_postprocess(), hydro_reduce(), HydraPriv::hydrokicks, slots_manager_type::info, is_timebin_active(), HydraPriv::MinEgySpec, DriftKickTimes::mintimebin, myfree, mymalloc, TreeWalk::ngbiter, TreeWalk::ngbiter_type_elsize, ActiveParticles::NumActiveParticle, P, DriftKickTimes::PM_kick, TreeWalk::postprocess, PressurePred(), TreeWalk::priv, TreeWalk::query_type_elsize, TreeWalk::reduce, TreeWalk::result_type_elsize, slot_info::size, SlotsManager, SPH_EOMDensity(), HydraPriv::SPH_predicted, SphP, DriftKickTimes::Ti_Current, DriftKickTimes::Ti_kick, DriftKickTimes::Ti_lastactivedrift, TIMEBINS, TreeWalk::timecommsumm1, TreeWalk::timecommsumm2, TreeWalk::timecomp1, TreeWalk::timecomp2, TreeWalk::timecomp3, HydraPriv::times, TreeWalk::timewait1, TreeWalk::timewait2, TreeWalk::tree, treewalk_run(), treewalk_visit_ngbiter(), TreeWalk::visit, walltime_add, WALLTIME_IGNORE, and walltime_measure.

Referenced by run().

Here is the call graph for this function:

Here is the caller graph for this function:

set_hydro_params()

void set_hydro_params

(

ParameterSet *

ps

)

Definition at line 37 of file hydra.c.

{
    int ThisTask;
    MPI_Comm_rank(MPI_COMM_WORLD, &ThisTask);
    if(ThisTask == 0) {
        HydroParams.ArtBulkViscConst = param_get_double(ps, "ArtBulkViscConst");
        HydroParams.DensityContrastLimit = param_get_double(ps, "DensityContrastLimit");
        HydroParams.DensityIndependentSphOn= param_get_int(ps, "DensityIndependentSphOn");
    }
    MPI_Bcast(&HydroParams, sizeof(struct hydro_params), MPI_BYTE, 0, MPI_COMM_WORLD);
}

References hydro_params::ArtBulkViscConst, hydro_params::DensityContrastLimit, hydro_params::DensityIndependentSphOn, HydroParams, 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: