diff --git a/src/cosmology.c b/src/cosmology.c
index 8b02f0cab6761694ff58188ad8d530ddbeaf7183..d304354f3e5068137096c971fed693104924cdb3 100644
--- a/src/cosmology.c
+++ b/src/cosmology.c
@@ -374,31 +374,45 @@ void cosmology_init_tables(struct cosmology *c) {
                       GSL_INTEG_GAUSS61, space, &result, &abserr);
   c->universe_age_at_present_day = result;
 
-  /* Inverse t(a) */
-  const double time_init = c->time_interp_table_offset;
+
+  /* Update the times */
+  c->time_begin = cosmology_get_time_since_big_bang(c, c->a_begin);
+  c->time_end = cosmology_get_time_since_big_bang(c, c->a_end);
+
+  /*
+   * Inverse t(a)
+   */
+
   const double delta_t =
-      (c->universe_age_at_present_day - time_init) / cosmology_table_length;
+    (c->time_end - c->time_begin) / cosmology_table_length;
 
-  int i_prev = 0;
-  for (int i = 0; i < cosmology_table_length; i++) {
-    /* Current time */
-    double time_interp = delta_t * i;
+  /* index in the time_interp_table */
+  int i_a = 0;
+
+  for (int i_time = 0; i_time < cosmology_table_length; i_time++) {
+    /* Current time
+     * time_interp_table = \int_a_begin^a => no need of time_begin */
+    double time_interp = delta_t * (i_time + 1);
 
     /* Find next time in time_interp_table */
-    while (i_prev < cosmology_table_length &&
-           c->time_interp_table[i_prev] <= time_interp) {
-      i_prev++;
+    while (i_a < cosmology_table_length &&
+           c->time_interp_table[i_a] <= time_interp) {
+      i_a++;
     }
 
     /* Find linear interpolation scaling */
-    double scale = time_interp - c->time_interp_table[i_prev - 1];
-    scale /= c->time_interp_table[i_prev] - c->time_interp_table[i_prev - 1];
-    scale += i_prev;
+    double scale = 0;
+    if (i_a != cosmology_table_length) {
+      scale = time_interp - c->time_interp_table[i_a - 1];
+      scale /= c->time_interp_table[i_a] - c->time_interp_table[i_a - 1];
+    }
+
+    scale += i_a;
 
     /* Compute interpolated scale factor */
     double log_a = c->log_a_begin + scale * (c->log_a_end - c->log_a_begin) /
-                                        cosmology_table_length;
-    c->scale_factor_interp_table[i] = exp(log_a) - c->a_begin;
+      cosmology_table_length;
+    c->scale_factor_interp_table[i_time] = exp(log_a) - c->a_begin;
   }
 
   /* Free the workspace and temp array */
@@ -672,9 +686,6 @@ double cosmology_get_delta_time(const struct cosmology *c,
 /**
  * @brief Compute scale factor from time since big bang (in internal units).
  *
- * WARNING: This method has a low accuracy at high redshift.
- * The relative error is around 1e-3 (testCosmology.c is measuring it).
- *
  * @param c The current #cosmology.
  * @param t time since the big bang
  * @return The scale factor.
diff --git a/tests/testCosmology.c b/tests/testCosmology.c
index 698351ad952e7d0b5f7d8e354c45a1a2dd53f968..bafad55471453f7308d1498daa15dbae3a76a6bc 100644
--- a/tests/testCosmology.c
+++ b/tests/testCosmology.c
@@ -24,7 +24,7 @@
 #include "swift.h"
 
 #define N_CHECK 20
-#define TOLERANCE 1e-3
+#define TOLERANCE 1e-7
 
 void test_params_init(struct swift_params *params) {
   parser_init("", params);
@@ -72,5 +72,6 @@ int main(int argc, char *argv[]) {
 
   message("Everything seems fine with cosmology.");
 
+  cosmology_clean(&cosmo);
   return 0;
 }