diff --git a/examples/CoolingBox/energy_plot.py b/examples/CoolingBox/energy_plot.py index 9bf055845a61bd204d94e55e8ad6cb761910e8ca..0533528d8aabfbf11e18da95911a5627db820b04 100644 --- a/examples/CoolingBox/energy_plot.py +++ b/examples/CoolingBox/energy_plot.py @@ -11,9 +11,9 @@ snap_filename = "uniformBox_000.hdf5" k_b = 1.38E-16 #boltzmann m_p = 1.67e-24 #proton mass #initial conditions set in makeIC.py -rho = 3.2e6 -P = 4.5e9 -n_H_cgs = 0.1 +rho = 3.2e3 +P = 4.5e6 +n_H_cgs = 0.0001 gamma = 5./3. T_init = 1.0e5 @@ -24,10 +24,10 @@ unit_mass = units.attrs["Unit mass in cgs (U_M)"] unit_length = units.attrs["Unit length in cgs (U_L)"] unit_time = units.attrs["Unit time in cgs (U_t)"] parameters = f["Parameters"] -cooling_lambda = float(parameters.attrs["CreaseyCooling:Lambda"]) -min_T = float(parameters.attrs["CreaseyCooling:minimum_temperature"]) -mu = float(parameters.attrs["CreaseyCooling:mean_molecular_weight"]) -X_H = float(parameters.attrs["CreaseyCooling:hydrogen_mass_abundance"]) +cooling_lambda = float(parameters.attrs["Cooling:lambda"]) +min_T = float(parameters.attrs["Cooling:minimum_temperature"]) +mu = float(parameters.attrs["Cooling:mean_molecular_weight"]) +X_H = float(parameters.attrs["Cooling:hydrogen_mass_abundance"]) #get number of particles header = f["Header"] @@ -38,23 +38,29 @@ time = array[:,0] total_energy = array[:,2] total_mass = array[:,1] +time = time[1:] +total_energy = total_energy[1:] +total_mass = total_mass[1:] + #conversions to cgs rho_cgs = rho * unit_mass / (unit_length)**3 time_cgs = time * unit_time u_init_cgs = total_energy[0]/(total_mass[0]) * unit_length**2 / (unit_time)**2 #find the energy floor +print min_T u_floor_cgs = k_b * min_T / (mu * m_p * (gamma - 1.)) #find analytic solution -analytic_time = np.linspace(0.0,time_cgs[-1],10000) +analytic_time = np.linspace(time_cgs[0],time_cgs[-1],1000) +print time_cgs[1] +print analytic_time[1] du_dt_cgs = -cooling_lambda * n_H_cgs**2 / rho_cgs -u_analytic = du_dt_cgs*analytic_time + u_init_cgs +u_analytic = du_dt_cgs*(analytic_time - analytic_time[0]) + u_init_cgs cooling_time = u_init_cgs/(-du_dt_cgs) #rescale energy to initial energy total_energy /= total_energy[0] u_analytic /= u_init_cgs u_floor_cgs /= u_init_cgs - # plot_title = r"$\Lambda \, = \, %1.1g \mathrm{erg}\mathrm{cm^3}\mathrm{s^{-1}} \, \, T_{init} = %1.1g\mathrm{K} \, \, T_{floor} = %1.1g\mathrm{K} \, \, n_H = %1.1g\mathrm{cm^{-3}}$" %(cooling_lambda,T_init,T_floor,n_H) # plot_filename = "energy_plot_creasey_no_cooling_T_init_1p0e5_n_H_0p1.png" #analytic_solution = np.zeros(n_snaps-1) @@ -64,11 +70,11 @@ for i in range(u_analytic.size): plt.plot(time_cgs,total_energy,'k',label = "Numerical solution") plt.plot(analytic_time,u_analytic,'--r',lw = 2.0,label = "Analytic Solution") plt.plot((cooling_time,cooling_time),(0,1),'b',label = "Cooling time") -plt.plot((time_cgs[1],time_cgs[1]),(0,1),'m',label = "First output") +plt.plot((time_cgs[0],time_cgs[0]),(0,1),'m',label = "First output") plt.title(r"$n_H = %1.1e \, \mathrm{cm}^{-3}$" %n_H_cgs) plt.xlabel("Time (seconds)") plt.ylabel("Energy/Initial energy") -plt.ylim(0,1) +plt.ylim(0.999,1.001) plt.xlim(0,min(10.0*cooling_time,time_cgs[-1])) plt.legend(loc = "upper right") if (int(sys.argv[1])==0): diff --git a/examples/CoolingBox/makeIC.py b/examples/CoolingBox/makeIC.py index e8e2e1f506bbc581fd87b5a3ef3de40cae42354c..e1bfc28f063ff1d34dbdcbadce33df468ccce4ce 100644 --- a/examples/CoolingBox/makeIC.py +++ b/examples/CoolingBox/makeIC.py @@ -29,8 +29,8 @@ from numpy import * periodic= 1 # 1 For periodic box boxSize = 1 #1 kiloparsec L = int(sys.argv[1]) # Number of particles along one axis -rho = 3.2e6 # Density in code units (0.1 hydrogen atoms per cm^3) -P = 4.5e9 # Pressure in code units (at 10^5K) +rho = 3.2e3 # Density in code units (0.01 hydrogen atoms per cm^3) +P = 4.5e6 # Pressure in code units (at 10^5K) gamma = 5./3. # Gas adiabatic index eta = 1.2349 # 48 ngbs with cubic spline kernel fileName = "uniformBox.hdf5" diff --git a/examples/CoolingBox/run.sh b/examples/CoolingBox/run.sh index 87bac7a3610d5cc5ebce009be5a4023d115e13f2..a5fd5e41aba67a3cfb19004e1c523e809616f300 100755 --- a/examples/CoolingBox/run.sh +++ b/examples/CoolingBox/run.sh @@ -5,6 +5,8 @@ echo "Generating initial conditions for the uniform box example..." python makeIC.py 10 -../swift -s -C -t 16 uniformBox.yml +../swift -s -t 16 uniformBox.yml python energy_plot.py 0 + +python add_energy_column diff --git a/examples/CoolingBox/uniformBox.yml b/examples/CoolingBox/uniformBox.yml index 995efbe14d4a05114085bc99ace72b9560ac02f7..61713f92ed2f37759570dfcefabd8c69ab0907eb 100644 --- a/examples/CoolingBox/uniformBox.yml +++ b/examples/CoolingBox/uniformBox.yml @@ -43,9 +43,8 @@ PointMass: mass: 1e10 # mass of external point mass in internal units # Cooling parameters (Creasey cooling) (always in cgs) - -CreaseyCooling: - Lambda: 1.0e-22 +Cooling: + lambda: 0.0 minimum_temperature: 1.0e4 mean_molecular_weight: 0.59 hydrogen_mass_abundance: 0.75 diff --git a/src/Makefile.am b/src/Makefile.am index df1491f32e70bf298ea2b18539f6f16e16ed6871..0a4b2c64e85b30d9d15b9216f6135517e4217cc5 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -72,7 +72,7 @@ nobase_noinst_HEADERS = approx_math.h atomic.h cycle.h error.h inline.h kernel_h hydro/Gizmo/hydro_debug.h hydro/Gizmo/hydro_part.h \ riemann.h \ riemann/riemann_hllc.h riemann/riemann_trrs.h riemann/riemann_exact.h \ - cooling/const/cooling.h + cooling/const_du/cooling.h cooling/const_lambda/cooling.h # Sources and flags for regular library libswiftsim_la_SOURCES = $(AM_SOURCES) diff --git a/src/const.h b/src/const.h index e83d71d8733bdc43d1d8735de0d3cf0112019a61..cbd814f4a9646e0709b6c4ebd9b4ac3c19be2f5b 100644 --- a/src/const.h +++ b/src/const.h @@ -95,8 +95,8 @@ //#define EXTERNAL_POTENTIAL_DISK_PATCH /* Cooling properties */ -#define COOLING_CONST_COOLING -//#define COOLING_CREASEY_COOLING +//#define COOLING_CONST_DU +#define COOLING_CONST_LAMBDA //#define COOLING_GRACKLE_COOLING /* Are we debugging ? */ diff --git a/src/cooling.h b/src/cooling.h index 4c4694c86b07e587aa4f6e9a48ff046d0b180f34..691b2fd7008a6cdd8c6688bd5fe644d5b48ac897 100644 --- a/src/cooling.h +++ b/src/cooling.h @@ -31,10 +31,10 @@ #include "const.h" /* Import the right cooling definition */ -#if defined(COOLING_CONST_COOLING) -#include "./cooling/const/cooling.h" -#elif defined(COOLING_CREASEY_COOLING) -#include "./cooling/creasey/cooling.h" +#if defined(COOLING_CONST_DU) +#include "./cooling/const_du/cooling.h" +#elif defined(COOLING_CONST_LAMBDA) +#include "./cooling/const_lambda/cooling.h" #elif defined(COOLING_GRACKLE_COOLING) #include "./cooling/grackle/cooling.h" #else diff --git a/src/cooling/const/cooling.h b/src/cooling/const/cooling.h deleted file mode 100644 index 41b1f8bd5d1c8d99aa0dfbadf64850a02db5df88..0000000000000000000000000000000000000000 --- a/src/cooling/const/cooling.h +++ /dev/null @@ -1,139 +0,0 @@ -/******************************************************************************* - * This file is part of SWIFT. - * Copyright (c) 2016 Matthieu Schaller (matthieu.schaller@durham.ac.uk) - * Richard Bower (r.g.bower@durham.ac.uk) - * Stefan Arridge (stefan.arridge@durham.ac.uk) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU Lesser General Public License as published - * by the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU Lesser General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - * - ******************************************************************************/ -#ifndef SWIFT_CONST_COOLING_H -#define SWIFT_CONST_COOLING_H - -/** - * @file src/cooling/const/cooling.h - * @brief Routines related to the "constant cooling" cooling function. - * - * This is the simplest possible cooling function. A constant cooling rate with - * a minimal energy floor is applied. - */ - -/* Some standard headers. */ -#include <math.h> - -/* Local includes. */ -#include "const.h" -#include "error.h" -#include "hydro.h" -#include "parser.h" -#include "part.h" -#include "physical_constants.h" -#include "units.h" - -/** - * @brief Properties of the cooling function. - */ -struct cooling_data { - - /*! Cooling rate in internal units */ - float lambda; - - /*! Minimally allowed internal energy of the particles */ - float min_energy; - - /*! Constant multiplication factor for time-step criterion */ - float cooling_tstep_mult; -}; - -/** - * @brief Apply the cooling function to a particle. - * - * @param phys_const The physical constants in internal units. - * @param us The internal system of units. - * @param cooling The #cooling_data used in the run. - * @param p Pointer to the particle data. - * @param dt The time-step of this particle. - */ -__attribute__((always_inline)) INLINE static void cooling_cool_part( - const struct phys_const* const phys_const, const struct UnitSystem* us, - const struct cooling_data* cooling, struct part* p, float dt) { - - /* Get current internal energy (dt=0) */ - const float u_old = hydro_get_internal_energy(p, 0.f); - - /* Get cooling function properties */ - const float du_dt = -cooling->lambda; - const float u_floor = cooling->min_energy; - - /* Constant cooling with a minimal floor */ - float u_new; - if (u_old - du_dt * dt > u_floor) { - u_new = u_old + du_dt * dt; - } else { - u_new = u_floor; - } - - /* Update the internal energy */ - hydro_set_internal_energy(p, u_new); -} - -/** - * @brief Computes the cooling time-step. - * - * @param cooling The #cooling_data used in the run. - * @param phys_const The physical constants in internal units. - * @param p Pointer to the particle data. - */ -__attribute__((always_inline)) INLINE static double cooling_timestep( - const struct cooling_data* cooling, - const struct phys_const* const phys_const, const struct part* const p) { - - const float cooling_rate = cooling->lambda; - const float internal_energy = - hydro_get_internal_energy(p, 0); // dt = 0 because using current entropy - return cooling->cooling_tstep_mult * internal_energy / cooling_rate; -} - -/** - * @brief Initialises the cooling properties. - * - * @param parameter_file The parsed parameter file. - * @param us The current internal system of units. - * @param phys_const The physical constants in internal units. - * @param cooling The cooling properties to initialize - */ -INLINE void cooling_init(const struct swift_params* parameter_file, - const struct UnitSystem* us, - const struct phys_const* phys_const, - struct cooling_data* cooling) { - - cooling->lambda = parser_get_param_double(parameter_file, "Cooling:lambda"); - cooling->min_energy = - parser_get_param_double(parameter_file, "Cooling:min_energy"); - cooling->cooling_tstep_mult = - parser_get_param_double(parameter_file, "Cooling:cooling_tstep_mult"); -} - -/** - * @brief Prints the properties of the cooling model to stdout. - * - * @param cooling The properties of the cooling function. - */ -INLINE void cooling_print(const struct cooling_data* cooling) { - - message("Cooling function is 'Constant cooling' with rate %f and floor %f", - cooling->lambda, cooling->min_energy); -} - -#endif /* SWIFT_CONST_COOLING_H */ diff --git a/src/cooling/creasey/cooling.c b/src/cooling/creasey/cooling.c deleted file mode 100644 index 4ee03aa9bcbe0fb616fe90c1844da98a3c1ccf79..0000000000000000000000000000000000000000 --- a/src/cooling/creasey/cooling.c +++ /dev/null @@ -1,175 +0,0 @@ - -/******************************************************************************* - * This file is part of SWIFT. - * Copyright (c) 2016 Tom Theuns (tom.theuns@durham.ac.uk) - * Matthieu Schaller (matthieu.schaller@durham.ac.uk) - * Richard Bower (r.g.bower@durham.ac.uk) - * Stefan Arridge (stefan.arridge@durham.ac.uk) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU Lesser General Public License as published - * by the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU Lesser General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - * - ******************************************************************************/ - -/* Config parameters. */ -#include "../config.h" - -/* This object's header. */ -#include "adiabatic_index.h" -#include "cooling.h" -#include "hydro.h" - -/** - * @brief Initialises the cooling properties in the internal system - * of units. - * - * @param parameter_file The parsed parameter file - * @param us The current internal system of units - * @param cooling The cooling properties to initialize - */ -void cooling_init(const struct swift_params* parameter_file, - struct UnitSystem* us, - const struct phys_const* const phys_const, - struct cooling_data* cooling) { - -#ifdef CONST_COOLING - cooling->const_cooling.lambda = - parser_get_param_double(parameter_file, "Cooling:lambda"); - cooling->const_cooling.min_energy = - parser_get_param_double(parameter_file, "Cooling:min_energy"); - cooling->const_cooling.cooling_tstep_mult = - parser_get_param_double(parameter_file, "Cooling:cooling_tstep_mult"); -#endif /* CONST_COOLING */ - -#ifdef CREASEY_COOLING - cooling->creasey_cooling.lambda = - parser_get_param_double(parameter_file, "CreaseyCooling:Lambda"); - cooling->creasey_cooling.min_temperature = parser_get_param_double( - parameter_file, "CreaseyCooling:minimum_temperature"); - cooling->creasey_cooling.mean_molecular_weight = parser_get_param_double( - parameter_file, "CreaseyCooling:mean_molecular_weight"); - cooling->creasey_cooling.hydrogen_mass_abundance = parser_get_param_double( - parameter_file, "CreaseyCooling:hydrogen_mass_abundance"); - cooling->creasey_cooling.cooling_tstep_mult = parser_get_param_double( - parameter_file, "CreaseyCooling:cooling_tstep_mult"); - - /*convert minimum temperature into minimum internal energy*/ - float u_floor = - phys_const->const_boltzmann_k * cooling->creasey_cooling.min_temperature / - (hydro_gamma_minus_one * cooling->creasey_cooling.mean_molecular_weight * - phys_const->const_proton_mass); - float u_floor_cgs = - u_floor * units_cgs_conversion_factor(us, UNIT_CONV_ENERGY_PER_UNIT_MASS); - - cooling->creasey_cooling.min_internal_energy = u_floor; - cooling->creasey_cooling.min_internal_energy_cgs = u_floor_cgs; -#endif /* CREASEY_COOLING */ -} - -/** - * @brief Prints the properties of the cooling model to stdout. - * - * @param cooling The cooling properties. - */ -void cooling_print(const struct cooling_data* cooling) { - -#ifdef CONST_COOLING - message( - "Cooling properties are (lambda, min_energy, tstep multiplier) %g %g %g ", - cooling->const_cooling.lambda, cooling->const_cooling.min_energy, - cooling->const_cooling.cooling_tstep_mult); -#endif /* CONST_COOLING */ - -#ifdef CREASEY_COOLING - message( - "Cooling properties for Creasey cooling are (lambda, min_temperature, " - "hydrogen_mass_abundance, mean_molecular_weight, tstep multiplier) %g %g " - "%g %g %g", - cooling->creasey_cooling.lambda, cooling->creasey_cooling.min_temperature, - cooling->creasey_cooling.hydrogen_mass_abundance, - cooling->creasey_cooling.mean_molecular_weight, - cooling->creasey_cooling.cooling_tstep_mult); -#endif /* CREASEY_COOLING */ -} - -void update_entropy(const struct phys_const* const phys_const, - const struct UnitSystem* us, - const struct cooling_data* cooling, struct part* p, - float dt) { - - /*updates the entropy of a particle after integrating the cooling equation*/ - float u_old; - float u_new; - float new_entropy; - // float old_entropy = p->entropy; - float rho = p->rho; - - // u_old = old_entropy/(GAMMA_MINUS1) * pow(rho,GAMMA_MINUS1); - u_old = - hydro_get_internal_energy(p, 0); // dt = 0 because using current entropy - u_new = calculate_new_thermal_energy(u_old, rho, dt, cooling, phys_const, us); - new_entropy = u_new * pow_minus_gamma_minus_one(rho) * hydro_gamma_minus_one; - p->entropy = new_entropy; -} - -/*This function integrates the cooling equation, given the initial - thermal energy, density and the timestep dt. Returns the final internal - energy*/ - -float calculate_new_thermal_energy(float u_old, float rho, float dt, - const struct cooling_data* cooling, - const struct phys_const* const phys_const, - const struct UnitSystem* us) { -#ifdef CONST_COOLING - /*du/dt = -lambda, independent of density*/ - float du_dt = -cooling->const_cooling.lambda; - float u_floor = cooling->const_cooling.min_energy; - float u_new; - if (u_old - du_dt * dt > u_floor) { - u_new = u_old + du_dt * dt; - } else { - u_new = u_floor; - } -#endif /*CONST_COOLING*/ - -#ifdef CREASEY_COOLING - /* rho*du/dt = -lambda*n_H^2 */ - float u_new; - float X_H = cooling->creasey_cooling.hydrogen_mass_abundance; - float lambda_cgs = cooling->creasey_cooling.lambda; // this is always in cgs - float u_floor_cgs = cooling->creasey_cooling.min_internal_energy_cgs; - - /*convert from internal code units to cgs*/ - float dt_cgs = dt * units_cgs_conversion_factor(us, UNIT_CONV_TIME); - float rho_cgs = rho * units_cgs_conversion_factor(us, UNIT_CONV_DENSITY); - float m_p_cgs = phys_const->const_proton_mass * - units_cgs_conversion_factor(us, UNIT_CONV_MASS); - float n_H_cgs = X_H * rho_cgs / m_p_cgs; - float u_old_cgs = - u_old * units_cgs_conversion_factor(us, UNIT_CONV_ENERGY_PER_UNIT_MASS); - float du_dt_cgs = -lambda_cgs * n_H_cgs * n_H_cgs / rho_cgs; - float u_new_cgs; - - if (u_old_cgs + du_dt_cgs * dt_cgs > u_floor_cgs) { - u_new_cgs = u_old_cgs + du_dt_cgs * dt_cgs; - } else { - u_new_cgs = u_floor_cgs; - } - /*convert back to internal code units when returning new internal energy*/ - - u_new = u_new_cgs / - units_cgs_conversion_factor(us, UNIT_CONV_ENERGY_PER_UNIT_MASS); - -#endif /*CREASEY_COOLING*/ - return u_new; -} diff --git a/src/cooling/creasey/cooling.h b/src/cooling/creasey/cooling.h deleted file mode 100644 index 62ecbd75d92d687c1f56124d6c0802ef6a1b6851..0000000000000000000000000000000000000000 --- a/src/cooling/creasey/cooling.h +++ /dev/null @@ -1,100 +0,0 @@ -/******************************************************************************* - * This file is part of SWIFT. - * Copyright (c) 2016 Tom Theuns (tom.theuns@durham.ac.uk) - * Matthieu Schaller (matthieu.schaller@durham.ac.uk) - * Richard Bower (r.g.bower@durham.ac.uk) - * Stefan Arridge (stefan.arridge@durham.ac.uk) - * - * This program is free software: you can redistribute it and/or modify - * it under the terms of the GNU Lesser General Public License as published - * by the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * This program is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU Lesser General Public License - * along with this program. If not, see <http://www.gnu.org/licenses/>. - * - ******************************************************************************/ - -#ifndef SWIFT_CREASEY_COOLING_H -#define SWIFT_CREASEY_COOLING_H - -/* Some standard headers. */ -#include <math.h> - -/* Local includes. */ -#include "const.h" -#include "error.h" -#include "hydro.h" -#include "parser.h" -#include "part.h" -#include "physical_constants.h" -#include "units.h" - -/* Cooling Properties */ -struct cooling_data { - -#ifdef CONST_COOLING - struct { - float lambda; - float min_energy; - float cooling_tstep_mult; - } const_cooling; -#endif - -#ifdef CREASEY_COOLING - struct { - float lambda; - float min_temperature; - float hydrogen_mass_abundance; - float mean_molecular_weight; - float min_internal_energy; - float min_internal_energy_cgs; - float cooling_tstep_mult; - } creasey_cooling; -#endif -}; - -/* Include Cooling */ -#ifdef CONST_COOLING - -/** - * @brief Computes the time-step due to cooling - * - * @param cooling The #cooling_data used in the run. - * @param phys_const The physical constants in internal units. - * @param Pointer to the particle data. - */ -__attribute__((always_inline)) INLINE static double cooling_timestep( - const struct cooling_data* cooling, - const struct phys_const* const phys_const, const struct part* const p) { - - const double cooling_rate = cooling->const_cooling.lambda; - const double internal_energy = - hydro_get_internal_energy(p, 0); // dt = 0 because using current entropy - return cooling->const_cooling.cooling_tstep_mult * internal_energy / - cooling_rate; -} - -#endif /* CONST_COOLING */ - -/* Now, some generic functions, defined in the source file */ -void cooling_init(const struct swift_params* parameter_file, - struct UnitSystem* us, - const struct phys_const* const phys_const, - struct cooling_data* cooling); - -void cooling_print(const struct cooling_data* cooling); -void update_entropy(const struct phys_const* const phys_const, - const struct UnitSystem* us, - const struct cooling_data* cooling, struct part* p, - float dt); -float calculate_new_thermal_energy(float u_old, float rho, float dt, - const struct cooling_data* cooling, - const struct phys_const* const phys_const, - const struct UnitSystem* us); -#endif /* SWIFT_CREASEY_COOLING_H */ diff --git a/src/runner.c b/src/runner.c index 004ec2e316fbe1e08a1df7a14d3b979df5d7387f..d3fc7bc9e04596bc2fb848b14c1c1cde56cd3fe7 100644 --- a/src/runner.c +++ b/src/runner.c @@ -1319,6 +1319,9 @@ void *runner_main(void *data) { case task_type_grav_external: runner_do_grav_external(r, t->ci, 1); break; + case task_type_cooling: + runner_do_cooling(r, t->ci, 1); + break; default: error("Unknown task type."); }