diff --git a/Tests/NumLib/FeTestData/TestFeHEX20.h b/Tests/NumLib/FeTestData/TestFeHEX20.h
index 2827347b9de8db643068c489fd3398bcdfd08f96..f2d274b3995ec1e4d674a9e16d5abecfedba97b7 100644
--- a/Tests/NumLib/FeTestData/TestFeHEX20.h
+++ b/Tests/NumLib/FeTestData/TestFeHEX20.h
@@ -34,28 +34,31 @@ public:
     /// create a mesh element
     MeshElementType* createMeshElement()
     {
-        // cubic
         auto** nodes = new MeshLib::Node*[e_nnodes];
+        // bottom
         nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
         nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
         nodes[2] = new MeshLib::Node(1.0, 1.0, 0.0);
         nodes[3] = new MeshLib::Node(0.0, 1.0, 0.0);
-
+        // top
         nodes[4] = new MeshLib::Node(0.0, 0.0, 1.0);
         nodes[5] = new MeshLib::Node(1.0, 0.0, 1.0);
         nodes[6] = new MeshLib::Node(1.0, 1.0, 1.0);
         nodes[7] = new MeshLib::Node(0.0, 1.0, 1.0);
 
+        // mid-edge nodes on bottom
         nodes[8] = new MeshLib::Node(0.5, 0.0, 0.0);
         nodes[9] = new MeshLib::Node(1.0, 0.5, 0.0);
         nodes[10] = new MeshLib::Node(0.5, 1.0, 0.0);
         nodes[11] = new MeshLib::Node(0.0, 0.5, 0.0);
 
+        // mid-edge nodes on top
         nodes[12] = new MeshLib::Node(0.5, 0.0, 1.0);
         nodes[13] = new MeshLib::Node(1.0, 0.5, 1.0);
         nodes[14] = new MeshLib::Node(0.5, 1.0, 1.0);
         nodes[15] = new MeshLib::Node(0.0, 0.5, 1.0);
 
+        // mid-edge nodes between the bottom and the top
         nodes[16] = new MeshLib::Node(0.0, 0.0, 0.5);
         nodes[17] = new MeshLib::Node(1.0, 0.0, 0.5);
         nodes[18] = new MeshLib::Node(1.0, 1.0, 0.5);
diff --git a/Tests/NumLib/FeTestData/TestFePRISM15.h b/Tests/NumLib/FeTestData/TestFePRISM15.h
index a1d25ee2030b51fcb362d3bf636a6a93c3328cfd..3e0e67cc7072f487da9e20bc2c932eb42ca40d64 100644
--- a/Tests/NumLib/FeTestData/TestFePRISM15.h
+++ b/Tests/NumLib/FeTestData/TestFePRISM15.h
@@ -34,24 +34,28 @@ public:
     /// create a mesh element
     MeshElementType* createMeshElement()
     {
-        // cubic
         auto** nodes = new MeshLib::Node*[e_nnodes];
+        // vertex nodes on bottom
         nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
         nodes[1] = new MeshLib::Node(0.5, 0.0, 0.0);
         nodes[2] = new MeshLib::Node(0.5, 0.5, 0.0);
 
+        // vertex nodes on top
         nodes[3] = new MeshLib::Node(0.0, 0.0, 0.5);
         nodes[4] = new MeshLib::Node(0.5, 0.0, 0.5);
         nodes[5] = new MeshLib::Node(0.5, 0.5, 0.5);
 
+        // mid-edge nodes on bottom
         nodes[6] = new MeshLib::Node(0.25, 0.0, 0.0);
         nodes[7] = new MeshLib::Node(0.5, 0.25, 0.0);
         nodes[8] = new MeshLib::Node(0.25, 0.25, 0.0);
 
+        // mid-edge nodes on top
         nodes[9] = new MeshLib::Node(0.25, 0.0, 0.5);
         nodes[10] = new MeshLib::Node(0.5, 0.25, 0.5);
         nodes[11] = new MeshLib::Node(0.25, 0.25, 0.5);
 
+        // mid-edge nodes between the bottom and the top
         nodes[12] = new MeshLib::Node(0.0, 0.0, 0.25);
         nodes[13] = new MeshLib::Node(0.5, 0.0, 0.25);
         nodes[14] = new MeshLib::Node(0.5, 0.5, 0.25);
diff --git a/Tests/NumLib/FeTestData/TestFePYRA13.h b/Tests/NumLib/FeTestData/TestFePYRA13.h
index 8f82b3a65a88eec2b5be5e60e73bdbed5ebb9b38..c45772228709029f0ee73eee43ed5d7f58e70e28 100644
--- a/Tests/NumLib/FeTestData/TestFePYRA13.h
+++ b/Tests/NumLib/FeTestData/TestFePYRA13.h
@@ -35,23 +35,26 @@ public:
     /// create a mesh element
     MeshElementType* createMeshElement()
     {
-        // cubic
         auto** nodes = new MeshLib::Node*[e_nnodes];
+        // nodes on base
         nodes[0] = new MeshLib::Node(0.0, 0.0, 0.5);
         nodes[1] = new MeshLib::Node(0.0, 0.0, 0.0);
         nodes[2] = new MeshLib::Node(0.0, 0.5, 0.0);
         nodes[3] = new MeshLib::Node(0.0, 0.5, 0.5);
+        // node at the top
         nodes[4] = new MeshLib::Node(0.25, 0.25, 0.25);
 
+        // mid-edge nodes on base
         nodes[5] = new MeshLib::Node(0.0, 0.0, 0.25);
         nodes[6] = new MeshLib::Node(0.0, 0.25, 0.0);
         nodes[7] = new MeshLib::Node(0.0, 0.5, 0.25);
         nodes[8] = new MeshLib::Node(0.0, 0.25, 0.5);
 
-        nodes[9] = new MeshLib::Node(0.125, 0.125, 0.3625);
+        // mid-edge nodes between the base and the top
+        nodes[9] = new MeshLib::Node(0.125, 0.125, 0.375);
         nodes[10] = new MeshLib::Node(0.125, 0.125, 0.125);
-        nodes[11] = new MeshLib::Node(0.125, 0.3625, 0.125);
-        nodes[12] = new MeshLib::Node(0.125, 0.3625, 0.3625);
+        nodes[11] = new MeshLib::Node(0.125, 0.375, 0.125);
+        nodes[12] = new MeshLib::Node(0.125, 0.375, 0.375);
 
         return new MeshElementType(nodes);
     }
diff --git a/Tests/NumLib/FeTestData/TestFeQUAD8.h b/Tests/NumLib/FeTestData/TestFeQUAD8.h
index ff7a037b20cfc3794a6553573b4b88f89830d7de..9362e6cb67d39936d83ff45b1a3b5f74c23697f1 100644
--- a/Tests/NumLib/FeTestData/TestFeQUAD8.h
+++ b/Tests/NumLib/FeTestData/TestFeQUAD8.h
@@ -15,6 +15,15 @@
 
 namespace FeTestData
 {
+/** Test the eight node iso parameter quadrilateral element
+ * Assuming that the edge points are translated outside of the square shape
+ * element with a distance of b (perturbation) to their owner edges, we can
+ * obtain the the area of the changed element as
+ * \f[
+ *   A = a^2 + \frac{8}{3}a b
+ * \f]
+ * where \f$a\f$ is the length of the edge of the original element.
+ */
 class TestFeQUAD8
 {
 public:
@@ -36,7 +45,6 @@ public:
     /// create a mesh element
     MeshElementType* createMeshElement()
     {
-        // square
         auto** nodes = new MeshLib::Node*[e_nnodes];
         nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
         nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
diff --git a/Tests/NumLib/FeTestData/TestFeQUAD9.h b/Tests/NumLib/FeTestData/TestFeQUAD9.h
index 79a89652835f8ae4b297d1b2a10199a2826a01f8..4e293c6ed8f04af79c3b577d206abc3db0b64e3a 100644
--- a/Tests/NumLib/FeTestData/TestFeQUAD9.h
+++ b/Tests/NumLib/FeTestData/TestFeQUAD9.h
@@ -17,6 +17,15 @@ namespace FeTestData
 {
 class TestFeQUAD9
 {
+    /** Test the nine node iso parameter quadrilateral element
+     * Assuming that the edge points are translated outside of the square shape
+     * element with a distance of b (perturbation) to their owner edges, we can
+     * obtain the the area of the changed element as
+     * \f[
+     *   A = a^2 + \frac{8}{3}a b
+     * \f]
+     * where \f$a\f$ is the length of the edge of the original element.
+     */
 public:
     using ShapeFunction = NumLib::ShapeQuad9;
 
@@ -36,7 +45,6 @@ public:
     /// create a mesh element
     MeshElementType* createMeshElement()
     {
-        // square
         auto** nodes = new MeshLib::Node*[e_nnodes];
         nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
         nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
@@ -48,7 +56,7 @@ public:
         nodes[6] = new MeshLib::Node(0.5, 1.0 - perturbation, 0.0);
         nodes[7] = new MeshLib::Node(perturbation, 0.5, 0.0);
 
-        nodes[8] = new MeshLib::Node(0.5, 0.5, 0.0);
+        nodes[8] = new MeshLib::Node(0.45, 0.25, 0.0);
         return new MeshElementType(nodes);
     }
 
diff --git a/Tests/NumLib/FeTestData/TestFeTET10.h b/Tests/NumLib/FeTestData/TestFeTET10.h
index 474f755710266be735cd1ce210dd3562a9dcec16..820f32febd56d4aea6adce3b5cd48cd66964d73e 100644
--- a/Tests/NumLib/FeTestData/TestFeTET10.h
+++ b/Tests/NumLib/FeTestData/TestFeTET10.h
@@ -36,15 +36,19 @@ public:
     MeshElementType* createMeshElement()
     {
         auto** nodes = new MeshLib::Node*[e_nnodes];
+        // nodes on the base
         nodes[0] = new MeshLib::Node(0.0, 0.0, 0.0);
         nodes[1] = new MeshLib::Node(0.5, 0.5, 0.5);
         nodes[2] = new MeshLib::Node(0.5, 0.0, 0.5);
+        // top node
         nodes[3] = new MeshLib::Node(1.0, 0.5, 0.5);
 
+        // mid-edge nodes on the base
         nodes[4] = new MeshLib::Node(0.25, 0.25, 0.25);
         nodes[5] = new MeshLib::Node(0.5, 0.25, 0.5);
         nodes[6] = new MeshLib::Node(0.25, 0.0, 0.25);
 
+        // mid-edge nodes between the base and the top
         nodes[7] = new MeshLib::Node(0.5, 0.25, 0.25);
         nodes[8] = new MeshLib::Node(0.75, 0.5, 0.5);
         nodes[9] = new MeshLib::Node(0.75, 0.25, 0.5);
diff --git a/Tests/NumLib/FeTestData/TestFeTRI6.h b/Tests/NumLib/FeTestData/TestFeTRI6.h
index 85740b04f130c75fac1ac3fd04196879358193ac..40f072bb48d41315068bf6b0fe1d9bbdd29b18ad 100644
--- a/Tests/NumLib/FeTestData/TestFeTRI6.h
+++ b/Tests/NumLib/FeTestData/TestFeTRI6.h
@@ -16,6 +16,16 @@
 
 namespace FeTestData
 {
+/** Test the six node iso parameter triangle element
+ * Assuming that the edge points are translated outside of the square shape
+ * element with a distance of b (perturbation) to their owner edges, we can
+ * obtain the the area of the changed element as
+ * \f[
+ *   A = A_0 + \sum_{i=1}^3\frac{2}{3}a_i b
+ * \f]
+ * where\f$A_0\f$ is the area of the original triangle, \f$a_i\f$ is the length
+ * of edge \f$i\f$ of the original element.
+ */
 class TestFeTRI6
 {
 public:
@@ -31,6 +41,7 @@ public:
     static const unsigned n_sample_pt_order2 = 3;
     static const unsigned n_sample_pt_order3 = 4;
     static const unsigned global_dim = MeshElementType::dimension;
+    const double perturbation = 0.1;
 
     /// create a mesh element
     MeshElementType* createMeshElement()
@@ -40,14 +51,23 @@ public:
         nodes[1] = new MeshLib::Node(1.0, 0.0, 0.0);
         nodes[2] = new MeshLib::Node(0.0, 1.0, 0.0);
 
-        nodes[3] = new MeshLib::Node(0.5, 0.0, 0.0);
-        nodes[4] = new MeshLib::Node(0.5, 0.5, 0.0);
-        nodes[5] = new MeshLib::Node(0.0, 0.5, 0.0);
+        nodes[3] = new MeshLib::Node(0.5, perturbation, 0.0);
+        const double perturbation_hypotenuse =
+            0.5 * perturbation * std::sqrt(2.0);
+        nodes[4] = new MeshLib::Node(0.5 - perturbation_hypotenuse,
+                                     0.5 - perturbation_hypotenuse, 0.0);
+        nodes[5] = new MeshLib::Node(perturbation, 0.5, 0.0);
 
         return new MeshElementType(nodes);
     }
 
-    double getVolume() const { return 0.5; }
+    double getVolume() const
+    {
+        // The length of hypotenuse is a_h=sqrt(2)
+        // Area = 0.5 - 2 * 2 * a * b /3 - 2 * a_h * b /3
+        // where a=1.0, b=perturbation
+        return 0.27238576250846036;
+    }
 };
 
 }  // namespace