Skip to content

Commit

Permalink
clean up variables.
Browse files Browse the repository at this point in the history
  • Loading branch information
@rainwoodman
rainwoodman committed Nov 25, 2014
1 parent 66aecba commit cd18953
Show file tree
Hide file tree
Showing 7 changed files with 79 additions and 138 deletions.
2 changes: 1 addition & 1 deletion allvars.h
View file
Expand Up @@ -786,7 +786,7 @@ extern struct global_data_all_processes
double a;
double a3inv;
double a2inv;
double afac;
double fac_egy;
double hubble;
double hubble_a2;
} cf;
Expand Down
34 changes: 9 additions & 25 deletions density.c
View file
Expand Up @@ -141,7 +141,6 @@ static void density_copy(int place, struct densdata_in * I);
* that one has to deal with substantially more than normal number of
* neighbours.)
*/
double a3inv, afac;

void density(void)
{
Expand Down Expand Up @@ -169,14 +168,6 @@ void density(void)

int64_t n_exported = 0;

if(All.ComovingIntegrationOn)
{
a3inv = 1 / (All.Time * All.Time * All.Time);
afac = pow(All.Time, 3 * GAMMA_MINUS1);
}
else
a3inv = afac = 1;

#if defined(EULERPOTENTIALS) || defined(VECT_PRO_CLEAN) || defined(TRACEDIVB) || defined(VECT_POTENTIAL)
double efak;

Expand Down Expand Up @@ -361,11 +352,7 @@ double density_decide_hsearch(int targettype, double h) {
* just like how it was done for grav smoothing.
* */
double rds;
if (All.ComovingIntegrationOn) {
rds = All.BlackHoleFeedbackRadiusMaxPhys / All.Time;
} else {
rds = All.BlackHoleFeedbackRadiusMaxPhys;
}
rds = All.BlackHoleFeedbackRadiusMaxPhys / All.cf.a;

if(rds > All.BlackHoleFeedbackRadius) {
rds = All.BlackHoleFeedbackRadius;
Expand Down Expand Up @@ -441,10 +428,7 @@ static void density_reduce(int place, struct densdata_out * remote, int mode) {
}
#ifdef HYDRO_COST_FACTOR
/* these will be added */
if(All.ComovingIntegrationOn)
P[place].GravCost += HYDRO_COST_FACTOR * All.Time * remote->Ninteractions;
else
P[place].GravCost += HYDRO_COST_FACTOR * remote->Ninteractions;
P[place].GravCost += HYDRO_COST_FACTOR * All.cf.a * remote->Ninteractions;
#endif

if(P[place].Type == 0)
Expand Down Expand Up @@ -837,9 +821,9 @@ static int density_evaluate(int target, int mode,
#pragma omp critical (_abundance_)
AbundanceRatios(DMAX(All.MinEgySpec,
SPHP(j).Entropy / GAMMA_MINUS1
* pow(SPHP(j).EOMDensity * a3inv,
* pow(SPHP(j).EOMDensity * All.cf.a3inv,
GAMMA_MINUS1)),
SPHP(j).d.Density * a3inv, &uvbg, &ne, &nh0, &nHeII);
SPHP(j).d.Density * All.cf.a3inv, &uvbg, &ne, &nh0, &nHeII);
#else
double nh0 = 1.0;
#endif
Expand Down Expand Up @@ -1090,29 +1074,29 @@ static void density_post_process(int i) {
#ifdef DENSITY_INDEPENDENT_SPH
#error pressure entropy incompatible with softereqs
/* use an intermediate EQS, between isothermal and the full multiphase model */
if(SPHP(i).d.Density * a3inv >= All.PhysDensThresh)
if(SPHP(i).d.Density * All.cf.a3inv >= All.PhysDensThresh)
SPHP(i).Pressure = All.FactorForSofterEQS *
(SPHP(i).Entropy + SPHP(i).e.DtEntropy * dt_entr) * pow(SPHP(i).d.Density, GAMMA) +
(1 -
All.FactorForSofterEQS) * afac * GAMMA_MINUS1 * SPHP(i).d.Density * All.InitGasU;
All.FactorForSofterEQS) * All.cf.fac_egy * GAMMA_MINUS1 * SPHP(i).d.Density * All.InitGasU;
else
SPHP(i).Pressure =
(SPHP(i).Entropy + SPHP(i).e.DtEntropy * dt_entr) * pow(SPHP(i).d.Density, GAMMA);
#else
/* use an intermediate EQS, between isothermal and the full multiphase model */
if(SPHP(i).d.Density * a3inv >= All.PhysDensThresh)
if(SPHP(i).d.Density * All.cf.a3inv >= All.PhysDensThresh)
SPHP(i).Pressure = All.FactorForSofterEQS *
(SPHP(i).Entropy + SPHP(i).e.DtEntropy * dt_entr) * pow(SPHP(i).d.Density, GAMMA) +
(1 -
All.FactorForSofterEQS) * afac * GAMMA_MINUS1 * SPHP(i).d.Density * All.InitGasU;
All.FactorForSofterEQS) * All.cf.fac_egy * GAMMA_MINUS1 * SPHP(i).d.Density * All.InitGasU;
else
SPHP(i).Pressure =
(SPHP(i).Entropy + SPHP(i).e.DtEntropy * dt_entr) * pow(SPHP(i).d.Density, GAMMA);
#endif // DENSITY_INDEPENDENT_SPH
#endif // SOFTEREQS
#else
/* Here we use an isothermal equation of state */
SPHP(i).Pressure = afac * GAMMA_MINUS1 * SPHP(i).d.Density * All.InitGasU;
SPHP(i).Pressure = All.cf.fac_egy * GAMMA_MINUS1 * SPHP(i).d.Density * All.InitGasU;
SPHP(i).Entropy = SPHP(i).Pressure / pow(SPHP(i).d.Density, GAMMA);
#endif
#else
Expand Down
2 changes: 1 addition & 1 deletion domain.c
View file
Expand Up @@ -985,7 +985,6 @@ int domain_fork_particle(int parent) {

NextActiveParticle[child] = FirstActiveParticle;
FirstActiveParticle = child;
NumForceUpdate++;

P[parent].Generation ++;
uint64_t g = P[parent].Generation;
Expand All @@ -1004,6 +1003,7 @@ int domain_fork_particle(int parent) {
NextInTimeBin[parent] = child;
if(LastInTimeBin[P[parent].TimeBin] == parent)
LastInTimeBin[P[parent].TimeBin] = child;
NumForceUpdate++;

/*! When a new additional star particle is created, we can put it into the
* tree at the position of the spawning gas particle. This is possible
Expand Down
10 changes: 1 addition & 9 deletions forcetree.c
View file
Expand Up @@ -719,14 +719,7 @@ void force_update_node_recursive(int no, int sib, int father)
#if defined(SFR) && defined(EDDINGTON_TENSOR_SFR)
if(pa->Type == 0)
{
double a;

if(All.ComovingIntegrationOn == 1)
a = All.Time;
else
a = 1.0;

if((SPHP(p).d.Density / (a * a * a)) >= All.PhysDensThresh)
if((SPHP(p).d.Density * All.cf.a3inv) >= All.PhysDensThresh)
{
stellar_s[0] += (pa->Mass * pa->Pos[0]);
stellar_s[1] += (pa->Mass * pa->Pos[1]);
Expand Down Expand Up @@ -1665,7 +1658,6 @@ static void real_force_drift_node(int no, int time1)
}
else
{

dt_drift_hmax = (time1 - Nodes[no].Ti_current) * All.Timebase_interval;
dt_drift = dt_drift_hmax;
}
Expand Down

0 comments on commit cd18953

Please sign in to comment.