12#include "gadgetconfig.h"
20#include "../data/allvars.h"
21#include "../data/dtypes.h"
22#include "../data/intposconvert.h"
23#include "../data/simparticles.h"
24#include "../lightcone/lightcone.h"
25#include "../logs/logs.h"
26#include "../logs/timer.h"
27#include "../main/main.h"
28#include "../main/simulation.h"
29#include "../mpi_utils/mpi_utils.h"
30#include "../system/system.h"
31#include "../time_integration/driftfac.h"
32#include "../time_integration/timestep.h"
47 for(
int bin = 0; bin <
TIMEBINS; bin++)
62 for(
int i = 0; i <
NumPart; i++)
81 if(
P[i].getType() == 0)
101 TimeBinSfr[bin] +=
SphP[i].Sfr;
133 ti_next_for_bin = (
All.
Ti_Current / dt_bin) * dt_bin + dt_bin;
140 if(ti_next_for_bin < ti_next_kick)
141 ti_next_kick = ti_next_for_bin;
146#ifdef ENLARGE_DYNAMIC_RANGE_IN_TIME
149 MPI_Allreduce(&ti_next_kick, &ti_next_kick_global, 1, MPI_INT, MPI_MIN,
Communicator);
152 return ti_next_kick_global;
157 int lowest_active_bin =
TIMEBINS, highest_active_bin = 0;
158 int lowest_occupied_bin =
TIMEBINS, highest_occupied_bin = 0;
159 int lowest_occupied_gravity_bin =
TIMEBINS, highest_occupied_gravity_bin = 0;
160 int highest_synchronized_bin = 0;
161 int nsynchronized_gravity = 0, nsynchronized_hydro = 0;
169 if(highest_occupied_gravity_bin < n)
170 highest_occupied_gravity_bin = n;
172 if(lowest_occupied_gravity_bin > n)
173 lowest_occupied_gravity_bin = n;
180 if(highest_occupied_bin < n)
181 highest_occupied_bin = n;
183 if(lowest_occupied_bin > n)
184 lowest_occupied_bin = n;
197 if(highest_synchronized_bin < n)
198 highest_synchronized_bin = n;
202 if(highest_active_bin < n)
203 highest_active_bin = n;
205 if(lowest_active_bin > n)
206 lowest_active_bin = n;
213 int lowest_in[3], lowest_out[3];
214 lowest_in[0] = lowest_occupied_bin;
215 lowest_in[1] = lowest_occupied_gravity_bin;
216 lowest_in[2] = lowest_active_bin;
217 MPI_Allreduce(lowest_in, lowest_out, 3, MPI_INT, MPI_MIN,
Communicator);
222 int highest_in[4], highest_out[4];
223 highest_in[0] = highest_occupied_bin;
224 highest_in[1] = highest_occupied_gravity_bin;
225 highest_in[2] = highest_active_bin;
226 highest_in[3] = highest_synchronized_bin;
227 MPI_Allreduce(highest_in, highest_out, 4, MPI_INT, MPI_MAX,
Communicator);
236 long long output_longs[2 + 2 *
TIMEBINS];
238 input_ints[0] = nsynchronized_hydro;
239 input_ints[1] = nsynchronized_gravity;
247 long long *tot_count_grav = output_longs + 2;
248 long long *tot_count_sph = output_longs + 2 +
TIMEBINS;
250 long long tot_grav = 0, tot_sph = 0;
254 tot_grav += tot_count_grav[n];
255 tot_sph += tot_count_sph[n];
259 tot_count_grav[n] += tot_count_grav[n - 1];
260 tot_count_sph[n] += tot_count_sph[n - 1];
278 for(
int i = 0; i <
NumPart; i++)
280#ifdef LIGHTCONE_MASSMAPS
285 MPI_Allreduce(MPI_IN_PLACE, &flag, 1, MPI_INT, MPI_MAX,
Communicator);
286 LightCone->lightcone_massmap_flush(0);
293#ifdef LIGHTCONE_MASSMAPS
298 MPI_Allreduce(MPI_IN_PLACE, &flag, 1, MPI_INT, MPI_MAX,
Communicator);
300 LightCone->lightcone_massmap_flush(0);
315 int buffer_full_flag = 0;
317 while(
P->
access.test_and_set(std::memory_order_acquire))
328 P->
access.clear(std::memory_order_release);
330 return buffer_full_flag;
334 Terminate(
"no prediction into past allowed: time0=%lld time1=%lld\n", (
long long)time0, (
long long)time1);
344 if(ignore_light_cone ==
false)
345 buffer_full_flag = LightCone->lightcone_test_for_particle_addition(
P, time0, time1, dt_drift);
349 for(
int j = 0; j < 3; j++)
350 posdiff[j] =
P->
Vel[j] * dt_drift;
355 for(
int j = 0; j < 3; j++)
362 double dt_hydrokick, dt_entr, dt_gravkick;
371 dt_gravkick = dt_entr = dt_hydrokick = dt_drift;
373 for(
int j = 0; j < 3; j++)
376#ifdef HIERARCHICAL_GRAVITY
381#if defined(PMGRID) && !defined(TREEPM_NOTIMESPLIT)
392#ifdef PRESSURE_ENTROPY_SPH
393 SphP->PressureSphDensity +=
SphP->DtPressureSphDensity * dt_drift;
404 P->
access.clear(std::memory_order_release);
407 return buffer_full_flag;
422 if(
P[i].getType() == 0)
424 if(
P[i].getTimeBinHydro() != n)
425 Terminate(
"P[i].TimeBinHydro=%d != timebin=%d",
P[i].getTimeBinHydro(), n);
427 if(
P[i].Ti_Current.load(std::memory_order_acquire) !=
All.
Ti_Current)
444 if(
P[i].Ti_Current.load(std::memory_order_acquire) !=
All.
Ti_Current)
global_data_all_processes All
double get_drift_factor(integertime time0, integertime time1)
double get_gravkick_factor(integertime time0, integertime time1)
double get_hydrokick_factor(integertime time0, integertime time1)
MySignedIntPosType pos_to_signedintpos(T posdiff)
void constrain_intpos(MyIntPosType *pos)
int drift_particle(particle_data *P, sph_particle_data *SphP, integertime time1, bool ignore_light_cone=false)
This function drifts a particle i to time1.
void reconstruct_timebins(void)
TimeBinData TimeBinsGravity
TimeBinData TimeBinsHydro
integertime find_next_sync_point(void)
This function finds the next synchronization point of the system. (i.e. the earliest point of time an...
int TimeBinSynchronized[TIMEBINS]
void make_list_of_active_particles(void)
void drift_all_particles(void)
void mark_active_timebins(void)
int32_t MySignedIntPosType
#define TIMER_START(counter)
Starts the timer counter.
#define TIMER_STOP(counter)
Stops the timer counter.
void sumup_large_ints(int n, int *src, long long *res, MPI_Comm comm)
void minimum_large_ints(int n, long long *src, long long *res, MPI_Comm comm)
expr pow(half base, half exp)
int TimeBinCount[TIMEBINS]
int FirstInTimeBin[TIMEBINS]
long long GlobalNActiveParticles
int LastInTimeBin[TIMEBINS]
long long GlobalNSynchronizedHydro
long long GlobalNSynchronizedGravity
int LowestOccupiedTimeBin
int HighestSynchronizedTimeBin
integertime Ti_begstep[TIMEBINS]
int ComovingIntegrationOn
int SmallestTimeBinWithDomainDecomposition
int LowestOccupiedGravTimeBin
double ActivePartFracForNewDomainDecomp
int HighestOccupiedGravTimeBin
int HighestOccupiedTimeBin
std::atomic< integertime > Ti_Current
unsigned char getTimeBinHydro(void)
unsigned char getType(void)
vector< MyFloat > GravAccel
void set_thermodynamic_variables(void)