diff --git a/src/adiabatic_index.h b/src/adiabatic_index.h
index f65f2dac13b9cf1a470ded155590cf5e443d0c76..cb3cabe357567eb9bc0a6f99dc8a7e1701eda094 100644
--- a/src/adiabatic_index.h
+++ b/src/adiabatic_index.h
@@ -108,7 +108,7 @@
  *
  * Computes \f$x^\gamma\f$.
  */
-__attribute__((always_inline)) INLINE static float pow_gamma(float x) {
+__attribute__((always_inline, const)) INLINE static float pow_gamma(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
 
@@ -141,7 +141,7 @@ __attribute__((always_inline)) INLINE static float pow_gamma(float x) {
  *
  * Computes \f$x^{(\gamma-1)}\f$.
  */
-__attribute__((always_inline)) INLINE static float pow_gamma_minus_one(
+__attribute__((always_inline, const)) INLINE static float pow_gamma_minus_one(
     float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
@@ -175,8 +175,8 @@ __attribute__((always_inline)) INLINE static float pow_gamma_minus_one(
  *
  * Computes \f$x^{-(\gamma-1)}\f$.
  */
-__attribute__((always_inline)) INLINE static float pow_minus_gamma_minus_one(
-    float x) {
+__attribute__((always_inline, const)) INLINE static float
+pow_minus_gamma_minus_one(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
 
@@ -212,7 +212,8 @@ __attribute__((always_inline)) INLINE static float pow_minus_gamma_minus_one(
  * @param x Argument
  * @return One over the argument to the power given by the adiabatic index
  */
-__attribute__((always_inline)) INLINE static float pow_minus_gamma(float x) {
+__attribute__((always_inline, const)) INLINE static float pow_minus_gamma(
+    float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
 
@@ -252,8 +253,8 @@ __attribute__((always_inline)) INLINE static float pow_minus_gamma(float x) {
  * @param x Argument
  * @return Argument to the power two divided by the adiabatic index minus one
  */
-__attribute__((always_inline)) INLINE static float pow_two_over_gamma_minus_one(
-    float x) {
+__attribute__((always_inline, const)) INLINE static float
+pow_two_over_gamma_minus_one(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
 
@@ -292,7 +293,7 @@ __attribute__((always_inline)) INLINE static float pow_two_over_gamma_minus_one(
  * @return Argument to the power two times the adiabatic index divided by the
  * adiabatic index minus one
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 pow_two_gamma_over_gamma_minus_one(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
@@ -336,7 +337,7 @@ pow_two_gamma_over_gamma_minus_one(float x) {
  * @return Argument to the power the adiabatic index minus one divided by two
  * times the adiabatic index
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 pow_gamma_minus_one_over_two_gamma(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
@@ -373,7 +374,7 @@ pow_gamma_minus_one_over_two_gamma(float x) {
  * @return Inverse argument to the power the adiabatic index plus one divided by
  * two times the adiabatic index
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 pow_minus_gamma_plus_one_over_two_gamma(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
@@ -408,7 +409,8 @@ pow_minus_gamma_plus_one_over_two_gamma(float x) {
  * @param x Argument
  * @return Argument to the power one over the adiabatic index
  */
-__attribute__((always_inline)) INLINE static float pow_one_over_gamma(float x) {
+__attribute__((always_inline, const)) INLINE static float pow_one_over_gamma(
+    float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
 
@@ -441,8 +443,8 @@ __attribute__((always_inline)) INLINE static float pow_one_over_gamma(float x) {
  *
  * @param x Argument
  */
-__attribute__((always_inline)) INLINE static float pow_three_gamma_minus_two(
-    float x) {
+__attribute__((always_inline, const)) INLINE static float
+pow_three_gamma_minus_two(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
 
@@ -476,7 +478,7 @@ __attribute__((always_inline)) INLINE static float pow_three_gamma_minus_two(
  *
  * @param x Argument
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 pow_three_gamma_minus_five_over_two(float x) {
 
 #if defined(HYDRO_GAMMA_5_3)
diff --git a/src/equation_of_state/ideal_gas/equation_of_state.h b/src/equation_of_state/ideal_gas/equation_of_state.h
index 0d57f6a5ce51091f82e79b009b0e85f1368d51cf..4da5f69f29a2ef3443443fa25fe0388fea29e495 100644
--- a/src/equation_of_state/ideal_gas/equation_of_state.h
+++ b/src/equation_of_state/ideal_gas/equation_of_state.h
@@ -45,7 +45,7 @@ struct eos_parameters {};
  * @param density The density \f$\rho\f$.
  * @param entropy The entropy \f$S\f$.
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 gas_internal_energy_from_entropy(float density, float entropy) {
 
   return entropy * pow_gamma_minus_one(density) *
@@ -60,8 +60,8 @@ gas_internal_energy_from_entropy(float density, float entropy) {
  * @param density The density \f$\rho\f$.
  * @param entropy The entropy \f$S\f$.
  */
-__attribute__((always_inline)) INLINE static float gas_pressure_from_entropy(
-    float density, float entropy) {
+__attribute__((always_inline, const)) INLINE static float
+gas_pressure_from_entropy(float density, float entropy) {
 
   return entropy * pow_gamma(density);
 }
@@ -75,8 +75,8 @@ __attribute__((always_inline)) INLINE static float gas_pressure_from_entropy(
  * @param pressure The pressure \f$P\f$.
  * @return The entropy \f$A\f$.
  */
-__attribute__((always_inline)) INLINE static float gas_entropy_from_pressure(
-    float density, float pressure) {
+__attribute__((always_inline, const)) INLINE static float
+gas_entropy_from_pressure(float density, float pressure) {
 
   return pressure * pow_minus_gamma(density);
 }
@@ -89,8 +89,8 @@ __attribute__((always_inline)) INLINE static float gas_entropy_from_pressure(
  * @param density The density \f$\rho\f$.
  * @param entropy The entropy \f$S\f$.
  */
-__attribute__((always_inline)) INLINE static float gas_soundspeed_from_entropy(
-    float density, float entropy) {
+__attribute__((always_inline, const)) INLINE static float
+gas_soundspeed_from_entropy(float density, float entropy) {
 
   return sqrtf(hydro_gamma * pow_gamma_minus_one(density) * entropy);
 }
@@ -103,7 +103,7 @@ __attribute__((always_inline)) INLINE static float gas_soundspeed_from_entropy(
  * @param density The density \f$\rho\f$
  * @param u The internal energy \f$u\f$
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 gas_entropy_from_internal_energy(float density, float u) {
 
   return hydro_gamma_minus_one * u * pow_minus_gamma_minus_one(density);
@@ -117,7 +117,7 @@ gas_entropy_from_internal_energy(float density, float u) {
  * @param density The density \f$\rho\f$
  * @param u The internal energy \f$u\f$
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 gas_pressure_from_internal_energy(float density, float u) {
 
   return hydro_gamma_minus_one * u * density;
@@ -132,7 +132,7 @@ gas_pressure_from_internal_energy(float density, float u) {
  * @param pressure The pressure \f$P\f$.
  * @return The internal energy \f$u\f$.
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 gas_internal_energy_from_pressure(float density, float pressure) {
   return hydro_one_over_gamma_minus_one * pressure / density;
 }
@@ -145,7 +145,7 @@ gas_internal_energy_from_pressure(float density, float pressure) {
  * @param density The density \f$\rho\f$
  * @param u The internal energy \f$u\f$
  */
-__attribute__((always_inline)) INLINE static float
+__attribute__((always_inline, const)) INLINE static float
 gas_soundspeed_from_internal_energy(float density, float u) {
 
   return sqrtf(u * hydro_gamma * hydro_gamma_minus_one);
@@ -159,8 +159,8 @@ gas_soundspeed_from_internal_energy(float density, float u) {
  * @param density The density \f$\rho\f$
  * @param P The pressure \f$P\f$
  */
-__attribute__((always_inline)) INLINE static float gas_soundspeed_from_pressure(
-    float density, float P) {
+__attribute__((always_inline, const)) INLINE static float
+gas_soundspeed_from_pressure(float density, float P) {
 
   const float density_inv = 1.f / density;
   return sqrtf(hydro_gamma * P * density_inv);
@@ -176,16 +176,16 @@ __attribute__((always_inline)) INLINE static float gas_soundspeed_from_pressure(
  * @param us The internal unit system.
  * @param params The parsed parameters.
  */
-__attribute__((always_inline)) INLINE static void eos_init(
-    struct eos_parameters *e, const struct phys_const *phys_const,
-    const struct unit_system *us, struct swift_params *params) {}
+INLINE static void eos_init(struct eos_parameters *e,
+                            const struct phys_const *phys_const,
+                            const struct unit_system *us,
+                            struct swift_params *params) {}
 /**
  * @brief Print the equation of state
  *
  * @param e The #eos_parameters
  */
-__attribute__((always_inline)) INLINE static void eos_print(
-    const struct eos_parameters *e) {
+INLINE static void eos_print(const struct eos_parameters *e) {
 
   message("Equation of state: Ideal gas.");
 
@@ -199,8 +199,8 @@ __attribute__((always_inline)) INLINE static void eos_print(
  * @param h_grpsph The HDF5 group in which to write
  * @param e The #eos_parameters
  */
-__attribute__((always_inline)) INLINE static void eos_print_snapshot(
-    hid_t h_grpsph, const struct eos_parameters *e) {
+INLINE static void eos_print_snapshot(hid_t h_grpsph,
+                                      const struct eos_parameters *e) {
 
   io_write_attribute_f(h_grpsph, "Adiabatic index", hydro_gamma);
 
diff --git a/src/periodic.h b/src/periodic.h
index bd46d7c8ef82f6585ffbb892f8c0392caeccb760..05cb2fe454feadf8fe75942294443a726fc94d89 100644
--- a/src/periodic.h
+++ b/src/periodic.h
@@ -49,7 +49,7 @@
  * Similarly for dx < -b.
  *
  */
-__attribute__((always_inline)) INLINE static double nearest(
+__attribute__((always_inline, const)) INLINE static double nearest(
     const double dx, const double box_size) {
 
   return ((dx > 0.5 * box_size)
@@ -67,7 +67,7 @@ __attribute__((always_inline)) INLINE static double nearest(
  * Similarly for dx < -b.
  *
  */
-__attribute__((always_inline)) INLINE static float nearestf(
+__attribute__((always_inline, const)) INLINE static float nearestf(
     const float dx, const float box_size) {
 
   return ((dx > 0.5f * box_size)
diff --git a/src/sign.h b/src/sign.h
index 04ecdb3a99c5baf35e40d6f8286d25102d22e5fd..a92875d5d6c5f319eb6019adb20aa34af86d620b 100644
--- a/src/sign.h
+++ b/src/sign.h
@@ -25,7 +25,7 @@
  * @param x The number of interest.
  * @return >0 if positive, 0 if negative.
  */
-__attribute__((always_inline)) INLINE static int signf(float x) {
+__attribute__((always_inline, const)) INLINE static int signf(float x) {
 #ifdef __GNUC__
   return !signbit(x);
 #else
@@ -39,7 +39,8 @@ __attribute__((always_inline)) INLINE static int signf(float x) {
  * @param x The first number
  * @param y The second number
  */
-__attribute__((always_inline)) INLINE static int same_signf(float x, float y) {
+__attribute__((always_inline, const)) INLINE static int same_signf(float x,
+                                                                   float y) {
   return signf(x) == signf(y);
 }
 
diff --git a/src/timeline.h b/src/timeline.h
index 66b2c9456fb9e2e9f8b36272a9bbf3b7764523fc..4078a904c3e9205b8ea6ae7090534bd3d3d0784f 100644
--- a/src/timeline.h
+++ b/src/timeline.h
@@ -51,7 +51,8 @@ typedef char timebin_t;
  *
  * @param bin The time bin of interest.
  */
-static INLINE integertime_t get_integer_timestep(timebin_t bin) {
+__attribute__((const)) static INLINE integertime_t
+get_integer_timestep(timebin_t bin) {
 
   if (bin <= 0) return 0;
   return 1LL << (bin + 1);
@@ -62,7 +63,8 @@ static INLINE integertime_t get_integer_timestep(timebin_t bin) {
  *
  * Assumes that integertime_t maps to an unsigned long long.
  */
-static INLINE timebin_t get_time_bin(integertime_t time_step) {
+__attribute__((const)) static INLINE timebin_t
+get_time_bin(integertime_t time_step) {
 
   /* ((int) log_2(time_step)) - 1 */
   return (timebin_t)(62 - intrinsics_clzll(time_step));
@@ -74,7 +76,8 @@ static INLINE timebin_t get_time_bin(integertime_t time_step) {
  * @param bin The time bin of interest.
  * @param timeBase the minimal time-step size of the simulation.
  */
-static INLINE double get_timestep(timebin_t bin, double timeBase) {
+__attribute__((const)) static INLINE double get_timestep(timebin_t bin,
+                                                         double timeBase) {
 
   return get_integer_timestep(bin) * timeBase;
 }
@@ -86,8 +89,8 @@ static INLINE double get_timestep(timebin_t bin, double timeBase) {
  * @param ti_current The current time on the integer time line.
  * @param bin The time bin of interest.
  */
-static INLINE integertime_t get_integer_time_begin(integertime_t ti_current,
-                                                   timebin_t bin) {
+__attribute__((const)) static INLINE integertime_t
+get_integer_time_begin(integertime_t ti_current, timebin_t bin) {
 
   const integertime_t dti = get_integer_timestep(bin);
   if (dti == 0)
@@ -103,8 +106,8 @@ static INLINE integertime_t get_integer_time_begin(integertime_t ti_current,
  * @param ti_current The current time on the integer time line.
  * @param bin The time bin of interest.
  */
-static INLINE integertime_t get_integer_time_end(integertime_t ti_current,
-                                                 timebin_t bin) {
+__attribute__((const)) static INLINE integertime_t
+get_integer_time_end(integertime_t ti_current, timebin_t bin) {
 
   const integertime_t dti = get_integer_timestep(bin);
   if (dti == 0)
@@ -118,7 +121,8 @@ static INLINE integertime_t get_integer_time_end(integertime_t ti_current,
  *
  * @param time The current point on the time line.
  */
-static INLINE timebin_t get_max_active_bin(integertime_t time) {
+__attribute__((const)) static INLINE timebin_t
+get_max_active_bin(integertime_t time) {
 
   if (time == 0) return num_time_bins;
 
@@ -134,8 +138,8 @@ static INLINE timebin_t get_max_active_bin(integertime_t time) {
  * @param ti_current The current point on the time line.
  * @param ti_old The last synchronisation point on the time line.
  */
-static INLINE timebin_t get_min_active_bin(integertime_t ti_current,
-                                           integertime_t ti_old) {
+__attribute__((const)) static INLINE timebin_t
+get_min_active_bin(integertime_t ti_current, integertime_t ti_old) {
 
   const timebin_t min_bin = get_max_active_bin(ti_current - ti_old);
   return (ti_old > 0) ? min_bin : (min_bin - 1);