diff --git a/harmonicmaps-skyrmions-eoc-testsuite/compute-disc-error-skyrmions-hexagon.parset b/harmonicmaps-skyrmions-eoc-testsuite/compute-disc-error-skyrmions-hexagon.parset
new file mode 100644
index 0000000000000000000000000000000000000000..0aae9f15f9ed1355e49518a16e900ccebdc5d15d
--- /dev/null
+++ b/harmonicmaps-skyrmions-eoc-testsuite/compute-disc-error-skyrmions-hexagon.parset
@@ -0,0 +1,23 @@
+#############################################
+#  Grid parameters
+#############################################
+
+structuredGrid = false
+
+path = /home/sander/data/liquidcrystals
+gridFile = hexagon.msh
+
+targetDim = 3
+targetSpace = UnitVector
+
+#numLevels = 1
+#numReferenceLevels = 2
+
+#order = 2
+
+#simulationData = cosserat-continuum-result-1.data
+
+#referenceData = cosserat-continuum-result-2.data
+
+discretizationErrorMode = discrete
+
diff --git a/harmonicmaps-skyrmions-eoc-testsuite/harmonicmaps-skyrmions-hexagon.parset b/harmonicmaps-skyrmions-eoc-testsuite/harmonicmaps-skyrmions-hexagon.parset
new file mode 100644
index 0000000000000000000000000000000000000000..a6c84913541c453153dc8ea05f3b9c39a9eca93f
--- /dev/null
+++ b/harmonicmaps-skyrmions-eoc-testsuite/harmonicmaps-skyrmions-hexagon.parset
@@ -0,0 +1,69 @@
+#############################################
+#  Grid parameters
+#############################################
+
+structuredGrid = false
+path = /home/sander/data/liquidcrystals
+gridFile = hexagon.msh
+
+# Number of grid levels
+numLevels = 3
+
+order = 1
+
+#############################################
+#  Solver parameters
+#############################################
+
+# Tolerance of the trust region solver
+tolerance = 1e-6
+
+# Max number of steps of the trust region solver
+maxTrustRegionSteps = 100
+
+# Initial trust-region radius
+initialTrustRegionRadius = 1
+
+# Number of multigrid iterations per trust-region step
+numIt = 200
+
+# Number of presmoothing steps
+nu1 = 3
+
+# Number of postsmoothing steps
+nu2 = 3
+
+# Number of coarse grid corrections
+mu = 1
+
+# Number of base solver iterations
+baseIt = 100
+
+# Tolerance of the multigrid solver
+mgTolerance = 1e-10
+
+# Tolerance of the base grid solver
+baseTolerance = 1e-8
+
+############################
+#   Problem specifications
+############################
+
+# Type of energy we are minimizing
+energy = chiral_skyrmion
+
+[energyParameters]
+
+h = 3
+kappa = 1
+
+[]
+
+# Inverse stereographic projection
+initialIterate = inverse-stereographic-projection
+
+dirichletValues = [0, 0, 1]
+
+startFromFile = yes
+initialIterateGridFilename = hexagon.msh
+initialIterateFilename = initial_skyrmion_hexagon_5_1.vtu
diff --git a/harmonicmaps-skyrmions-eoc-testsuite/testsuite.sh b/harmonicmaps-skyrmions-eoc-testsuite/testsuite.sh
new file mode 100755
index 0000000000000000000000000000000000000000..dbef4aed319527ee224818477734d219451a817c
--- /dev/null
+++ b/harmonicmaps-skyrmions-eoc-testsuite/testsuite.sh
@@ -0,0 +1,76 @@
+#!/bin/bash
+
+set -e
+
+runComputation(){
+
+numLevels=$1
+order=$2
+
+LOGFILE="./harmonicmaps_${order}_${numLevels}.log"
+
+#  run the actual simulation
+./harmonicmaps-${order} harmonicmaps-skyrmions-hexagon.parset -numLevels ${numLevels} | tee ${LOGFILE}
+
+}
+
+
+MAXPROCS=1
+
+order=1
+numReferenceLevels=10
+
+for numLevels in $(seq 5 ${numReferenceLevels}); do
+
+  # Do one simulation run
+  runComputation $numLevels ${order}
+
+  # Measure the discretization errors against the solution on the finest grid
+  LOGFILE="./compute-disc-error_${order}_${numLevels}.log"
+
+  ../build-cmake/src/compute-disc-error compute-disc-error-skyrmions-hexagon.parset \
+                                        -order ${order} \
+                                        -numLevels ${numLevels} \
+                                        -numReferenceLevels ${numReferenceLevels} \
+                                        -simulationData harmonicmaps-result-${order}-${numLevels}.data \
+                                        -referenceData  harmonicmaps-result-${order}-${numReferenceLevels}.data | tee ${LOGFILE}
+done
+
+order=2
+numReferenceLevels=9
+
+for numLevels in $(seq 4 ${numReferenceLevels}); do
+
+  # Do one simulation run
+  runComputation $numLevels $order
+
+  # Measure the discretization errors against the solution on the finest grid
+  LOGFILE="./compute-disc-error_${order}_${numLevels}.log"
+
+  ../build-cmake/src/compute-disc-error compute-disc-error-skyrmions-hexagon.parset \
+                                        -order ${order} \
+                                        -numLevels ${numLevels} \
+                                        -numReferenceLevels ${numReferenceLevels} \
+                                        -simulationData harmonicmaps-result-${order}-${numLevels}.data \
+                                        -referenceData  harmonicmaps-result-${order}-${numReferenceLevels}.data | tee ${LOGFILE}
+done
+
+order=3
+numReferenceLevels=8
+
+for numLevels in $(seq 3 $numReferenceLevels); do
+
+  # Do one simulation run
+  runComputation $numLevels $order
+
+  # Measure the discretization errors against the solution on the finest grid
+  LOGFILE="./compute-disc-error_${order}_${numLevels}.log"
+
+  ../build-cmake/src/compute-disc-error compute-disc-error-skyrmions-hexagon.parset \
+                                        -order ${order} \
+                                        -numLevels ${numLevels} \
+                                        -numReferenceLevels ${numReferenceLevels} \
+                                        -simulationData harmonicmaps-result-${order}-${numLevels}.data \
+                                        -referenceData  harmonicmaps-result-${order}-${numReferenceLevels}.data | tee ${LOGFILE}
+done
+