diff --git a/.gitignore b/.gitignore
index 28a830818af36faad3f4278c6adcba5562b59ee7..a8778a3b82b2cf89f2169f2580d16b0975bb41ba 100644
--- a/.gitignore
+++ b/.gitignore
@@ -35,11 +35,16 @@ examples/*/*/used_parameters.yml
 examples/*/gravity_checks_*.dat
 
 tests/testPair
+tests/brute_force_periodic_BC_standard.dat
+tests/swift_periodic_BC_standard.dat
+tests/brute_force_periodic_BC_pertrubed.dat
+tests/swift_periodic_BC_perturbed.dat
 tests/brute_force_standard.dat
 tests/swift_dopair_standard.dat
 tests/brute_force_perturbed.dat
 tests/swift_dopair_perturbed.dat
 tests/test27cells
+tests/testPeriodicBC
 tests/test125cells
 tests/brute_force_27_standard.dat
 tests/swift_dopair_27_standard.dat
@@ -64,6 +69,8 @@ tests/testMaths
 tests/testThreadpool
 tests/testParser
 tests/parser_output.yml
+tests/testPeriodicBC.sh
+tests/testPeriodicBCPerturbed.sh
 tests/test27cells.sh
 tests/test27cellsPerturbed.sh
 tests/test125cells.sh
diff --git a/configure.ac b/configure.ac
index 788bb57eed801c1a1dff2204b57b34c4fadf3b58..75a333a3d8f6554752dc70978c3ca45aecb09dc9 100644
--- a/configure.ac
+++ b/configure.ac
@@ -854,13 +854,23 @@ AC_CONFIG_FILES([tests/test27cells.sh], [chmod +x tests/test27cells.sh])
 AC_CONFIG_FILES([tests/test27cellsPerturbed.sh], [chmod +x tests/test27cellsPerturbed.sh])
 AC_CONFIG_FILES([tests/test125cells.sh], [chmod +x tests/test125cells.sh])
 AC_CONFIG_FILES([tests/test125cellsPerturbed.sh], [chmod +x tests/test125cellsPerturbed.sh])
+AC_CONFIG_FILES([tests/testPeriodicBC.sh], [chmod +x tests/testPeriodicBC.sh])
+AC_CONFIG_FILES([tests/testPeriodicBCPerturbed.sh], [chmod +x tests/testPeriodicBCPerturbed.sh])
 AC_CONFIG_FILES([tests/testParser.sh], [chmod +x tests/testParser.sh])
 
 # Save the compilation options
 AC_DEFINE_UNQUOTED([SWIFT_CONFIG_FLAGS],["$swift_config_flags"],[Flags passed to configure])
 
+# Make sure the latest git revision string gets included
+touch src/version.c
+
+# Generate output.
+AC_OUTPUT
+
 # Report general configuration.
-AC_MSG_RESULT([
+AC_MSG_RESULT([ 
+ ------- Summary --------
+   
    Compiler        : $CC
     - vendor       : $ax_cv_c_compiler_vendor
     - version      : $ax_cv_c_compiler_version
@@ -891,10 +901,5 @@ AC_MSG_RESULT([
    Task debugging    : $enable_task_debugging
    Debugging checks  : $enable_debugging_checks
    Gravity checks    : $gravity_force_checks
-])
-
-# Make sure the latest git revision string gets included
-touch src/version.c
 
-# Generate output.
-AC_OUTPUT
+ ------------------------])
diff --git a/src/runner.c b/src/runner.c
index 36074144e14ec31741ee5964d8e2b98c9673192b..c4922e3e389a5c5b3a8d93d904d706938bbe4524 100644
--- a/src/runner.c
+++ b/src/runner.c
@@ -1822,7 +1822,7 @@ void *runner_main(void *data) {
             runner_dopair1_branch_density(r, ci, cj);
 #ifdef EXTRA_HYDRO_LOOP
           else if (t->subtype == task_subtype_gradient)
-            runner_dopair1_gradient(r, ci, cj);
+            runner_dopair1_branch_gradient(r, ci, cj);
 #endif
           else if (t->subtype == task_subtype_force)
             runner_dopair2_force(r, ci, cj);
diff --git a/src/tools.c b/src/tools.c
index 73684c82662870d368f7dd360c84635654f06434..b2932f3eb09bfe2acaf435a18e25b0421984542a 100644
--- a/src/tools.c
+++ b/src/tools.c
@@ -37,6 +37,7 @@
 #include "gravity.h"
 #include "hydro.h"
 #include "part.h"
+#include "periodic.h"
 #include "runner.h"
 
 /**
@@ -181,6 +182,7 @@ void pairs_all_density(struct runner *r, struct cell *ci, struct cell *cj) {
 
   float r2, hi, hj, hig2, hjg2, dx[3];
   struct part *pi, *pj;
+  const double dim[3] = {r->e->s->dim[0], r->e->s->dim[1], r->e->s->dim[2]};
 
   /* Implements a double-for loop and checks every interaction */
   for (int i = 0; i < ci->count; ++i) {
@@ -197,6 +199,7 @@ void pairs_all_density(struct runner *r, struct cell *ci, struct cell *cj) {
       r2 = 0.0f;
       for (int k = 0; k < 3; k++) {
         dx[k] = ci->parts[i].x[k] - cj->parts[j].x[k];
+        dx[k] = nearest(dx[k], dim[k]);
         r2 += dx[k] * dx[k];
       }
 
@@ -224,6 +227,7 @@ void pairs_all_density(struct runner *r, struct cell *ci, struct cell *cj) {
       r2 = 0.0f;
       for (int k = 0; k < 3; k++) {
         dx[k] = cj->parts[j].x[k] - ci->parts[i].x[k];
+        dx[k] = nearest(dx[k], dim[k]);
         r2 += dx[k] * dx[k];
       }
 
@@ -241,6 +245,7 @@ void pairs_all_force(struct runner *r, struct cell *ci, struct cell *cj) {
 
   float r2, hi, hj, hig2, hjg2, dx[3];
   struct part *pi, *pj;
+  const double dim[3] = {r->e->s->dim[0], r->e->s->dim[1], r->e->s->dim[2]};
 
   /* Implements a double-for loop and checks every interaction */
   for (int i = 0; i < ci->count; ++i) {
@@ -259,6 +264,7 @@ void pairs_all_force(struct runner *r, struct cell *ci, struct cell *cj) {
       r2 = 0.0f;
       for (int k = 0; k < 3; k++) {
         dx[k] = ci->parts[i].x[k] - cj->parts[j].x[k];
+        dx[k] = nearest(dx[k], dim[k]);
         r2 += dx[k] * dx[k];
       }
 
@@ -288,6 +294,7 @@ void pairs_all_force(struct runner *r, struct cell *ci, struct cell *cj) {
       r2 = 0.0f;
       for (int k = 0; k < 3; k++) {
         dx[k] = cj->parts[j].x[k] - ci->parts[i].x[k];
+        dx[k] = nearest(dx[k], dim[k]);
         r2 += dx[k] * dx[k];
       }
 
diff --git a/tests/Makefile.am b/tests/Makefile.am
index 92f60384eb7539d3a1cd357105e58df419cc739a..b356086b54b6f4f07c8bedafe1fccb543121871e 100644
--- a/tests/Makefile.am
+++ b/tests/Makefile.am
@@ -25,7 +25,8 @@ TESTS = testGreetings testMaths testReading.sh testSingle testKernel testSymmetr
         testParser.sh testSPHStep test125cells.sh test125cellsPerturbed.sh testFFT \
         testAdiabaticIndex testRiemannExact testRiemannTRRS testRiemannHLLC \
         testMatrixInversion testThreadpool testDump testLogger \
-        testVoronoi1D testVoronoi2D testVoronoi3D
+        testVoronoi1D testVoronoi2D testVoronoi3D \
+	testPeriodicBC.sh testPeriodicBCPerturbed.sh
 
 # List of test programs to compile
 check_PROGRAMS = testGreetings testReading testSingle testTimeIntegration \
@@ -34,7 +35,7 @@ check_PROGRAMS = testGreetings testReading testSingle testTimeIntegration \
                  testSymmetry testThreadpool benchmarkInteractions \
                  testAdiabaticIndex testRiemannExact testRiemannTRRS \
                  testRiemannHLLC testMatrixInversion testDump testLogger \
-		 testVoronoi1D testVoronoi2D testVoronoi3D
+		 testVoronoi1D testVoronoi2D testVoronoi3D testPeriodicBC
 
 # Sources for the individual programs
 testGreetings_SOURCES = testGreetings.c
@@ -55,6 +56,8 @@ testPair_SOURCES = testPair.c
 
 test27cells_SOURCES = test27cells.c
 
+testPeriodicBC_SOURCES = testPeriodicBC.c
+
 test125cells_SOURCES = test125cells.c
 
 testParser_SOURCES = testParser.c
@@ -91,9 +94,9 @@ testLogger_SOURCES = testLogger.c
 
 # Files necessary for distribution
 EXTRA_DIST = testReading.sh makeInput.py testPair.sh testPairPerturbed.sh \
-	     test27cells.sh test27cellsPerturbed.sh testParser.sh \
-	     test125cells.sh test125cellsPerturbed.sh testParserInput.yaml difffloat.py \
-	     tolerance_125_normal.dat tolerance_125_perturbed.dat \
+	     test27cells.sh test27cellsPerturbed.sh testParser.sh testPeriodicBC.sh \
+	     testPeriodicBCPerturbed.sh test125cells.sh test125cellsPerturbed.sh testParserInput.yaml \
+	     difffloat.py tolerance_125_normal.dat tolerance_125_perturbed.dat \
              tolerance_27_normal.dat tolerance_27_perturbed.dat tolerance_27_perturbed_h.dat \
 	     tolerance_pair_normal.dat tolerance_pair_perturbed.dat \
-	     fft_params.yml
+	     fft_params.yml tolerance_periodic_BC_normal.dat tolerance_periodic_BC_perturbed.dat
diff --git a/tests/difffloat.py b/tests/difffloat.py
index 8c8713a8d2968005a86a9c1e7fdce34bf3e624f7..ddcf7bcb29758afa3429dea8bcf50e1c5c0477dc 100644
--- a/tests/difffloat.py
+++ b/tests/difffloat.py
@@ -42,11 +42,6 @@ if len(sys.argv) >= 4:
 if len(sys.argv) >= 5:
     number_to_check = int(sys.argv[4])
 
-if len(sys.argv) == 6:
-    ignoreSmallRhoDh = int(sys.argv[5])
-else:
-    ignoreSmallRhoDh = 0
-
 # Get the particle properties being compared from the header.
 with open(file1, 'r') as f:
   line = f.readline()
@@ -69,7 +64,7 @@ n_lines = shape(data1)[0]
 n_columns = shape(data1)[1]
 
 if fileTol != "":
-    if n_linesTol != 2:
+    if n_linesTol != 3:
         print "Incorrect number of lines in tolerance file '%s'."%fileTol
     if n_columnsTol != n_columns:
         print "Incorrect number of columns in tolerance file '%s'."%fileTol
@@ -79,10 +74,12 @@ if fileTol == "":
     print "Relative difference tolerance:", rel_tol
     absTol = ones(n_columns) * abs_tol
     relTol = ones(n_columns) * rel_tol
+    limTol = zeros(n_columns)
 else:
     print "Tolerances read from file"
     absTol = dataTol[0,:]
     relTol = dataTol[1,:]
+    limTol = dataTol[2,:]
 
 n_lines_to_check = 0
 if number_to_check > 0:
@@ -113,11 +110,8 @@ for i in range(n_lines_to_check):
             print ""
             error = True
 
-        if abs(data1[i,j]) < 4e-6 and abs(data2[i,j]) < 4e-6 : continue
+        if abs(data1[i,j]) + abs(data2[i,j]) < limTol[j] : continue
 
-        # Ignore pathological cases with rho_dh
-        if ignoreSmallRhoDh and j == 8 and abs(data1[i,j]) < 2e-4: continue
-        
         if( rel_diff > 1.1*relTol[j]):
             print "Relative difference larger than tolerance (%e) for particle %d, column %s:"%(relTol[j], data1[i,0], part_props[j])
             print "%10s:           a = %e"%("File 1", data1[i,j])
diff --git a/tests/test125cells.c b/tests/test125cells.c
index 9c113f98c89e6213d8c01fcd4515da076787e93e..0f3939c834cd3c552cea7882211321b62cc7d966 100644
--- a/tests/test125cells.c
+++ b/tests/test125cells.c
@@ -623,7 +623,7 @@ int main(int argc, char *argv[]) {
 /* Do the density calculation */
 #if !(defined(MINIMAL_SPH) && defined(WITH_VECTORIZATION))
 
-    /* Initialise the particle cache. */
+/* Initialise the particle cache. */
 #ifdef WITH_VECTORIZATION
     runner.ci_cache.count = 0;
     cache_init(&runner.ci_cache, 512);
diff --git a/tests/testPair.c b/tests/testPair.c
index 2ca2092cf1bb41c4910ee2639c22677e2fb1dbd2..92068e286fce1570010cf83862a6cdc9a5ad4ffc 100644
--- a/tests/testPair.c
+++ b/tests/testPair.c
@@ -45,7 +45,6 @@
 #define DOPAIR1_NAME "runner_dopair1_density"
 #endif
 
-
 /* n is both particles per axis and box size:
  * particles are generated on a mesh with unit spacing
  */
diff --git a/tests/testPeriodicBC.c b/tests/testPeriodicBC.c
new file mode 100644
index 0000000000000000000000000000000000000000..19faf252745e18624f9ed5d2a3cbfc42098ef49e
--- /dev/null
+++ b/tests/testPeriodicBC.c
@@ -0,0 +1,587 @@
+/*******************************************************************************
+ * This file is part of SWIFT.
+ * Copyright (C) 2015 Matthieu Schaller (matthieu.schaller@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"
+
+/* Some standard headers. */
+#include <fenv.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <unistd.h>
+
+/* Local headers. */
+#include "swift.h"
+
+#define ACC_THRESHOLD 1e-5
+
+#if defined(WITH_VECTORIZATION)
+#define DOSELF1 runner_doself1_density_vec
+#define DOPAIR1 runner_dopair1_branch_density
+#define DOSELF1_NAME "runner_doself1_density_vec"
+#define DOPAIR1_NAME "runner_dopair1_density_vec"
+#endif
+
+#ifndef DOSELF1
+#define DOSELF1 runner_doself1_density
+#define DOSELF1_NAME "runner_doself1_density"
+#endif
+
+#ifndef DOPAIR1
+#define DOPAIR1 runner_dopair1_branch_density
+#define DOPAIR1_NAME "runner_dopair1_density"
+#endif
+
+enum velocity_types {
+  velocity_zero,
+  velocity_random,
+  velocity_divergent,
+  velocity_rotating
+};
+
+/**
+ * @brief Constructs a cell and all of its particle in a valid state prior to
+ * a DOPAIR or DOSELF calcuation.
+ *
+ * @param n The cube root of the number of particles.
+ * @param offset The position of the cell offset from (0,0,0).
+ * @param size The cell size.
+ * @param h The smoothing length of the particles in units of the inter-particle
+ *separation.
+ * @param density The density of the fluid.
+ * @param partId The running counter of IDs.
+ * @param pert The perturbation to apply to the particles in the cell in units
+ *of the inter-particle separation.
+ * @param vel The type of velocity field (0, random, divergent, rotating)
+ */
+struct cell *make_cell(size_t n, double *offset, double size, double h,
+                       double density, long long *partId, double pert,
+                       enum velocity_types vel) {
+  const size_t count = n * n * n;
+  const double volume = size * size * size;
+  struct cell *cell = malloc(sizeof(struct cell));
+  bzero(cell, sizeof(struct cell));
+
+  if (posix_memalign((void **)&cell->parts, part_align,
+                     count * sizeof(struct part)) != 0) {
+    error("couldn't allocate particles, no. of particles: %d", (int)count);
+  }
+  bzero(cell->parts, count * sizeof(struct part));
+
+  float h_max = 0.f;
+
+  /* Construct the parts */
+  struct part *part = cell->parts;
+  for (size_t x = 0; x < n; ++x) {
+    for (size_t y = 0; y < n; ++y) {
+      for (size_t z = 0; z < n; ++z) {
+        part->x[0] =
+            offset[0] +
+            size * (x + 0.5 + random_uniform(-0.5, 0.5) * pert) / (float)n;
+        part->x[1] =
+            offset[1] +
+            size * (y + 0.5 + random_uniform(-0.5, 0.5) * pert) / (float)n;
+        part->x[2] =
+            offset[2] +
+            size * (z + 0.5 + random_uniform(-0.5, 0.5) * pert) / (float)n;
+        switch (vel) {
+          case velocity_zero:
+            part->v[0] = 0.f;
+            part->v[1] = 0.f;
+            part->v[2] = 0.f;
+            break;
+          case velocity_random:
+            part->v[0] = random_uniform(-0.05, 0.05);
+            part->v[1] = random_uniform(-0.05, 0.05);
+            part->v[2] = random_uniform(-0.05, 0.05);
+            break;
+          case velocity_divergent:
+            part->v[0] = part->x[0] - 1.5 * size;
+            part->v[1] = part->x[1] - 1.5 * size;
+            part->v[2] = part->x[2] - 1.5 * size;
+            break;
+          case velocity_rotating:
+            part->v[0] = part->x[1];
+            part->v[1] = -part->x[0];
+            part->v[2] = 0.f;
+            break;
+        }
+        part->h = size * h / (float)n;
+        h_max = fmax(h_max, part->h);
+        part->id = ++(*partId);
+
+#if defined(GIZMO_SPH) || defined(SHADOWFAX_SPH)
+        part->conserved.mass = density * volume / count;
+
+#ifdef SHADOWFAX_SPH
+        double anchor[3] = {0., 0., 0.};
+        double side[3] = {1., 1., 1.};
+        voronoi_cell_init(&part->cell, part->x, anchor, side);
+#endif
+
+#else
+        part->mass = density * volume / count;
+#endif
+
+#if defined(HOPKINS_PE_SPH)
+        part->entropy = 1.f;
+        part->entropy_one_over_gamma = 1.f;
+#endif
+
+        part->time_bin = 1;
+
+#ifdef SWIFT_DEBUG_CHECKS
+        part->ti_drift = 8;
+        part->ti_kick = 8;
+#endif
+
+        ++part;
+      }
+    }
+  }
+
+  /* Cell properties */
+  cell->split = 0;
+  cell->h_max = h_max;
+  cell->count = count;
+  cell->dx_max_part = 0.;
+  cell->dx_max_sort = 0.;
+  cell->width[0] = size;
+  cell->width[1] = size;
+  cell->width[2] = size;
+  cell->loc[0] = offset[0];
+  cell->loc[1] = offset[1];
+  cell->loc[2] = offset[2];
+
+  cell->ti_old_part = 8;
+  cell->ti_end_min = 8;
+  cell->ti_end_max = 8;
+  cell->ti_sort = 8;
+
+  shuffle_particles(cell->parts, cell->count);
+
+  cell->sorted = 0;
+  cell->sort = NULL;
+  cell->sortsize = 0;
+
+  return cell;
+}
+
+void clean_up(struct cell *ci) {
+  free(ci->parts);
+  free(ci->sort);
+  free(ci);
+}
+
+/**
+ * @brief Initializes all particles field to be ready for a density calculation
+ */
+void zero_particle_fields(struct cell *c) {
+  for (int pid = 0; pid < c->count; pid++) {
+    hydro_init_part(&c->parts[pid], NULL);
+  }
+}
+
+/**
+ * @brief Ends the loop by adding the appropriate coefficients
+ */
+void end_calculation(struct cell *c) {
+  for (int pid = 0; pid < c->count; pid++) {
+    hydro_end_density(&c->parts[pid]);
+  }
+}
+
+/**
+ * @brief Dump all the particles to a file
+ */
+void dump_particle_fields(char *fileName, struct cell *main_cell, int i, int j,
+                          int k) {
+  FILE *file = fopen(fileName, "a");
+
+  /* Write header */
+  fprintf(file,
+          "# %4s %10s %10s %10s %10s %10s %10s %13s %13s %13s %13s %13s "
+          "%13s %13s %13s\n",
+          "ID", "pos_x", "pos_y", "pos_z", "v_x", "v_y", "v_z", "rho", "rho_dh",
+          "wcount", "wcount_dh", "div_v", "curl_vx", "curl_vy", "curl_vz");
+
+  fprintf(file, "# Centre cell at (i,j,k)=(%d, %d, %d) ---------------------\n",
+          i, j, k);
+
+  /* Write main cell */
+  for (int pid = 0; pid < main_cell->count; pid++) {
+    fprintf(file,
+            "%6llu %10f %10f %10f %10f %10f %10f %13e %13e %13e %13e %13e "
+            "%13e %13e %13e\n",
+            main_cell->parts[pid].id, main_cell->parts[pid].x[0],
+            main_cell->parts[pid].x[1], main_cell->parts[pid].x[2],
+            main_cell->parts[pid].v[0], main_cell->parts[pid].v[1],
+            main_cell->parts[pid].v[2],
+            hydro_get_density(&main_cell->parts[pid]),
+#if defined(GIZMO_SPH) || defined(SHADOWFAX_SPH)
+            0.f,
+#else
+            main_cell->parts[pid].density.rho_dh,
+#endif
+            main_cell->parts[pid].density.wcount,
+            main_cell->parts[pid].density.wcount_dh,
+#if defined(GADGET2_SPH) || defined(DEFAULT_SPH) || defined(HOPKINS_PE_SPH)
+            main_cell->parts[pid].density.div_v,
+            main_cell->parts[pid].density.rot_v[0],
+            main_cell->parts[pid].density.rot_v[1],
+            main_cell->parts[pid].density.rot_v[2]
+#else
+            0., 0., 0., 0.
+#endif
+            );
+  }
+  fclose(file);
+}
+
+/**
+ * @brief Compares the vectorised result against
+ * the serial result of the interaction.
+ *
+ * @param serial_parts Particle array that has been interacted serially
+ * @param vec_parts Particle array to be interacted using vectors
+ * @param count No. of particles that have been interacted
+ * @param threshold Level of accuracy needed
+ *
+ * @return Non-zero value if difference found, 0 otherwise
+ */
+int check_results(struct part *serial_parts, struct part *vec_parts, int count,
+                  double threshold) {
+  int result = 0;
+
+  for (int i = 0; i < count; i++)
+    result += compare_particles(serial_parts[i], vec_parts[i], threshold);
+
+  return result;
+}
+
+/* Just a forward declaration... */
+void runner_doself1_density(struct runner *r, struct cell *ci);
+void runner_doself1_density_vec(struct runner *r, struct cell *ci);
+void runner_dopair1_branch_density(struct runner *r, struct cell *ci,
+                                   struct cell *cj);
+
+void test_boundary_conditions(struct cell **cells, struct runner runner,
+                              const int loc_i, const int loc_j, const int loc_k,
+                              const int dim, char *swiftOutputFileName,
+                              char *bruteForceOutputFileName) {
+
+  /* Store the main cell for future use */
+  struct cell *main_cell = cells[loc_i * (dim * dim) + loc_j * dim + loc_k];
+
+  /* Zero the fields */
+  for (int j = 0; j < 512; ++j) zero_particle_fields(cells[j]);
+
+/* Run all the pairs */
+#if !(defined(MINIMAL_SPH) && defined(WITH_VECTORIZATION))
+
+#ifdef WITH_VECTORIZATION
+  runner.ci_cache.count = 0;
+  cache_init(&runner.ci_cache, 512);
+  runner.cj_cache.count = 0;
+  cache_init(&runner.cj_cache, 512);
+#endif
+
+  /* Now loop over all the neighbours of this cell
+   * and perform the pair interactions. */
+  for (int ii = -1; ii < 2; ii++) {
+    int iii = loc_i + ii;
+    iii = (iii + dim) % dim;
+    for (int jj = -1; jj < 2; jj++) {
+      int jjj = loc_j + jj;
+      jjj = (jjj + dim) % dim;
+      for (int kk = -1; kk < 2; kk++) {
+        int kkk = loc_k + kk;
+        kkk = (kkk + dim) % dim;
+
+        /* Get the neighbouring cell */
+        struct cell *cj = cells[iii * (dim * dim) + jjj * dim + kkk];
+
+        if (cj != main_cell) DOPAIR1(&runner, main_cell, cj);
+      }
+    }
+  }
+
+  /* And now the self-interaction */
+
+  DOSELF1(&runner, main_cell);
+
+#endif
+
+  /* Let's get physical ! */
+  end_calculation(main_cell);
+
+  /* Dump particles from the main cell. */
+  dump_particle_fields(swiftOutputFileName, main_cell, loc_i, loc_j, loc_k);
+
+  /* Now perform a brute-force version for accuracy tests */
+
+  /* Zero the fields */
+  for (int i = 0; i < 512; ++i) zero_particle_fields(cells[i]);
+
+#if !(defined(MINIMAL_SPH) && defined(WITH_VECTORIZATION))
+
+  /* Now loop over all the neighbours of this cell
+   * and perform the pair interactions. */
+  for (int ii = -1; ii < 2; ii++) {
+    int iii = loc_i + ii;
+    iii = (iii + dim) % dim;
+    for (int jj = -1; jj < 2; jj++) {
+      int jjj = loc_j + jj;
+      jjj = (jjj + dim) % dim;
+      for (int kk = -1; kk < 2; kk++) {
+        int kkk = loc_k + kk;
+        kkk = (kkk + dim) % dim;
+
+        /* Get the neighbouring cell */
+        struct cell *cj = cells[iii * (dim * dim) + jjj * dim + kkk];
+
+        if (cj != main_cell) pairs_all_density(&runner, main_cell, cj);
+      }
+    }
+  }
+
+  /* And now the self-interaction */
+  self_all_density(&runner, main_cell);
+
+#endif
+
+  /* Let's get physical ! */
+  end_calculation(main_cell);
+
+  /* Dump */
+  dump_particle_fields(bruteForceOutputFileName, main_cell, loc_i, loc_j,
+                       loc_k);
+}
+
+/* And go... */
+int main(int argc, char *argv[]) {
+
+  engine_pin();
+  size_t runs = 0, particles = 0;
+  double h = 1.23485, size = 1., rho = 1.;
+  double perturbation = 0.;
+  double threshold = ACC_THRESHOLD;
+  char outputFileNameExtension[200] = "";
+  char swiftOutputFileName[200] = "";
+  char bruteForceOutputFileName[200] = "";
+  enum velocity_types vel = velocity_zero;
+
+  /* Initialize CPU frequency, this also starts time. */
+  unsigned long long cpufreq = 0;
+  clocks_set_cpufreq(cpufreq);
+
+  /* Choke on FP-exceptions */
+  feenableexcept(FE_DIVBYZERO | FE_INVALID | FE_OVERFLOW);
+
+  /* Get some randomness going */
+  srand(0);
+
+  char c;
+  while ((c = getopt(argc, argv, "m:s:h:n:r:t:d:f:v:a:")) != -1) {
+    switch (c) {
+      case 'h':
+        sscanf(optarg, "%lf", &h);
+        break;
+      case 's':
+        sscanf(optarg, "%lf", &size);
+        break;
+      case 'n':
+        sscanf(optarg, "%zu", &particles);
+        break;
+      case 'r':
+        sscanf(optarg, "%zu", &runs);
+        break;
+      case 'd':
+        sscanf(optarg, "%lf", &perturbation);
+        break;
+      case 'm':
+        sscanf(optarg, "%lf", &rho);
+        break;
+      case 'f':
+        strcpy(outputFileNameExtension, optarg);
+        break;
+      case 'v':
+        sscanf(optarg, "%d", (int *)&vel);
+        break;
+      case 'a':
+        sscanf(optarg, "%lf", &threshold);
+        break;
+      case '?':
+        error("Unknown option.");
+        break;
+    }
+  }
+
+  if (h < 0 || particles == 0 || runs == 0) {
+    printf(
+        "\nUsage: %s -n PARTICLES_PER_AXIS -r NUMBER_OF_RUNS [OPTIONS...]\n"
+        "\nGenerates 27 cells, filled with particles on a Cartesian grid."
+        "\nThese are then interacted using runner_dopair1_density() and "
+        "runner_doself1_density()."
+        "\n\nOptions:"
+        "\n-h DISTANCE=1.2348 - Smoothing length in units of <x>"
+        "\n-m rho             - Physical density in the cell"
+        "\n-s size            - Physical size of the cell"
+        "\n-d pert            - Perturbation to apply to the particles [0,1["
+        "\n-v type (0,1,2,3)  - Velocity field: (zero, random, divergent, "
+        "rotating)"
+        "\n-f fileName        - Part of the file name used to save the dumps\n",
+        argv[0]);
+    exit(1);
+  }
+
+  /* Help users... */
+  message("DOSELF1 function called: %s", DOSELF1_NAME);
+  message("DOPAIR1 function called: %s", DOPAIR1_NAME);
+  message("Vector size: %d", VEC_SIZE);
+  message("Adiabatic index: ga = %f", hydro_gamma);
+  message("Hydro implementation: %s", SPH_IMPLEMENTATION);
+  message("Smoothing length: h = %f", h * size);
+  message("Kernel:               %s", kernel_name);
+  message("Neighbour target: N = %f", pow_dimension(h) * kernel_norm);
+  message("Density target: rho = %f", rho);
+  message("div_v target:   div = %f", vel == 2 ? 3.f : 0.f);
+  message("curl_v target: curl = [0., 0., %f]", vel == 3 ? -2.f : 0.f);
+
+  printf("\n");
+
+  /* Build the infrastructure */
+  struct space space;
+  space.periodic = 1;
+  space.dim[0] = 8.;
+  space.dim[1] = 8.;
+  space.dim[2] = 8.;
+
+  struct hydro_props hp;
+  hp.h_max = FLT_MAX;
+
+  struct engine engine;
+  engine.s = &space;
+  engine.time = 0.1f;
+  engine.ti_current = 8;
+  engine.max_active_bin = num_time_bins;
+  engine.hydro_properties = &hp;
+
+  struct runner runner;
+  runner.e = &engine;
+
+  /* Construct some cells */
+  struct cell *cells[512];
+  const int dim = 8;
+  static long long partId = 0;
+  for (int i = 0; i < dim; ++i) {
+    for (int j = 0; j < dim; ++j) {
+      for (int k = 0; k < dim; ++k) {
+        double offset[3] = {i * size, j * size, k * size};
+        cells[i * (dim * dim) + j * dim + k] = make_cell(
+            particles, offset, size, h, rho, &partId, perturbation, vel);
+
+        runner_do_drift_part(&runner, cells[i * (dim * dim) + j * dim + k], 0);
+
+        runner_do_sort(&runner, cells[i * (dim * dim) + j * dim + k], 0x1FFF,
+                       0);
+      }
+    }
+  }
+
+  /* Create output file names. */
+  sprintf(swiftOutputFileName, "swift_periodic_BC_%s.dat",
+          outputFileNameExtension);
+  sprintf(bruteForceOutputFileName, "brute_force_periodic_BC_%s.dat",
+          outputFileNameExtension);
+
+  /* Delete files if they already exist. */
+  remove(swiftOutputFileName);
+  remove(bruteForceOutputFileName);
+
+  const int half_dim = (dim - 1) / 2;
+
+  /* Test the periodic boundary conditions for each of the 8 corners. */
+  test_boundary_conditions(cells, runner, 0, 0, 0, dim, swiftOutputFileName,
+                           bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, 0, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, 0, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, 0, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, dim - 1, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, dim - 1, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, dim - 1, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, dim - 1, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+
+  /* Test the boundary conditions for cells at the centre of each face of the
+   * box. */
+  test_boundary_conditions(cells, runner, half_dim, half_dim, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, half_dim, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, half_dim, half_dim, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, half_dim, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, half_dim, 0, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, half_dim, dim - 1, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+
+  /* Test the boundary conditions for cells at the centre of each edge of the
+   * box. */
+  test_boundary_conditions(cells, runner, half_dim, dim - 1, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, dim - 1, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, half_dim, dim - 1, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, dim - 1, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+
+  test_boundary_conditions(cells, runner, 0, half_dim, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, half_dim, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, half_dim, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, half_dim, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+
+  test_boundary_conditions(cells, runner, half_dim, 0, 0, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, dim - 1, 0, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, half_dim, 0, dim - 1, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+  test_boundary_conditions(cells, runner, 0, 0, half_dim, dim,
+                           swiftOutputFileName, bruteForceOutputFileName);
+
+  /* Clean things to make the sanitizer happy ... */
+  for (int i = 0; i < 512; ++i) clean_up(cells[i]);
+
+  return 0;
+}
diff --git a/tests/testPeriodicBC.sh.in b/tests/testPeriodicBC.sh.in
new file mode 100755
index 0000000000000000000000000000000000000000..075acc0b68686bd2f418cf457140b3d6b93093d5
--- /dev/null
+++ b/tests/testPeriodicBC.sh.in
@@ -0,0 +1,30 @@
+#!/bin/bash
+
+for v in {0..3}
+do
+    echo ""
+	
+    rm -f brute_force_periodic_BC_standard.dat swift_periodic_BC_standard.dat
+
+    echo "Running ./testPeriodicBC -n 6 -r 1 -d 0 -f standard -v $v"
+    ./testPeriodicBC -n 6 -r 1 -d 0 -f standard -v $v
+
+    if [ -e brute_force_periodic_BC_standard.dat ]
+    then
+	if python @srcdir@/difffloat.py brute_force_periodic_BC_standard.dat swift_periodic_BC_standard.dat @srcdir@/tolerance_periodic_BC_normal.dat
+  then
+	    echo "Accuracy test passed"
+	else
+	    echo "Accuracy test failed"
+	    exit 1
+	fi
+    else
+	echo "Error Missing test output file"
+	exit 1
+    fi
+
+    echo "------------"
+    
+done
+
+exit $?
diff --git a/tests/testPeriodicBCPerturbed.sh.in b/tests/testPeriodicBCPerturbed.sh.in
new file mode 100755
index 0000000000000000000000000000000000000000..ac190d5a80654154dcd329e69c1c9cc9fe45833a
--- /dev/null
+++ b/tests/testPeriodicBCPerturbed.sh.in
@@ -0,0 +1,30 @@
+#!/bin/bash
+
+for v in {0..3}
+do
+    echo ""
+	
+    rm -f brute_force_periodic_BC_perturbed.dat swift_periodic_BC_perturbed.dat
+
+    echo "Running ./testPeriodicBC -n 6 -r 1 -d 0.1 -f perturbed -v $v"
+    ./testPeriodicBC -n 6 -r 1 -d 0.1 -f perturbed -v $v
+
+    if [ -e brute_force_periodic_BC_perturbed.dat ]
+    then
+	if python @srcdir@/difffloat.py brute_force_periodic_BC_perturbed.dat swift_periodic_BC_perturbed.dat @srcdir@/tolerance_periodic_BC_perturbed.dat
+	then
+	    echo "Accuracy test passed"
+	else
+	    echo "Accuracy test failed"
+	    exit 1
+	fi
+    else
+	echo "Error Missing test output file"
+	exit 1
+    fi
+
+    echo "------------"
+    
+done
+
+exit $?
diff --git a/tests/tolerance_125_normal.dat b/tests/tolerance_125_normal.dat
index c9ad23d4472c46e64e8418e46c5fe71f813b23b5..0f11d03507b23c76b5703e118eede1359fe2afba 100644
--- a/tests/tolerance_125_normal.dat
+++ b/tests/tolerance_125_normal.dat
@@ -1,3 +1,4 @@
 #   ID    pos_x    pos_y    pos_z      v_x      v_y      v_z        h      rho    div_v        S        u        P        c      a_x      a_y      a_z     h_dt    v_sig    dS/dt    du/dt
     0	  1e-4	   1e-4	    1e-4       1e-4	1e-4	 1e-4	    1e-4   1e-4	  1e-4	       1e-4	1e-4	 1e-4	  1e-4	 1e-4	  1e-4	   1e-4	   1e-4	   1e-4	    1e-4     1e-4
     0	  1e-4	   1e-4	    1e-4       1e-4	1e-4	 1e-4	    1e-4   1e-4	  1e-4	       1e-4	1e-4	 1e-4	  1e-4	 1e-4	  1e-4	   1e-4	   1e-4	   1e-4	    1e-4     1e-4
+    0	  1e-6	   1e-6	    1e-6       1e-6	1e-6	 1e-6	    1e-6   1e-6	  1e-6	       1e-6	1e-6	 1e-6	  1e-6	 1e-5	  1e-5	   1e-5	   1e-5	   1e-5	    1e-5     1e-5
diff --git a/tests/tolerance_125_perturbed.dat b/tests/tolerance_125_perturbed.dat
index a45a419c1b484be226591ddf5ce16c8ef328e9a5..71922dca49da8fb38df78cacab4bbe1e61e912bd 100644
--- a/tests/tolerance_125_perturbed.dat
+++ b/tests/tolerance_125_perturbed.dat
@@ -1,3 +1,4 @@
 #   ID    pos_x    pos_y    pos_z      v_x      v_y      v_z        h      rho    div_v        S        u        P        c      a_x      a_y      a_z     h_dt    v_sig    dS/dt    du/dt
     0	  1e-4	   1e-4	    1e-4       1e-4	1e-4	 1e-4	    1e-4   1e-4	  1e-4	       1e-4	1e-4	 1e-4	  1e-4	 1e-4	  1e-4	   1e-4	   1e-4	   1e-4	    1e-4     1e-4
-    0	  1e-4	   1e-4	    1e-4       1e-4	1e-4	 1e-4	    1e-4   1e-4	  1e-4	       1e-4	1e-4	 1e-4	  1e-4	 1e-2	  1e-2	   1e-2	   1e-4	   1e-4	    1e-4     1e-4
+    0	  1e-4	   1e-4	    1e-4       1e-4	1e-4	 1e-4	    1e-4   1e-4	  1e-4	       1e-4	1e-4	 1e-4	  1e-4	 1e-3	  1e-3	   1e-3	   1e-4	   1e-4	    1e-4     1e-4
+    0	  1e-6	   1e-6	    1e-6       1e-6	1e-6	 1e-6	    1e-6   1e-6	  1e-6	       1e-6	1e-6	 1e-6	  1e-6	 2e-4	  2e-4	   2e-4	   1e-6	   1e-6	    1e-6     1e-6
diff --git a/tests/tolerance_27_normal.dat b/tests/tolerance_27_normal.dat
index 31ee002bb9c73ff8d74cce545aff715476b33507..3ec58173c9c5a0465382de9094b5d15b81b88a2a 100644
--- a/tests/tolerance_27_normal.dat
+++ b/tests/tolerance_27_normal.dat
@@ -1,3 +1,4 @@
 #   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   2e-6	      4e-5	    2e-4       2e-3		 1e-5	     6e-6	   6e-6		 6e-6
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1.2e-4	    1e-4       1e-4		 2e-4	     1e-4	   1e-4	 	 1e-4
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1e-6	    1e-6       1e-6		 1e-6	     1e-6	   1e-6		 1e-6
diff --git a/tests/tolerance_27_perturbed.dat b/tests/tolerance_27_perturbed.dat
index cd31b4c5aa25059e5abc012ba3d84a1c88a2c746..ed7df2ca9f55dd6211ed6321ffc40faa37d98092 100644
--- a/tests/tolerance_27_perturbed.dat
+++ b/tests/tolerance_27_perturbed.dat
@@ -1,3 +1,4 @@
 #   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   2e-6       1e-4	    2e-4       3e-3		 1e-5	     3e-6	   3e-6		 7e-6
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      2e-3	    1e-5       1e-4		 6e-5	     2e-3	   2e-3	 	 2e-3
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      4e-4	    1e-6       1e-6		 1e-6	     2e-6	   2e-6		 2e-6
diff --git a/tests/tolerance_27_perturbed_h.dat b/tests/tolerance_27_perturbed_h.dat
index 43521a94c4339adbcf43f2e46ae578a101c23bc6..5e70ef2a173a0347714a057b5ae7110e511337c2 100644
--- a/tests/tolerance_27_perturbed_h.dat
+++ b/tests/tolerance_27_perturbed_h.dat
@@ -1,3 +1,4 @@
 #   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   2e-6       1e-4	    2e-4       4e-3		 1e-5	     3e-6	   3e-6		 7e-6
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1.4e-2	    1e-5       1e-4		 2.5e-4	     3e-3	   3e-3	 	 3e-3
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1e-6	    1e-6       1e-6		 1e-6	     3e-6	   3e-6		 3e-6
diff --git a/tests/tolerance_pair_normal.dat b/tests/tolerance_pair_normal.dat
index f5031c5f47dfa203300ebcc9a47fbac42f854d26..3fca591dce8156a08c970e4a6579b56b2de12553 100644
--- a/tests/tolerance_pair_normal.dat
+++ b/tests/tolerance_pair_normal.dat
@@ -1,3 +1,4 @@
 #   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-5	      1e-5	    2e-5       3e-2		 1e-5	     1e-5	   1e-5		 1e-5
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-5	      1.2e-5	    1e-5       1e-2		 1e-4	     1e-4	   1e-4		 1e-4
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1e-6	    1e-6       1e-6		 1e-6	     1e-6	   1e-6		 1e-6
diff --git a/tests/tolerance_pair_perturbed.dat b/tests/tolerance_pair_perturbed.dat
index ca58ff45995158e031eca6b60eec498aa6c627ef..59a80a7b1fa0db373b7cb43b9793146330b17d53 100644
--- a/tests/tolerance_pair_perturbed.dat
+++ b/tests/tolerance_pair_perturbed.dat
@@ -1,3 +1,4 @@
 #   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-5	      1e-5	    2e-5       3e-2		 1e-5	     1e-5	   1e-5		 1e-5
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-3	      4e-4	    8e-3       2e-2		 1e-4	     1.6e-4	   1.6e-4	 1.6e-4
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1e-6	    1e-6       1e-6		 1e-6	     1e-6	   1e-6		 1e-6
diff --git a/tests/tolerance_periodic_BC_normal.dat b/tests/tolerance_periodic_BC_normal.dat
new file mode 100644
index 0000000000000000000000000000000000000000..a7e8bfdaef6d3b4e49708e624042ef975b9e8e42
--- /dev/null
+++ b/tests/tolerance_periodic_BC_normal.dat
@@ -0,0 +1,4 @@
+#   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   4e-6	      4e-5	    3e-4       3e-3		 2e-4	     1e-4	   1e-4		 1e-4
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   2e-6	      2e-4	    1e-4       1e-4		 6e-4	     2e-3	   2e-3	 	 2e-3
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1e-4	    1e-6       1e-4		 5e-4	     2e-4	   2e-4	 	 2e-4
diff --git a/tests/tolerance_periodic_BC_perturbed.dat b/tests/tolerance_periodic_BC_perturbed.dat
new file mode 100644
index 0000000000000000000000000000000000000000..78e573ee01811bfdc5c72f59b20f17fe00a54348
--- /dev/null
+++ b/tests/tolerance_periodic_BC_perturbed.dat
@@ -0,0 +1,4 @@
+#   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   3e-6	      4e-5	    2e-4       3e-3		 2e-4	     1e-4	   1e-4		 1e-4
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   2e-6	      6e-3	    1e-4       1e-4		 1e-2	     6e-3	   6e-3	 	 6e-3
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      2e-3	    1e-6       1e-4		 5e-4	     3e-3	   3e-3 	 3e-3