diff --git a/examples/main.c b/examples/main.c
index d4a920db0a56066c7f89e172d15094dbfbfd4d28..b6a016761989da0fa909a4ebe951c6e46aef3f99 100644
--- a/examples/main.c
+++ b/examples/main.c
@@ -73,8 +73,7 @@ void print_help_message() {
"Overwrite the CPU frequency (Hz) to be used for time measurements");
printf(" %2s %8s %s\n", "-g", "",
"Run with an external gravitational potential");
- printf(" %2s %8s %s\n", "-F", "",
- "Run with feedback ");
+ printf(" %2s %8s %s\n", "-F", "", "Run with feedback ");
printf(" %2s %8s %s\n", "-G", "", "Run with self-gravity");
printf(" %2s %8s %s\n", "-n", "{int}",
"Execute a fixed number of time steps");
@@ -182,8 +181,8 @@ int main(int argc, char *argv[]) {
case 'G':
with_self_gravity = 1;
break;
- case 'F':
- with_sourceterms = 1;
+ case 'F':
+ with_sourceterms = 1;
break;
case 'h':
if (myrank == 0) print_help_message();
diff --git a/src/Makefile.am b/src/Makefile.am
index e1f0235cab7bb6111ab42bf69292fabc88f63fd8..731ce5687522b491756a22267c24e0a4001ff03a 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -42,8 +42,8 @@ endif
include_HEADERS = space.h runner.h queue.h task.h lock.h cell.h part.h const.h \
engine.h swift.h serial_io.h timers.h debug.h scheduler.h proxy.h parallel_io.h \
common_io.h single_io.h multipole.h map.h tools.h partition.h clocks.h parser.h \
- physical_constants.h physical_constants_cgs.h potentials.h version.h hydro_properties.h \
- sourceterms.h
+ physical_constants.h physical_constants_cgs.h potentials.h version.h \
+ hydro_properties.h sourceterms.h threadpool.h
# Common source files
AM_SOURCES = space.c runner.c queue.c task.c cell.c engine.c \
@@ -51,7 +51,7 @@ AM_SOURCES = space.c runner.c queue.c task.c cell.c engine.c \
units.c common_io.c single_io.c multipole.c version.c map.c \
kernel_hydro.c tools.c part.c partition.c clocks.c parser.c \
physical_constants.c potentials.c hydro_properties.c \
- runner_doiact_fft.c sourceterms.c
+ runner_doiact_fft.c sourceterms.c threadpool.c
# Include files for distribution, not installation.
nobase_noinst_HEADERS = approx_math.h atomic.h cycle.h error.h inline.h kernel_hydro.h kernel_gravity.h \
diff --git a/src/const.h b/src/const.h
index 7f882f828b5df491f6d2f88467b2acc383d8684d..22200fb88c37ea678d07c7288ad7107e4d1e1951 100644
--- a/src/const.h
+++ b/src/const.h
@@ -36,13 +36,17 @@
/* Time integration constants. */
#define const_max_u_change 0.1f
+/* Thermal energy per unit mass used as a constant for the isothermal EoS */
+#define const_isothermal_internal_energy 20.2615290634f
+
/* Dimensionality of the problem */
-#define HYDRO_DIMENSION_3D
-//#define HYDRO_DIMENSION_2D
+//#define HYDRO_DIMENSION_3D
+#define HYDRO_DIMENSION_2D
//#define HYDRO_DIMENSION_1D
/* Hydrodynamical adiabatic index. */
#define HYDRO_GAMMA_5_3
+//#define HYDRO_GAMMA_7_5
//#define HYDRO_GAMMA_4_3
//#define HYDRO_GAMMA_2_1
@@ -60,8 +64,24 @@
/* SPH variant to use */
//#define MINIMAL_SPH
-#define GADGET2_SPH
+//#define GADGET2_SPH
//#define DEFAULT_SPH
+#define GIZMO_SPH
+
+/* Riemann solver to use (GIZMO_SPH only) */
+#define RIEMANN_SOLVER_EXACT
+//#define RIEMANN_SOLVER_TRRS
+//#define RIEMANN_SOLVER_HLLC
+
+/* Type of gradients to use (GIZMO_SPH only) */
+/* If no option is chosen, no gradients are used (first order scheme) */
+//#define GRADIENTS_SPH
+#define GRADIENTS_GIZMO
+
+/* Types of slope limiter to use (GIZMO_SPH only) */
+/* Different slope limiters can be combined */
+#define SLOPE_LIMITER_PER_FACE
+#define SLOPE_LIMITER_CELL_WIDE
/* Self gravity stuff. */
#define const_gravity_multipole_order 2
@@ -70,8 +90,9 @@
#define const_gravity_eta 0.025f
/* External gravity properties */
-#define EXTERNAL_POTENTIAL_POINTMASS
+//#define EXTERNAL_POTENTIAL_POINTMASS
//#define EXTERNAL_POTENTIAL_ISOTHERMALPOTENTIAL
+#define EXTERNAL_POTENTIAL_DISK_PATCH
/* Source terms */
#define SN_FEEDBACK
diff --git a/src/engine.c b/src/engine.c
index 91a82bab529d250811bf270c5c3d6f587d108e4d..2a010c1eb82a1761444b15579cc0c793a07702b9 100644
--- a/src/engine.c
+++ b/src/engine.c
@@ -228,7 +228,11 @@ void engine_make_hydro_hierarchical_tasks(struct engine *e, struct cell *c,
/* Generate the ghost task. */
c->ghost = scheduler_addtask(s, task_type_ghost, task_subtype_none, 0, 0,
c, NULL, 0);
-
+#ifdef EXTRA_HYDRO_LOOP
+ /* Generate the extra ghost task. */
+ c->extra_ghost = scheduler_addtask(s, task_type_extra_ghost,
+ task_subtype_none, 0, 0, c, NULL, 0);
+#endif
/* add source terms */
if (is_with_sourceterms)
c->sourceterms = scheduler_addtask(s, task_type_sourceterms, task_subtype_none, 0, 0,
diff --git a/src/engine.h b/src/engine.h
index 164b0f1202f8760be92993fc45c280b8e6006d4b..5605131f491ed9d3823bb3f3c68ad41cdcf377a2 100644
--- a/src/engine.h
+++ b/src/engine.h
@@ -41,9 +41,9 @@
#include "parser.h"
#include "partition.h"
#include "potentials.h"
-#include "sourceterms.h"
#include "runner.h"
#include "scheduler.h"
+#include "sourceterms.h"
#include "space.h"
#include "task.h"
#include "units.h"
@@ -224,7 +224,7 @@ void engine_init(struct engine *e, struct space *s,
const struct phys_const *physical_constants,
const struct hydro_props *hydro,
const struct external_potential *potential,
- const struct sourceterms *sourceterms);
+ const struct sourceterms *sourceterms);
void engine_launch(struct engine *e, int nr_runners, unsigned int mask,
unsigned int submask);
void engine_prepare(struct engine *e);
diff --git a/src/hydro/Gizmo/hydro.h b/src/hydro/Gizmo/hydro.h
index f69dc3f1798f014e895c4a63760805b1739cec94..1858d622228f0d93f940f06ce03913285ab14ca3 100644
--- a/src/hydro/Gizmo/hydro.h
+++ b/src/hydro/Gizmo/hydro.h
@@ -618,3 +618,18 @@ __attribute__((always_inline)) INLINE static float hydro_get_internal_energy(
struct part* p) {
return 0.f;
}
+/**
+ * @brief Modifies the thermal state of a particle to the imposed internal
+ * energy
+ *
+ * This overrides the current state of the particle but does *not* change its
+ * time-derivatives
+ *
+ * @param p The particle
+ * @param u The new internal energy
+ */
+__attribute__((always_inline)) INLINE static void hydro_set_internal_energy(
+ struct part *restrict p, float u) {
+
+ p->conserved.energy = u;
+}
diff --git a/src/runner.c b/src/runner.c
index f8e03fba5c5ba674ca9c911c5fa4dde512dc00dc..56deb59a50094a45668f29afbb1c1ab1587ac462 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -1226,19 +1226,6 @@ void *runner_main(void *data) {
case task_type_sourceterms:
runner_do_sourceterms(r, t->ci, 1);
break;
- case task_type_part_sort:
- space_do_parts_sort();
- break;
- case task_type_gpart_sort:
- space_do_gparts_sort();
- break;
- case task_type_split_cell:
- space_do_split(e->s, t->ci);
- break;
- case task_type_rewait:
- scheduler_do_rewait((struct task *)t->ci, (struct task *)t->cj,
- t->flags, t->rank);
- break;
default:
error("Unknown task type.");
}
diff --git a/src/sourceterms.c b/src/sourceterms.c
index a3d32d829a02910ce9d98c6dd14ed905b21e2626..1ad903509f80e9192e4011748800a2207065e390 100644
--- a/src/sourceterms.c
+++ b/src/sourceterms.c
@@ -26,11 +26,11 @@
/* Local includes. */
#include "const.h"
+#include "equation_of_state.h"
#include "error.h"
#include "parser.h"
#include "part.h"
#include "physical_constants.h"
-#include "equation_of_state.h"
#include "units.h"
/* This object's header. */
@@ -44,15 +44,15 @@
* @param sourceterms: the structure that has all the source term properties
*/
void source_terms_init(const struct swift_params* parameter_file,
- struct UnitSystem* us,
- struct sourceterms* source) {
+ struct UnitSystem* us, struct sourceterms* source) {
#ifdef SN_FEEDBACK
- source->supernova.time = parser_get_param_double(parameter_file, "SN:time");
- source->supernova.energy = parser_get_param_double(parameter_file, "SN:energy");
- source->supernova.x = parser_get_param_double(parameter_file, "SN:x");
- source->supernova.y = parser_get_param_double(parameter_file, "SN:y");
- source->supernova.z = parser_get_param_double(parameter_file, "SN:z");
+ source->supernova.time = parser_get_param_double(parameter_file, "SN:time");
+ source->supernova.energy =
+ parser_get_param_double(parameter_file, "SN:energy");
+ source->supernova.x = parser_get_param_double(parameter_file, "SN:x");
+ source->supernova.y = parser_get_param_double(parameter_file, "SN:y");
+ source->supernova.z = parser_get_param_double(parameter_file, "SN:z");
#endif /* SN_FEEDBCK */
};
@@ -64,13 +64,18 @@ void source_terms_init(const struct swift_params* parameter_file,
void source_terms_print(const struct sourceterms* source) {
#ifdef SN_FEEDBACK
- message(" Single SNe of energy= %e will explode at time= %e at location (%e,%e,%e)",
- source->supernova.energy, source->supernova.time, source->supernova.x, source->supernova.y, source->supernova.z);
+ message(
+ " Single SNe of energy= %e will explode at time= %e at location "
+ "(%e,%e,%e)",
+ source->supernova.energy, source->supernova.time, source->supernova.x,
+ source->supernova.y, source->supernova.z);
#endif /* SN_FEEDBACK */
-
};
-__attribute__((always_inline)) INLINE float calculate_entropy(const float entropy_old, const float density, const float u_old, const float u_new)
-{
- return entropy_old + hydro_gamma_minus_one * (u_new-u_old) * pow_minus_gamma_minus_one(density);
-};
+/* __attribute__((always_inline)) INLINE float calculate_entropy( */
+/* const float entropy_old, const float density, const float u_old, */
+/* const float u_new) { */
+/* return entropy_old + */
+/* hydro_gamma_minus_one * (u_new - u_old) * */
+/* pow_minus_gamma_minus_one(density); */
+/* }; */
diff --git a/src/sourceterms.h b/src/sourceterms.h
index 506aadde82e66d9bc67c5bfe5bccb020d1dc7358..10166a0bd612010a9a46b7b25100ab27e2a8ba88 100644
--- a/src/sourceterms.h
+++ b/src/sourceterms.h
@@ -21,18 +21,17 @@
#include "./const.h"
/* So far only one model here */
-struct sourceterms
-{
+struct sourceterms {
#ifdef SN_FEEDBACK
#include "sourceterms/sn_feedback/sn_feedback_struct.h"
#endif
};
void source_terms_init(const struct swift_params* parameter_file,
- struct UnitSystem* us,
- struct sourceterms* source);
+ struct UnitSystem* us, struct sourceterms* source);
void source_terms_print(const struct sourceterms* source);
-float calculate_entropy(const float entropy_old, const float density, const float u_old, const float u_new);
+/* float calculate_entropy(const float entropy_old, const float density, */
+/* const float u_old, const float u_new); */
#ifdef SN_FEEDBACK
#include "sourceterms/sn_feedback/sn_feedback.h"
#endif
diff --git a/src/sourceterms/Default/sourceterms.h b/src/sourceterms/Default/sourceterms.h
index a31a675fbe36629fe79849691f8d0377626cb5b5..ea860f9cd17fda2c9d7767d6f88618f1090c9c8c 100644
--- a/src/sourceterms/Default/sourceterms.h
+++ b/src/sourceterms/Default/sourceterms.h
@@ -22,7 +22,7 @@
#include "feedback.h"
/**
- * @brief Computes the feedback time-step of a given particle due to
+ * @brief Computes the feedback time-step of a given particle due to
* a source term
*
* This function only branches towards the potential chosen by the user.
@@ -31,21 +31,18 @@
* @param phys_const The physical constants in internal units.
* @param g Pointer to the particle data.
*/
-__attribute__((always_inline)) INLINE static float
-sourceterms_compute_timestep(const struct sourceterms *feedback,
- const struct phys_const* const phys_const,
- const struct part* const p){
+__attribute__((always_inline)) INLINE static float sourceterms_compute_timestep(
+ const struct sourceterms* feedback,
+ const struct phys_const* const phys_const, const struct part* const p) {
float dt = FLT_MAX;
#ifdef SN_FEEDBACK
- dt =
- fmin(dt, sn_feedback_timestep(feedback, phys_const, p));
+ dt = fmin(dt, sn_feedback_timestep(feedback, phys_const, p));
#endif
}
-__attribute__((always_inline)) INLINE static void feedback(const struct sourceterms *feedback,
- const struct phys_const* const phys_const,
- struct part* p)
+__attribute__((always_inline)) INLINE static void feedback(
+ const struct sourceterms* feedback,
+ const struct phys_const* const phys_const, struct part* p)
#ifdef SN_FEEDBACK
- sn_feedback(feedback, phys_const, p);
+ sn_feedback(feedback, phys_const, p);
#endif
}
-
diff --git a/src/sourceterms/sn_feedback/sn_feedback.c b/src/sourceterms/sn_feedback/sn_feedback.c
index ebc9f69c55eb39bf300cef7b10ae157d97da97cb..6e6d3bade4ea7197355946e8b67e66016a7cdbca 100644
--- a/src/sourceterms/sn_feedback/sn_feedback.c
+++ b/src/sourceterms/sn_feedback/sn_feedback.c
@@ -1,4 +1,4 @@
-#error: this file is no longer in use!
+#error : this file is no longer in use!
/*******************************************************************************
* This file is part of SWIFT.
* Copyright (c) 2016 Tom Theuns (tom.theuns@durham.ac.uk)
@@ -26,26 +26,26 @@
#include "../sourceterms.h"
void supernova_feedback_init(const struct swift_params* parameter_file,
- struct UnitSystem* us,
- struct supernova* sn){
- sn.time = parser_get_param_double(parameter_file, "SN:time");
+ struct UnitSystem* us, struct supernova* sn) {
+ sn.time = parser_get_param_double(parameter_file, "SN:time");
sn.energy = parser_get_param_double(parameter_file, "SN:energy");
- sn.x = parser_get_param_double(parameter_file, "SN:x");
- sn.y = parser_get_param_double(parameter_file, "SN:y");
- sn.z = parser_get_param_double(parameter_file, "SN:z");
+ sn.x = parser_get_param_double(parameter_file, "SN:x");
+ sn.y = parser_get_param_double(parameter_file, "SN:y");
+ sn.z = parser_get_param_double(parameter_file, "SN:z");
}
-void supernova_feedback_print(const struct supernova* sn){
- message(" Single SNe of energy= %e will explode at time= %e at location (%e,%e,%e)",
- sn.energy, sn.time, sn.x, sn.y, sn.z);
+void supernova_feedback_print(const struct supernova* sn) {
+ message(
+ " Single SNe of energy= %e will explode at time= %e at location "
+ "(%e,%e,%e)",
+ sn.energy, sn.time, sn.x, sn.y, sn.z);
};
-__attribute__((always_inline)) INLINE static void do_supernova_feedback(const struct sourceterms* sourceterms, struct part* p)
-{
-};
+__attribute__((always_inline)) INLINE static void do_supernova_feedback(
+ const struct sourceterms* sourceterms, struct part* p){};
-__attribute__((always_inline)) INLINE void update_entropy(const sourceterms *sourceterms, struct part* p)
-{
+__attribute__((always_inline)) INLINE void update_entropy(
+ const sourceterms* sourceterms, struct part* p) {
/*updates the entropy of a particle due to feedback */
float u_old;
@@ -55,10 +55,10 @@ __attribute__((always_inline)) INLINE void update_entropy(const sourceterms *sou
float rho = p->rho;
// u_old = old_entropy/(GAMMA_MINUS1) * pow(rho,GAMMA_MINUS1);
- const float u_old = hydro_get_internal_energy(p,0); // dt = 0 because using current entropy
+ const float u_old =
+ hydro_get_internal_energy(p, 0); // dt = 0 because using current entropy
const float u_new = u_old + sourceterms->supernova.energy;
- const float new_entropy = u_new*pow_minus_gamma_minus_one(-p>rho) * hydro_gamma_minus_one;
+ const float new_entropy =
+ u_new * pow_minus_gamma_minus_one(-p > rho) * hydro_gamma_minus_one;
p->entropy = new_entropy;
}
-
-
diff --git a/src/sourceterms/sn_feedback/sn_feedback.h b/src/sourceterms/sn_feedback/sn_feedback.h
index 16a60c7b6f7bc10af00755281d660df0233e3696..2cb05d5a43feb5f474ba9ed70c5d6d398a36e3c3 100644
--- a/src/sourceterms/sn_feedback/sn_feedback.h
+++ b/src/sourceterms/sn_feedback/sn_feedback.h
@@ -24,36 +24,25 @@
#include <float.h>
#include "equation_of_state.h"
#include "hydro.h"
-// void do_supernova_feedback(const struct source_terms* source_terms, const struct phys_const* const constants, const struct part* p);
-
+// void do_supernova_feedback(const struct source_terms* source_terms, const
+// struct phys_const* const constants, const struct part* p);
/* determine whether location is in cell */
-__attribute__((always_inline)) INLINE static int is_in_cell(const double cell_min[], const double cell_width[], const double location[], const int dimen)
-{
- for(int i=0; i<dimen; i++){
- if (cell_min[i] > location[i])
- return 0;
- if ( (cell_min[i] + cell_width[i]) <= location[i])
- return 0;
+__attribute__((always_inline)) INLINE static int is_in_cell(
+ const double cell_min[], const double cell_width[], const double location[],
+ const int dimen) {
+ for (int i = 0; i < dimen; i++) {
+ if (cell_min[i] > location[i]) return 0;
+ if ((cell_min[i] + cell_width[i]) <= location[i]) return 0;
}
return 1;
};
-__attribute__((always_inline)) INLINE static void do_supernova_feedback(const struct sourceterms* sourceterms, struct part* p)
-{
- const float density = p->rho;
- const float u_old = hydro_get_internal_energy(p,0);
- const float u_new = u_old + sourceterms->supernova.energy / p->mass;
- const float entropy_old = gas_entropy_from_internal_energy(density, u_old);
- const float entropy_new = calculate_entropy(entropy_old, density, u_old, u_new);
+__attribute__((always_inline)) INLINE static void do_supernova_feedback(
+ const struct sourceterms* sourceterms, struct part* p) {
+ const float u_old = hydro_get_internal_energy(p, 0);
+ const float u_new = u_old + sourceterms->supernova.energy / hydro_get_mass(p);
hydro_set_internal_energy(p, u_new);
- hydro_set_entropy(p,entropy_new);
- message(" u_old= %e u_new= %e S_old= %e S_new= %e rho= %e", u_old, u_new, entropy_old, entropy_new, density);
- message(" p->s %e", p->entropy);
+
};
#endif /* SWIFT_SN_FEEDBACK_H */
-
-
-
-
-
diff --git a/src/sourceterms/sn_feedback/sn_feedback_struct.h b/src/sourceterms/sn_feedback/sn_feedback_struct.h
index fe662e49c77f6870e20a9dc8b34815b4f2ac1f8d..9531477b30e5bbb58885d92dfef322cea27539d2 100644
--- a/src/sourceterms/sn_feedback/sn_feedback_struct.h
+++ b/src/sourceterms/sn_feedback/sn_feedback_struct.h
@@ -1,7 +1,6 @@
/* sn feedback */
struct {
- double time;
- double energy;
- double x, y, z;
- } supernova;
-
+ double time;
+ double energy;
+ double x, y, z;
+} supernova;
diff --git a/src/swift.h b/src/swift.h
index 6d8d57968f1fd8283044f2dbb753d2943f0d43eb..fc6168cdda7cf47658d4f9f4616b717441b297a7 100644
--- a/src/swift.h
+++ b/src/swift.h
@@ -43,12 +43,12 @@
#include "partition.h"
#include "physical_constants.h"
#include "potentials.h"
-#include "sourceterms.h"
#include "queue.h"
#include "runner.h"
#include "scheduler.h"
#include "serial_io.h"
#include "single_io.h"
+#include "sourceterms.h"
#include "space.h"
#include "task.h"
#include "timers.h"
diff --git a/src/task.c b/src/task.c
index 54b9ec5fae4e5d6040eb0f65194e1b4280a0f1a1..1aa02893b151bdc3b9a7c472fda2004bf9b83109 100644
--- a/src/task.c
+++ b/src/task.c
@@ -48,11 +48,10 @@
/* Task type names. */
const char *taskID_names[task_type_count] = {
- "none", "sort", "self", "pair", "sub_self",
- "sub_pair", "init", "ghost", "drift", "kick",
- "kick_fixdt", "send", "recv", "grav_gather_m", "grav_fft",
- "grav_mm", "grav_up", "grav_external", "part_sort", "gpart_sort",
- "split_cell", "rewait"};
+ "none", "sort", "self", "pair", "sub_self",
+ "sub_pair", "init", "ghost", "extra_ghost", "kick",
+ "kick_fixdt", "send", "recv", "grav_gather_m", "grav_fft",
+ "grav_mm", "grav_up", "grav_external", "sourceterms"};
const char *subtaskID_names[task_subtype_count] = {"none", "density", "force",
"grav", "tend"};
@@ -111,6 +110,7 @@ __attribute__((always_inline)) INLINE static enum task_actions task_acts_on(
case task_type_sort:
case task_type_ghost:
+ case task_type_extra_ghost:
case task_type_sourceterms:
return task_action_part;
break;
diff --git a/src/task.h b/src/task.h
index a5cf89c9da62a20e3f0e356e81beee26d756255b..4c98e9ab9903ff6746acfc883e7e43223dd2e11a 100644
--- a/src/task.h
+++ b/src/task.h
@@ -52,10 +52,6 @@ enum task_types {
task_type_grav_up,
task_type_grav_external,
task_type_sourceterms,
- task_type_part_sort,
- task_type_gpart_sort,
- task_type_split_cell,
- task_type_rewait,
task_type_count
};