diff --git a/src/tools.c b/src/tools.c
index d25b7401a1e0515c650333b41193d54b5e155d39..8fc393144ffec7071968e20d3d1f8db0d19cc0bb 100644
--- a/src/tools.c
+++ b/src/tools.c
@@ -471,3 +471,33 @@ void engine_single_force(double *dim, long long int pid,
           p.a_hydro[1], p.a_hydro[2]);
   fflush(stdout);
 }
+
+/**
+ * Returns a random number (uniformly distributed) in [a,b[
+ */
+double random_uniform(double a, double b) {
+  return (rand() / (double)RAND_MAX) * (b - a) + a;
+}
+
+/**
+ * @brief Randomly shuffle an array of particles.
+ */
+void shuffle_particles(struct part *parts, const int count) {
+
+  if(count > 1) {
+    
+    for(int i=0; i<count - 1; i++) {
+      
+      int j = i + random_uniform(0.,(double)(count - 1 - i));
+      
+      struct part particle = parts[j];
+
+      parts[j] = parts[i];
+
+      parts[i] = particle;
+    }
+
+  }
+  else error("Array not big enough to shuffle!");
+
+}
diff --git a/src/tools.h b/src/tools.h
index ccffc77ceb8a967fd40c3737651ba75d529eee0f..8e0212652922fb4afd1ac89a64d4a01c71ecd4a9 100644
--- a/src/tools.h
+++ b/src/tools.h
@@ -38,4 +38,7 @@ void self_all_density(struct runner *r, struct cell *ci);
 void pairs_n2(double *dim, struct part *__restrict__ parts, int N,
               int periodic);
 
+double random_uniform(double a, double b);
+void shuffle_particles(struct part *parts, const int count);
+
 #endif /* SWIFT_TOOL_H */