diff --git a/examples/SubgridTests/StellarEvolution/plot_box_evolution.py b/examples/SubgridTests/StellarEvolution/plot_box_evolution.py
index 15fa363ca8e606f1d842b5703f62955b141df987..9636c822f5469ddee8f8f1a45752072b23472f12 100644
--- a/examples/SubgridTests/StellarEvolution/plot_box_evolution.py
+++ b/examples/SubgridTests/StellarEvolution/plot_box_evolution.py
@@ -42,11 +42,11 @@ params = {'axes.labelsize': 10,
 'ytick.labelsize': 10,
 'text.usetex': True,
  'figure.figsize' : (9.90,6.45),
-'figure.subplot.left'    : 0.1,
+'figure.subplot.left'    : 0.05,
 'figure.subplot.right'   : 0.99,
 'figure.subplot.bottom'  : 0.1,
 'figure.subplot.top'     : 0.95,
-'figure.subplot.wspace'  : 0.2,
+'figure.subplot.wspace'  : 0.25,
 'figure.subplot.hspace'  : 0.2,
 'lines.markersize' : 6,
 'lines.linewidth' : 3.,
@@ -183,7 +183,7 @@ xlabel("${\\rm{Time}} (Gyr)$", labelpad=0)
 ylabel("Change in element mass of gas particles (Msun)", labelpad=2)
 xscale("log")
 yscale("log")
-legend(bbox_to_anchor=(2.1, 0.5), ncol=3)
+legend(bbox_to_anchor=(1.005, 1.), ncol=1, fontsize=8, handlelength=1)
 
 # Box gas metal mass --------------------------------
 subplot(224)
diff --git a/examples/SubgridTests/StellarEvolution/stellar_evolution.yml b/examples/SubgridTests/StellarEvolution/stellar_evolution.yml
index 36279011f842863940d3bb5e6a6b8dc5c07f0b72..431dce995ef8908bf171c03268a9ad06f6c9c3a1 100644
--- a/examples/SubgridTests/StellarEvolution/stellar_evolution.yml
+++ b/examples/SubgridTests/StellarEvolution/stellar_evolution.yml
@@ -17,17 +17,16 @@ Cosmology:
 
 # Parameters governing the time integration
 TimeIntegration:
-  time_begin: 0     # The starting time of the simulation (in internal units).
-  time_end:   5e-3  # The end time of the simulation (in internal units).
-  dt_min:     1e-12 # The minimal time-step size of the simulation (in internal units).
-  dt_max:     8e-6  # The maximal time-step size of the simulation (in internal units).
+  time_begin: 0      # The starting time of the simulation (in internal units).
+  time_end:   1.3e-2 # The end time of the simulation (in internal units).
+  dt_min:     1e-10  # The minimal time-step size of the simulation (in internal units).
+  dt_max:     8e-6   # The maximal time-step size of the simulation (in internal units).
 
 # Parameters governing the snapshots
 Snapshots:
   basename:            stellar_evolution # Common part of the name of output files
   time_first:          0.       # Time of the first output (in internal units)
-  delta_time:          1.e-5     # Time difference between consecutive outputs without cosmology (internal units)
-  compression:         4
+  delta_time:          3.e-5    # Time difference between consecutive outputs without cosmology (internal units)
 
 # Parameters governing the conserved quantities statistics
 Statistics:
diff --git a/src/feedback/EAGLE/feedback.c b/src/feedback/EAGLE/feedback.c
index 8b97bb510a540df9093bcf7a3863c71e4f9300fe..206af3dd120faa30775762c8e837a7eeff348d51 100644
--- a/src/feedback/EAGLE/feedback.c
+++ b/src/feedback/EAGLE/feedback.c
@@ -350,6 +350,8 @@ inline static void evolve_SNIa(const float log10_min_mass,
  * IMF weighted by the yields read from tables for each of the quantities of
  * interest.
  *
+ * Note for Matthieu: This function is poorly written and needs improving.
+ *
  * @param log10_min_mass log10 mass at the end of step
  * @param log10_max_mass log10 mass at the beginning of step
  * @param stellar_yields array to store calculated yields for passing to
@@ -504,6 +506,8 @@ inline static void evolve_SNII(float log10_min_mass, float log10_max_mass,
  * IMF weighted by the yields read from tables for each of the quantities of
  * interest.
  *
+ * Note for Matthieu: This function is poorly written and needs improving.
+ *
  * @param log10_min_mass log10 mass at the end of step
  * @param log10_max_mass log10 mass at the beginning of step
  * @param stellar_yields array to store calculated yields for passing to
@@ -646,6 +650,10 @@ inline static void evolve_AGB(const float log10_min_mass, float log10_max_mass,
   }
 }
 
+/**
+ * @brief Prepares a star's feedback field before computing what
+ * needs to be distributed.
+ */
 inline static void feedback_init_to_distribute(struct spart* sp) {
 
   /* Zero the amount of mass that is distributed */
diff --git a/src/feedback/EAGLE/feedback_iact.h b/src/feedback/EAGLE/feedback_iact.h
index 63cecf017126484535572cac5192729af12b349e..25c21ef29ae2f92a3c64afad0c2a95abf590dc45 100644
--- a/src/feedback/EAGLE/feedback_iact.h
+++ b/src/feedback/EAGLE/feedback_iact.h
@@ -140,77 +140,82 @@ runner_iact_nonsym_feedback_apply(const float r2, const float *dx,
     pj->chemistry_data.metal_mass_fraction[elem] = new_metal_mass / new_mass;
   }
 
-  /* /\* Update iron mass fraction from SNIa  *\/ */
-  /* const float current_iron_from_SNIa_mass = */
-  /*     pj->chemistry_data.iron_mass_fraction_from_SNIa * current_mass; */
-  /* const float new_iron_from_SNIa_mass = */
-  /*     current_iron_from_SNIa_mass + */
-  /*     si->feedback_data.to_distribute.Fe_mass_from_SNIa *
-   * Omega_frac; */
-  /* pj->chemistry_data.iron_mass_fraction_from_SNIa = */
-  /*     new_iron_from_SNIa_mass / new_mass; */
-
-  /* /\* Update mass fraction from SNIa  *\/ */
-  /* const float current_mass_from_SNIa = */
-  /*     pj->chemistry_data.mass_from_SNIa * current_mass; */
-  /* const float new_mass_from_SNIa = */
-  /*     current_mass_from_SNIa + */
-  /*     si->feedback_data.to_distribute.mass_from_SNIa * Omega_frac;
-   */
-  /* pj->chemistry_data.mass_from_SNIa = new_mass_from_SNIa / new_mass; */
-
-  /* /\* Update metal mass fraction from SNIa *\/ */
-  /* const float current_metal_mass_fraction_from_SNIa = */
-  /*     pj->chemistry_data.metal_mass_fraction_from_SNIa * current_mass; */
-  /* const float new_metal_mass_fraction_from_SNIa = */
-  /*     current_metal_mass_fraction_from_SNIa + */
-  /*     si->feedback_data.to_distribute.metal_mass_from_SNIa * */
-  /*         Omega_frac; */
-  /* pj->chemistry_data.metal_mass_fraction_from_SNIa = */
-  /*     new_metal_mass_fraction_from_SNIa / new_mass; */
-
-  /* /\* Update mass fraction from SNII  *\/ */
-  /* const float current_mass_from_SNII = */
-  /*     pj->chemistry_data.mass_from_SNII * current_mass; */
-  /* const float new_mass_from_SNII = */
-  /*     current_mass_from_SNII + */
-  /*     si->feedback_data.to_distribute.mass_from_SNII * Omega_frac;
-   */
-  /* pj->chemistry_data.mass_from_SNII = new_mass_from_SNII / new_mass; */
-
-  /* /\* Update metal mass fraction from SNII *\/ */
-  /* const float current_metal_mass_fraction_from_SNII = */
-  /*     pj->chemistry_data.metal_mass_fraction_from_SNII * current_mass; */
-  /* const float new_metal_mass_fraction_from_SNII = */
-  /*     current_metal_mass_fraction_from_SNII + */
-  /*     si->feedback_data.to_distribute.metal_mass_from_SNII * */
-  /*         Omega_frac; */
-  /* pj->chemistry_data.metal_mass_fraction_from_SNII = */
-  /*     new_metal_mass_fraction_from_SNII / new_mass; */
-
-  /* /\* Update mass fraction from AGB  *\/ */
-  /* const float current_mass_from_AGB = */
-  /*     pj->chemistry_data.mass_from_AGB * current_mass; */
-  /* const float new_mass_from_AGB = */
-  /*     current_mass_from_AGB + */
-  /*     si->feedback_data.to_distribute.mass_from_AGB * Omega_frac;
-   */
-  /* pj->chemistry_data.mass_from_AGB = new_mass_from_AGB / new_mass; */
-
-  /* /\* Update metal mass fraction from AGB *\/ */
-  /* const float current_metal_mass_fraction_from_AGB = */
-  /*     pj->chemistry_data.metal_mass_fraction_from_AGB * current_mass; */
-  /* const float new_metal_mass_fraction_from_AGB = */
-  /*     current_metal_mass_fraction_from_AGB + */
-  /*     si->feedback_data.to_distribute.metal_mass_from_AGB * */
-  /*         Omega_frac; */
-  /* pj->chemistry_data.metal_mass_fraction_from_AGB = */
-  /*     new_metal_mass_fraction_from_AGB / new_mass; */
+  /* Update iron mass fraction from SNIa  */
+  const double current_iron_from_SNIa_mass =
+      pj->chemistry_data.iron_mass_fraction_from_SNIa * current_mass;
+  const double delta_iron_from_SNIa_mass =
+      si->feedback_data.to_distribute.Fe_mass_from_SNIa * Omega_frac;
+  const double new_iron_from_SNIa_mass =
+      current_iron_from_SNIa_mass + delta_iron_from_SNIa_mass;
+
+  pj->chemistry_data.iron_mass_fraction_from_SNIa =
+      new_iron_from_SNIa_mass / new_mass;
+
+  /* Update mass fraction from SNIa  */
+  const double current_mass_from_SNIa =
+      pj->chemistry_data.mass_from_SNIa * current_mass;
+  const double delta_mass_from_SNIa =
+      si->feedback_data.to_distribute.mass_from_SNIa * Omega_frac;
+  const double new_mass_from_SNIa =
+      current_mass_from_SNIa + delta_mass_from_SNIa;
+
+  pj->chemistry_data.mass_from_SNIa = new_mass_from_SNIa / new_mass;
+
+  /* Update metal mass fraction from SNIa */
+  const double current_metal_mass_from_SNIa =
+      pj->chemistry_data.metal_mass_fraction_from_SNIa * current_mass;
+  const double delta_metal_mass_from_SNIa =
+      si->feedback_data.to_distribute.metal_mass_from_SNIa * Omega_frac;
+  const double new_metal_mass_from_SNIa =
+      current_metal_mass_from_SNIa + delta_metal_mass_from_SNIa;
+
+  pj->chemistry_data.metal_mass_fraction_from_SNIa =
+      new_metal_mass_from_SNIa / new_mass;
+
+  /* Update mass fraction from SNII  */
+  const double current_mass_from_SNII =
+      pj->chemistry_data.mass_from_SNII * current_mass;
+  const double delta_mass_from_SNII =
+      si->feedback_data.to_distribute.mass_from_SNII * Omega_frac;
+  const double new_mass_from_SNII =
+      current_mass_from_SNII + delta_mass_from_SNII;
+
+  pj->chemistry_data.mass_from_SNII = new_mass_from_SNII / new_mass;
+
+  /* Update metal mass fraction from SNII */
+  const double current_metal_mass_from_SNII =
+      pj->chemistry_data.metal_mass_fraction_from_SNII * current_mass;
+  const double delta_metal_mass_from_SNII =
+      si->feedback_data.to_distribute.metal_mass_from_SNII * Omega_frac;
+  const double new_metal_mass_from_SNII =
+      current_metal_mass_from_SNII + delta_metal_mass_from_SNII;
+
+  pj->chemistry_data.metal_mass_fraction_from_SNII =
+      new_metal_mass_from_SNII / new_mass;
+
+  /* Update mass fraction from AGB  */
+  const double current_mass_from_AGB =
+      pj->chemistry_data.mass_from_AGB * current_mass;
+  const double delta_mass_from_AGB =
+      si->feedback_data.to_distribute.mass_from_AGB * Omega_frac;
+  const double new_mass_from_AGB = current_mass_from_AGB + delta_mass_from_AGB;
+
+  pj->chemistry_data.mass_from_AGB = new_mass_from_AGB / new_mass;
+
+  /* Update metal mass fraction from AGB */
+  const double current_metal_mass_from_AGB =
+      pj->chemistry_data.metal_mass_fraction_from_AGB * current_mass;
+  const double delta_metal_mass_from_AGB =
+      si->feedback_data.to_distribute.metal_mass_from_AGB * Omega_frac;
+  const double new_metal_mass_from_AGB =
+      current_metal_mass_from_AGB + delta_metal_mass_from_AGB;
+
+  pj->chemistry_data.metal_mass_fraction_from_AGB =
+      new_metal_mass_from_AGB / new_mass;
 
   /* /\* Update momentum *\/ */
   /* for (int i = 0; i < 3; i++) { */
-  /*   pj->v[i] += si->feedback_data.to_distribute.mass * Omega_frac
-   * * */
+  /*   pj->v[i] += si->feedback_data.to_distribute.mass * Omega_frac * */
   /*               (si->v[i] - pj->v[i]); */
   /* } */
 
diff --git a/src/feedback/EAGLE/feedback_struct.h b/src/feedback/EAGLE/feedback_struct.h
index f5950772a6f61671fef205693edb43852c820741..42ee43f3cdecfd6865c28cd7edf683f0c81c9504 100644
--- a/src/feedback/EAGLE/feedback_struct.h
+++ b/src/feedback/EAGLE/feedback_struct.h
@@ -19,6 +19,8 @@
 #ifndef SWIFT_FEEDBACK_STRUCT_EAGLE_H
 #define SWIFT_FEEDBACK_STRUCT_EAGLE_H
 
+#include "chemistry_struct.h"
+
 /**
  * @brief Feedback fields carried by each star particles
  */
@@ -26,6 +28,9 @@ struct feedback_spart_data {
 
   union {
 
+    /**
+     * @brief Values collected from the gas neighbours.
+     */
     struct {
 
       /*! Inverse of normalisation factor used for the enrichment */
@@ -36,6 +41,9 @@ struct feedback_spart_data {
 
     } to_collect;
 
+    /**
+     * @brief Values to be distributed to the gas neighbours.
+     */
     struct {
 
       /*! Normalisation factor used for the enrichment */