Skip to content
Snippets Groups Projects
Commit 0ef836d1 authored by Stefan Arridge's avatar Stefan Arridge
Browse files

Added cooling to runner.c

parent f06e74fa
No related branches found
No related tags found
1 merge request!242Add the cooling infrastructure and the const_du and const_lambda
...@@ -64,7 +64,7 @@ void cooling_print(const struct cooling_data* cooling) { ...@@ -64,7 +64,7 @@ void cooling_print(const struct cooling_data* cooling) {
void update_entropy(const struct cooling_data* cooling, void update_entropy(const struct cooling_data* cooling,
const struct phys_const* const phys_const, struct part* p, const struct phys_const* const phys_const, struct part* p,
float dt){ double dt){
/*updates the entropy of a particle after integrating the cooling equation*/ /*updates the entropy of a particle after integrating the cooling equation*/
float u_old; float u_old;
...@@ -80,7 +80,7 @@ void update_entropy(const struct cooling_data* cooling, ...@@ -80,7 +80,7 @@ void update_entropy(const struct cooling_data* cooling,
} }
float calculate_new_thermal_energy(float u_old, float dt, const struct cooling_data* cooling){ float calculate_new_thermal_energy(float u_old, double dt, const struct cooling_data* cooling){
#ifdef CONST_COOLING #ifdef CONST_COOLING
//This function integrates the cooling equation, given the initial thermal energy and the timestep dt. //This function integrates the cooling equation, given the initial thermal energy and the timestep dt.
//Returns 0 if successful and 1 if not //Returns 0 if successful and 1 if not
......
...@@ -97,5 +97,5 @@ void cooling_print(const struct cooling_data* cooling); ...@@ -97,5 +97,5 @@ void cooling_print(const struct cooling_data* cooling);
float calculate_new_thermal_energy(float u_old, float dt, const struct cooling_data* cooling); float calculate_new_thermal_energy(float u_old, float dt, const struct cooling_data* cooling);
void update_entropy(const struct cooling_data* cooling, void update_entropy(const struct cooling_data* cooling,
const struct phys_const* const phys_const, struct part* p, const struct phys_const* const phys_const, struct part* p,
float dt); double dt);
#endif /* SWIFT_COOLING_H */ #endif /* SWIFT_COOLING_H */
...@@ -135,6 +135,53 @@ void runner_do_grav_external(struct runner *r, struct cell *c, int timer) { ...@@ -135,6 +135,53 @@ void runner_do_grav_external(struct runner *r, struct cell *c, int timer) {
if (timer) TIMER_TOC(timer_dograv_external); if (timer) TIMER_TOC(timer_dograv_external);
} }
/**
* @brief Calculate change in entropy from cooling
*
* @param r runner task
* @param c cell
* @param timer 1 if the time is to be recorded.
*/
void runner_do_cooling(struct runner *r, struct cell *c, int timer) {
struct part *restrict parts = c->parts;
const int count = c->count;
const int ti_current = r->e->ti_current;
const struct cooling_data *cooling = r->e->cooling;
const struct phys_const *constants = r->e->physical_constants;
const double timeBase = r->e->timeBase;
double dt;
TIMER_TIC;
/* Recurse? */
if (c->split) {
for (int k = 0; k < 8; k++)
if (c->progeny[k] != NULL) runner_do_cooling(r, c->progeny[k], 0);
return;
}
#ifdef TASK_VERBOSE
OUT;
#endif
/* Loop over the gparts in this cell. */
for (int i = 0; i < count; i++) {
/* Get a direct pointer on the part. */
struct part *restrict p = &parts[i];
/* Is this part within the time step? */
if (p->ti_end <= ti_current) {
dt = (p->ti_end - p->ti_begin)*timeBase;
update_entropy(cooling, constants, p, dt);
}
}
if (timer) TIMER_TOC(timer_do_cooling);
}
/** /**
* @brief Sort the entries in ascending order using QuickSort. * @brief Sort the entries in ascending order using QuickSort.
* *
...@@ -1158,6 +1205,9 @@ void *runner_main(void *data) { ...@@ -1158,6 +1205,9 @@ void *runner_main(void *data) {
case task_type_grav_external: case task_type_grav_external:
runner_do_grav_external(r, t->ci, 1); runner_do_grav_external(r, t->ci, 1);
break; break;
case task_type_cooling:
runner_do_cooling(r, t->ci, 1);
break;
case task_type_part_sort: case task_type_part_sort:
space_do_parts_sort(); space_do_parts_sort();
break; break;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please to comment