diff --git a/.gitignore b/.gitignore
index 5dc1aacf..ec633d61 100644
--- a/.gitignore
+++ b/.gitignore
@@ -27,3 +27,4 @@ testMesh.txt
 *.synctex.gz
 compile_commands.json
 build/
+.vscode/*
\ No newline at end of file
diff --git a/CMakeLists.txt b/CMakeLists.txt
index deadc296..b0245b1b 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -23,7 +23,7 @@ find_package (Threads REQUIRED)
 
 find_package(precice REQUIRED CONFIG)
 
-find_package(Boost 1.65.1 REQUIRED COMPONENTS system program_options filesystem)
+find_package(Boost 1.65.1 REQUIRED COMPONENTS system program_options filesystem unit_test_framework)
 
 find_package(VTK REQUIRED)
 if (VTK_FOUND)
@@ -62,6 +62,15 @@ if(METIS_FOUND)
   target_link_libraries(preciceMap metisAPI)
 endif()
 
+add_executable(testing src/tests/testing.cpp src/tests/read_test.cpp src/tests/write_test.cpp src/mesh.cpp)
+target_include_directories(testing PRIVATE src)
+target_link_libraries(testing
+  Boost::filesystem
+  Boost::boost
+  Boost::unit_test_framework
+  ${VTK_LIBRARIES}
+)
+
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/make_mesh.py           DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/visualize_partition.py DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/partition_mesh.py      DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
@@ -70,3 +79,8 @@ file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/eval_mesh.py           DESTINATION ${C
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/vtk_calculator.py      DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/mesh_io.py             DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
 file(COPY ${CMAKE_CURRENT_SOURCE_DIR}/src/mesh.py                DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
+
+enable_testing()
+
+add_test(read_write_test testing)
+
diff --git a/src/mesh.cpp b/src/mesh.cpp
index b6b5075c..e558bed3 100644
--- a/src/mesh.cpp
+++ b/src/mesh.cpp
@@ -9,14 +9,6 @@
 #include <sstream>
 #include <stdexcept>
 
-#include <vtkDoubleArray.h>
-#include <vtkGenericDataObjectReader.h>
-#include <vtkPointData.h>
-#include <vtkPoints.h>
-#include <vtkSmartPointer.h>
-#include <vtkUnstructuredGrid.h>
-#include <vtkUnstructuredGridWriter.h>
-
 #include <vtkCell.h>
 #include <vtkCellArray.h>
 #include <vtkDoubleArray.h>
@@ -41,11 +33,21 @@ std::string MeshName::filename() const
   return _mname + ".vtk";
 }
 
+namespace aste {
+
+}
+namespace {
+Mesh::VID vtkToPos(vtkIdType id)
+{
+  assert(id >= 0);
+  return static_cast<Mesh::VID>(id);
+}
+} // namespace
+
 namespace {
 // Reads the main file containing the vertices and data
 void readMainFile(Mesh &mesh, const std::string &filename, const std::string &dataname, const int &dim)
 {
-
   if (!fs::is_regular_file(filename)) {
     throw std::invalid_argument{"The mesh file does not exist: " + filename};
   }
@@ -116,19 +118,20 @@ void readMainFile(Mesh &mesh, const std::string &filename, const std::string &da
   }
 
   for (int i = 0; i < reader->GetUnstructuredGridOutput()->GetNumberOfCells(); i++) {
-    int cellType = reader->GetUnstructuredGridOutput()->GetCell(1)->GetCellType();
+    int cellType = reader->GetUnstructuredGridOutput()->GetCell(i)->GetCellType();
 
+    //Here we use static cast since VTK library returns a long long unsigned int however preCICE uses int for PointId's
     if (cellType == VTK_TRIANGLE) {
-      vtkCell *             cell = reader->GetUnstructuredGridOutput()->GetCell(i);
-      std::array<size_t, 3> elem{cell->GetPointId(0), cell->GetPointId(1), cell->GetPointId(2)};
+      vtkCell *                cell = reader->GetUnstructuredGridOutput()->GetCell(i);
+      std::array<Mesh::VID, 3> elem{vtkToPos(cell->GetPointId(0)), vtkToPos(cell->GetPointId(1)), vtkToPos(cell->GetPointId(2))};
       mesh.triangles.push_back(elem);
     } else if (cellType == VTK_LINE) {
-      vtkCell *             cell = reader->GetUnstructuredGridOutput()->GetCell(i);
-      std::array<size_t, 2> elem{cell->GetPointId(0), cell->GetPointId(1)};
+      vtkCell *                cell = reader->GetUnstructuredGridOutput()->GetCell(i);
+      std::array<Mesh::VID, 2> elem{vtkToPos(cell->GetPointId(0)), vtkToPos(cell->GetPointId(1))};
       mesh.edges.push_back(elem);
     } else if (cellType == VTK_QUAD) {
-      vtkCell *             cell = reader->GetUnstructuredGridOutput()->GetCell(i);
-      std::array<size_t, 4> elem{cell->GetPointId(0), cell->GetPointId(1), cell->GetPointId(2), cell->GetPointId(3)};
+      vtkCell *                cell = reader->GetUnstructuredGridOutput()->GetCell(i);
+      std::array<Mesh::VID, 4> elem{vtkToPos(cell->GetPointId(0)), vtkToPos(cell->GetPointId(1)), vtkToPos(cell->GetPointId(2)), vtkToPos(cell->GetPointId(3))};
       mesh.quadrilaterals.push_back(elem);
     } else {
       throw std::runtime_error{
@@ -155,7 +158,7 @@ void MeshName::createDirectories() const
 
 void MeshName::save(const Mesh &mesh, const std::string &dataname) const
 {
-  const int                            numComp          = mesh.positions.size() / mesh.data.size();
+  const int                            numComp          = mesh.data.size() / mesh.positions.size();
   vtkSmartPointer<vtkUnstructuredGrid> unstructuredGrid = vtkSmartPointer<vtkUnstructuredGrid>::New();
   vtkSmartPointer<vtkPoints>           points           = vtkSmartPointer<vtkPoints>::New();
   vtkSmartPointer<vtkDoubleArray>      data             = vtkDoubleArray::New();
@@ -164,59 +167,63 @@ void MeshName::save(const Mesh &mesh, const std::string &dataname) const
   data->SetNumberOfComponents(numComp);
 
   // Insert Points and Point Data
-  for (size_t i = 0; i < mesh.positions.size() - 1; i++) {
+  for (size_t i = 0; i < mesh.positions.size(); i++) {
     points->InsertNextPoint(mesh.positions[i][0], mesh.positions[i][1], mesh.positions[i][2]);
-    for (int j = 0; j < numComp; j++)
-      data->InsertNextTuple(&mesh.data[i * numComp + j]);
+    std::vector<double> pointData;
+    for (int j = 0; j < numComp; j++) {
+      pointData.push_back(mesh.data[i * numComp + j]);
+    }
+    data->InsertNextTuple(pointData.data());
   }
 
   unstructuredGrid->SetPoints(points);
   unstructuredGrid->GetPointData()->AddArray(data);
 
   // Connectivity Information
+  vtkSmartPointer<vtkCellArray> cellArray = vtkSmartPointer<vtkCellArray>::New();
+
+  std::vector<int> cellTypes;
+  cellTypes.reserve(mesh.quadrilaterals.size() + mesh.triangles.size() + mesh.edges.size());
 
   if (mesh.quadrilaterals.size() > 0) {
-    vtkSmartPointer<vtkCellArray> quadArray = vtkSmartPointer<vtkCellArray>::New();
+
     for (size_t i = 0; i < mesh.quadrilaterals.size(); i++) {
       vtkSmartPointer<vtkQuad> quadrilateral = vtkSmartPointer<vtkQuad>::New();
-
       quadrilateral->GetPointIds()->SetId(0, mesh.quadrilaterals[i][0]);
-      quadrilateral->GetPointIds()->SetId(0, mesh.quadrilaterals[i][1]);
-      quadrilateral->GetPointIds()->SetId(0, mesh.quadrilaterals[i][2]);
-      quadrilateral->GetPointIds()->SetId(0, mesh.quadrilaterals[i][3]);
+      quadrilateral->GetPointIds()->SetId(1, mesh.quadrilaterals[i][1]);
+      quadrilateral->GetPointIds()->SetId(2, mesh.quadrilaterals[i][2]);
+      quadrilateral->GetPointIds()->SetId(3, mesh.quadrilaterals[i][3]);
 
-      quadArray->InsertNextCell(quadrilateral);
+      cellArray->InsertNextCell(quadrilateral);
+      cellTypes.push_back(VTK_QUAD);
     }
-    unstructuredGrid->SetCells(VTK_QUAD, quadArray);
   }
 
   if (mesh.triangles.size() > 0) {
-    vtkSmartPointer<vtkCellArray> triArray = vtkSmartPointer<vtkCellArray>::New();
     for (size_t i = 0; i < mesh.triangles.size(); i++) {
       vtkSmartPointer<vtkTriangle> triangle = vtkSmartPointer<vtkTriangle>::New();
-
       triangle->GetPointIds()->SetId(0, mesh.triangles[i][0]);
       triangle->GetPointIds()->SetId(1, mesh.triangles[i][1]);
       triangle->GetPointIds()->SetId(2, mesh.triangles[i][2]);
 
-      triArray->InsertNextCell(triangle);
+      cellArray->InsertNextCell(triangle);
+      cellTypes.push_back(VTK_TRIANGLE);
     }
-    unstructuredGrid->SetCells(VTK_TRIANGLE, triArray);
   }
 
   if (mesh.edges.size() > 0) {
-    vtkSmartPointer<vtkCellArray> lineArray = vtkSmartPointer<vtkCellArray>::New();
     for (size_t i = 0; i < mesh.edges.size(); i++) {
       vtkSmartPointer<vtkLine> line = vtkSmartPointer<vtkLine>::New();
-
       line->GetPointIds()->SetId(0, mesh.edges[i][0]);
       line->GetPointIds()->SetId(1, mesh.edges[i][1]);
 
-      lineArray->InsertNextCell(line);
+      cellArray->InsertNextCell(line);
+      cellTypes.push_back(VTK_LINE);
     }
-    unstructuredGrid->SetCells(VTK_LINE, lineArray);
   }
 
+  unstructuredGrid->SetCells(cellTypes.data(), cellArray);
+
   // Write file
   vtkSmartPointer<vtkUnstructuredGridWriter> writer =
       vtkSmartPointer<vtkUnstructuredGridWriter>::New();
@@ -253,12 +260,22 @@ std::vector<MeshName> BaseName::findAll(const ExecutionContext &context) const
 
     // Check multiple timesteps
     std::vector<MeshName> meshNames;
-    for (int t = 0; true; ++t) {
+    for (int t = 1; true; ++t) {
       std::string stepMeshName = _bname + ".dt" + std::to_string(t);
       if (!fs::is_regular_file(stepMeshName + ".vtk"))
         break;
       meshNames.push_back(MeshName{stepMeshName});
     }
+    {
+      auto initMeshName = std::string{_bname + ".init" + ".vtk"};
+      if (fs::is_regular_file(initMeshName))
+        meshNames.push_back(MeshName{initMeshName});
+    }
+    {
+      auto finalMeshName = std::string{_bname + ".final" + ".vtk"};
+      if (fs::is_regular_file(finalMeshName))
+        meshNames.push_back(MeshName{finalMeshName});
+    }
     std::cerr << "Names: " << meshNames.size() << '\n';
     return meshNames;
   } else {
diff --git a/src/mesh.hpp b/src/mesh.hpp
index db716f3a..1ea84fb2 100644
--- a/src/mesh.hpp
+++ b/src/mesh.hpp
@@ -71,11 +71,15 @@ class BaseName {
   std::string _bname;
 };
 
+namespace aste {
+
+} // namespace aste
 struct Mesh {
   using Vertex   = std::array<double, 3>;
-  using Edge     = std::array<size_t, 2>;
-  using Triangle = std::array<size_t, 3>;
-  using Quad     = std::array<size_t, 4>;
+  using VID      = std::vector<Vertex>::size_type;
+  using Edge     = std::array<VID, 2>;
+  using Triangle = std::array<VID, 3>;
+  using Quad     = std::array<VID, 4>;
   std::vector<Vertex>   positions;
   std::vector<Edge>     edges;
   std::vector<Triangle> triangles;
diff --git a/src/tests/edges.conn.txt b/src/tests/edges.conn.txt
deleted file mode 100644
index 88e903c0..00000000
--- a/src/tests/edges.conn.txt
+++ /dev/null
@@ -1,4 +0,0 @@
-0 1
-0 2
-2 3
-1 3
diff --git a/src/tests/edges.txt b/src/tests/edges.txt
deleted file mode 100644
index a61c92ba..00000000
--- a/src/tests/edges.txt
+++ /dev/null
@@ -1,4 +0,0 @@
-0 0.0 0.0 0
-0 1.0 0.0 0
-0 0.0 1.0 0
-0 1.0 1.0 0
diff --git a/src/tests/read_test.cpp b/src/tests/read_test.cpp
new file mode 100644
index 00000000..578f0fba
--- /dev/null
+++ b/src/tests/read_test.cpp
@@ -0,0 +1,48 @@
+#include "testing.hpp"
+
+void readtest(const Case &current_case)
+{
+
+  auto read_test = aste::BaseName{current_case.fname}.with(aste::ExecutionContext());
+  auto mesh      = read_test.load(current_case.dim, current_case.dataname);
+
+  BOOST_TEST(mesh.positions.size() == 12);
+  BOOST_TEST(mesh.edges.size() == 2);
+  BOOST_TEST(mesh.quadrilaterals.size() == 2);
+  BOOST_TEST(mesh.triangles.size() == 4);
+  //std::cout << "Number of mesh elements are correctly loaded\n";
+
+  BOOST_TEST(mesh.edges[1][0] == 10);
+  BOOST_TEST(mesh.edges[1][1] == 11);
+  //std::cout << "Edges loaded correctly\n";
+  BOOST_TEST(mesh.triangles[0][0] == 0);
+  BOOST_TEST(mesh.triangles[0][1] == 1);
+  BOOST_TEST(mesh.triangles[0][2] == 3);
+  //std::cout << "Triangles loaded correctly\n";
+  BOOST_TEST(mesh.quadrilaterals[1][0] == 4);
+  BOOST_TEST(mesh.quadrilaterals[1][1] == 5);
+  BOOST_TEST(mesh.quadrilaterals[1][2] == 8);
+  BOOST_TEST(mesh.quadrilaterals[1][3] == 7);
+  //std::cout << "Quads loaded correctly\n";
+
+  switch (current_case.dim) {
+  case 1:
+    BOOST_TEST(mesh.data.size() == 12);
+    break;
+  case 2:
+    BOOST_TEST(mesh.data.size() == 24);
+    break;
+  case 3:
+    BOOST_TEST(mesh.data.size() == 36);
+    break;
+  }
+
+  std::vector<double> testdata;
+  testdata.resize(mesh.data.size());
+
+  std::iota(testdata.begin(), testdata.end(), 0);
+
+  for (size_t i = 0; i < mesh.data.size(); ++i) {
+    BOOST_TEST(mesh.data[i] == testdata[i]);
+  }
+}
\ No newline at end of file
diff --git a/src/tests/reference_scalars.vtk b/src/tests/reference_scalars.vtk
new file mode 100644
index 00000000..86475fa4
--- /dev/null
+++ b/src/tests/reference_scalars.vtk
@@ -0,0 +1,35 @@
+# vtk DataFile Version 4.1
+vtk output
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 12 float
+0 0 0 1 0 0 2 0 0 
+0 1 0 1 1 0 2 1 0 
+0 2 0 1 2 0 2 2 0 
+4 0 0 5 0 0 6 0 0 
+
+CELLS 8 32
+4 3 4 7 6 
+4 4 5 8 7 
+3 0 1 3 
+3 1 4 3 
+3 1 2 4 
+3 2 4 5 
+2 9 10 
+2 10 11 
+
+CELL_TYPES 8
+9
+9
+5
+5
+5
+5
+3
+3
+
+POINT_DATA 12
+FIELD FieldData 1
+Scalars 1 12 double
+0 1 2 3 4 5 6 7 8 
+9 10 11 
diff --git a/src/tests/reference_vector2d.vtk b/src/tests/reference_vector2d.vtk
new file mode 100644
index 00000000..09085d4b
--- /dev/null
+++ b/src/tests/reference_vector2d.vtk
@@ -0,0 +1,36 @@
+# vtk DataFile Version 4.1
+vtk output
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 12 float
+0 0 0 1 0 0 2 0 0 
+0 1 0 1 1 0 2 1 0 
+0 2 0 1 2 0 2 2 0 
+4 0 0 5 0 0 6 0 0 
+
+CELLS 8 32
+4 3 4 7 6 
+4 4 5 8 7 
+3 0 1 3 
+3 1 4 3 
+3 1 2 4 
+3 2 4 5 
+2 9 10 
+2 10 11 
+
+CELL_TYPES 8
+9
+9
+5
+5
+5
+5
+3
+3
+
+POINT_DATA 12
+FIELD FieldData 1
+Vector2D 2 12 double
+0 1 2 3 4 5 6 7 8 
+9 10 11 12 13 14 15 16 17 
+18 19 20 21 22 23 
diff --git a/src/tests/reference_vector3d.vtk b/src/tests/reference_vector3d.vtk
new file mode 100644
index 00000000..34e65562
--- /dev/null
+++ b/src/tests/reference_vector3d.vtk
@@ -0,0 +1,38 @@
+# vtk DataFile Version 4.1
+vtk output
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 12 float
+0 0 0 1 0 0 2 0 0 
+0 1 0 1 1 0 2 1 0 
+0 2 0 1 2 0 2 2 0 
+4 0 0 5 0 0 6 0 0 
+
+CELLS 8 32
+4 3 4 7 6 
+4 4 5 8 7 
+3 0 1 3 
+3 1 4 3 
+3 1 2 4 
+3 2 4 5 
+2 9 10 
+2 10 11 
+
+CELL_TYPES 8
+9
+9
+5
+5
+5
+5
+3
+3
+
+POINT_DATA 12
+FIELD FieldData 1
+Vector3D 3 12 double
+0 1 2 3 4 5 6 7 8 
+9 10 11 12 13 14 15 16 17 
+18 19 20 21 22 23 24 25 26 
+27 28 29 30 31 32 33 34 35 
+
diff --git a/src/tests/simple.txt b/src/tests/simple.txt
deleted file mode 100644
index e66c2bff..00000000
--- a/src/tests/simple.txt
+++ /dev/null
@@ -1,4 +0,0 @@
-0 0.0 0.0 0
-0 1.0 0.0 0.1
-0 0.0 1.0 -0.2
-0 1.0 1.0 5.2
diff --git a/src/tests/testing.cpp b/src/tests/testing.cpp
new file mode 100644
index 00000000..37ba81e5
--- /dev/null
+++ b/src/tests/testing.cpp
@@ -0,0 +1,37 @@
+#define BOOST_TEST_MODULE ASTE
+#define BOOST_TEST_DYN_LINK
+#include "testing.hpp"
+
+BOOST_AUTO_TEST_SUITE(read_write_tests)
+
+BOOST_AUTO_TEST_CASE(read_test_scalar)
+{
+  readtest(Case{"../src/tests/reference_scalars", "Scalars", 1});
+}
+
+BOOST_AUTO_TEST_CASE(read_test_2dvector)
+{
+  readtest(Case{"../src/tests/reference_vector2d", "Vector2D", 2});
+}
+
+BOOST_AUTO_TEST_CASE(read_test_3dvector)
+{
+  readtest(Case{"../src/tests/reference_vector3d", "Vector3D", 3});
+}
+
+BOOST_AUTO_TEST_CASE(write_test_scalar)
+{
+  writetest(Case{"./scalar_write", "Scalars", 1});
+}
+
+BOOST_AUTO_TEST_CASE(write_test_2dvector)
+{
+  writetest(Case{"./vector2d_write", "Vector2D", 2});
+}
+
+BOOST_AUTO_TEST_CASE(write_test_3dvector)
+{
+  writetest(Case{"./vector3d_write", "Vector3D", 3});
+}
+
+BOOST_AUTO_TEST_SUITE_END() //read_write_tests
\ No newline at end of file
diff --git a/src/tests/testing.hpp b/src/tests/testing.hpp
new file mode 100644
index 00000000..60c3aa33
--- /dev/null
+++ b/src/tests/testing.hpp
@@ -0,0 +1,13 @@
+#include <boost/test/unit_test.hpp>
+#include <mesh.hpp>
+#include <numeric>
+#include <string>
+
+struct Case {
+  std::string fname{};
+  std::string dataname{};
+  int         dim{};
+};
+
+void writetest(const Case &current_case);
+void readtest(const Case &current_case);
diff --git a/src/tests/triangles.conn.txt b/src/tests/triangles.conn.txt
deleted file mode 100644
index 249f4cd1..00000000
--- a/src/tests/triangles.conn.txt
+++ /dev/null
@@ -1,2 +0,0 @@
-0 1 2
-1 3 2
diff --git a/src/tests/triangles.txt b/src/tests/triangles.txt
deleted file mode 100644
index a61c92ba..00000000
--- a/src/tests/triangles.txt
+++ /dev/null
@@ -1,4 +0,0 @@
-0 0.0 0.0 0
-0 1.0 0.0 0
-0 0.0 1.0 0
-0 1.0 1.0 0
diff --git a/src/tests/write_test.cpp b/src/tests/write_test.cpp
new file mode 100644
index 00000000..1d051b74
--- /dev/null
+++ b/src/tests/write_test.cpp
@@ -0,0 +1,77 @@
+#include "testing.hpp"
+
+void writetest(const Case &current_case)
+{
+  using VID = aste::Mesh::VID;
+
+  aste::Mesh testMesh;
+  testMesh.positions.reserve(12);
+
+  // Points for Quads and Tri Elements
+  for (double y = 0; y < 3.0; ++y) {
+    for (double x = 0; x < 3.0; ++x) {
+      const std::array<double, 3> pos{x, y, 0.0};
+      testMesh.positions.push_back(pos);
+    }
+  }
+  //Points for Line Elements
+  for (double x = 4; x < 7.0; ++x) {
+    const std::array<double, 3> pos{x, 0.0, 0.0};
+    testMesh.positions.push_back(pos);
+  }
+
+  //Create Lines
+  std::array<VID, 2> line1{9, 10};
+  std::array<VID, 2> line2{10, 11};
+  testMesh.edges.push_back(line1);
+  testMesh.edges.push_back(line2);
+
+  // Create Triangles
+  std::array<VID, 3> tri1{0, 1, 3};
+  std::array<VID, 3> tri2{1, 4, 3};
+  std::array<VID, 3> tri3{1, 2, 4};
+  std::array<VID, 3> tri4{2, 4, 5};
+  testMesh.triangles.push_back(tri1);
+  testMesh.triangles.push_back(tri2);
+  testMesh.triangles.push_back(tri3);
+  testMesh.triangles.push_back(tri4);
+
+  //Create Quad Elements
+  std::array<VID, 4> quad1{3, 4, 7, 6};
+  std::array<VID, 4> quad2{4, 5, 8, 7};
+  testMesh.quadrilaterals.push_back(quad1);
+  testMesh.quadrilaterals.push_back(quad2);
+
+  testMesh.data.resize(testMesh.positions.size() * current_case.dim);
+  std::iota(testMesh.data.begin(), testMesh.data.end(), 0);
+  auto scalar_test = aste::BaseName{current_case.fname}.with(aste::ExecutionContext());
+  scalar_test.save(testMesh, current_case.dataname);
+
+  // Read written data and compare with created data
+  auto read       = aste::BaseName{current_case.fname}.with(aste::ExecutionContext());
+  auto loadedMesh = read.load(current_case.dim, current_case.dataname);
+
+  //Check Elements are correctly written
+  BOOST_TEST(loadedMesh.positions.size() == testMesh.positions.size());
+  BOOST_TEST(loadedMesh.edges.size() == testMesh.edges.size());
+  BOOST_TEST(loadedMesh.quadrilaterals.size() == testMesh.quadrilaterals.size());
+  BOOST_TEST(loadedMesh.triangles.size() == testMesh.triangles.size());
+  //Check Edges
+  BOOST_TEST(loadedMesh.edges[1][0] == testMesh.edges[1][0]);
+  BOOST_TEST(loadedMesh.edges[1][1] == testMesh.edges[1][1]);
+  //Check Triangles
+  BOOST_TEST(loadedMesh.triangles[0][0] == testMesh.triangles[0][0]);
+  BOOST_TEST(loadedMesh.triangles[0][1] == testMesh.triangles[0][1]);
+  BOOST_TEST(loadedMesh.triangles[0][2] == testMesh.triangles[0][2]);
+  //Check Quads
+  BOOST_TEST(loadedMesh.quadrilaterals[1][0] == testMesh.quadrilaterals[1][0]);
+  BOOST_TEST(loadedMesh.quadrilaterals[1][1] == testMesh.quadrilaterals[1][1]);
+  BOOST_TEST(loadedMesh.quadrilaterals[1][2] == testMesh.quadrilaterals[1][2]);
+  BOOST_TEST(loadedMesh.quadrilaterals[1][3] == testMesh.quadrilaterals[1][3]);
+  //Check Datasize
+  BOOST_TEST(loadedMesh.data.size() == testMesh.data.size());
+  //Check Data Values
+  for (size_t i = 0; i < loadedMesh.data.size(); ++i) {
+    BOOST_TEST(loadedMesh.data[i] == testMesh.data[i]);
+  }
+}