diff --git a/src/eam_utils.c b/src/eam_utils.c
index 6437d39..9a872d9 100644
--- a/src/eam_utils.c
+++ b/src/eam_utils.c
@@ -44,11 +44,13 @@ void initEam(Eam* eam, Parameter* param) {
 
 void coeff(Eam* eam, Parameter* param) {
     read_file(&eam->file, param->input_file);
+    param->mass = eam->file.mass;
     param->cutforce = eam->file.cut;
     param->cutneigh = param->cutforce + 1.0;
     param->temp = 600.0;
     param->dt = 0.001;
     param->rho = 0.07041125;
+    param->dtforce = 0.5 * param->dt / param->mass;
 }
 
 void init_style(Eam* eam, Parameter* param) {
diff --git a/src/includes/parameter.h b/src/includes/parameter.h
index ac6c737..f23e0e8 100644
--- a/src/includes/parameter.h
+++ b/src/includes/parameter.h
@@ -35,6 +35,7 @@
 typedef struct {
     int force_field;
     char* input_file;
+    char* vtk_file;
     MD_FLOAT epsilon;
     MD_FLOAT sigma6;
     MD_FLOAT temp;
diff --git a/src/includes/vtk.h b/src/includes/vtk.h
new file mode 100644
index 0000000..0f03a74
--- /dev/null
+++ b/src/includes/vtk.h
@@ -0,0 +1,28 @@
+/*
+ * =======================================================================================
+ *
+ *   Author:   Jan Eitzinger (je), jan.eitzinger@fau.de
+ *   Copyright (c) 2020 RRZE, University Erlangen-Nuremberg
+ *
+ *   This file is part of MD-Bench.
+ *
+ *   MD-Bench 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.
+ *
+ *   MD-Bench 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 Lesser General Public License for more
+ *   details.
+ *
+ *   You should have received a copy of the GNU Lesser General Public License along
+ *   with MD-Bench.  If not, see <https://www.gnu.org/licenses/>.
+ * =======================================================================================
+ */
+#include <atom.h>
+
+#ifndef __VTK_H_
+#define __VTK_H_
+extern int write_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep);
+#endif
diff --git a/src/main.c b/src/main.c
index 362dee9..6dbb648 100644
--- a/src/main.c
+++ b/src/main.c
@@ -40,6 +40,7 @@
 #include <pbc.h>
 #include <timers.h>
 #include <eam.h>
+#include <vtk.h>
 
 #define HLINE "----------------------------------------------------------------------------\n"
 
@@ -50,6 +51,7 @@ extern double computeForceEam(Eam* eam, Atom *atom, Neighbor *neighbor, Stats *s
 void init(Parameter *param)
 {
     param->input_file = NULL;
+    param->vtk_file = NULL;
     param->force_field = FF_LJ;
     param->epsilon = 1.0;
     param->sigma6 = 1.0;
@@ -77,6 +79,7 @@ double setup(
         Neighbor *neighbor,
         Stats *stats)
 {
+    if(param->force_field == FF_EAM) { initEam(eam, param); }
     double S, E;
     param->lattice = pow((4.0 / param->rho), (1.0 / 3.0));
     param->xprd = param->nx * param->lattice;
@@ -84,7 +87,6 @@ double setup(
     param->zprd = param->nz * param->lattice;
 
     S = getTimeStamp();
-    if(param->force_field == FF_EAM) { initEam(eam, param); }
     initAtom(atom);
     initNeighbor(neighbor, param);
     initPbc();
@@ -232,6 +234,11 @@ int main(int argc, char** argv)
             param.proc_freq = atof(argv[++i]);
             continue;
         }
+        if((strcmp(argv[i], "--vtk") == 0))
+        {
+            param.vtk_file = strdup(argv[++i]);
+            continue;
+        }
         if((strcmp(argv[i], "-h") == 0) || (strcmp(argv[i], "--help") == 0))
         {
             printf("MD Bench: A minimalistic re-implementation of miniMD\n");
@@ -241,6 +248,7 @@ int main(int argc, char** argv)
             printf("-n / --nsteps <int>:  set number of timesteps for simulation\n");
             printf("-nx/-ny/-nz <int>:    set linear dimension of systembox in x/y/z direction\n");
             printf("--freq <real>:        processor frequency (GHz)\n");
+            printf("--vtk <string>:       VTK file for visualization\n");
             printf(HLINE);
             exit(EXIT_SUCCESS);
         }
@@ -262,6 +270,10 @@ int main(int argc, char** argv)
     timer[NEIGH] = 0.0;
     timer[TOTAL] = getTimeStamp();
 
+    if(param.vtk_file != NULL) {
+        write_atoms_to_vtk_file(param.vtk_file, &atom, 0);
+    }
+
     for(int n = 0; n < param.ntimes; n++) {
         initialIntegrate(&param, &atom);
 
@@ -285,6 +297,10 @@ int main(int argc, char** argv)
         if(!((n + 1) % param.nstat) && (n+1) < param.ntimes) {
             computeThermo(n + 1, &param, &atom);
         }
+
+        if(param.vtk_file != NULL) {
+            write_atoms_to_vtk_file(param.vtk_file, &atom, n + 1);
+        }
     }
 
     timer[TOTAL] = getTimeStamp() - timer[TOTAL];
diff --git a/src/thermo.c b/src/thermo.c
index 2361579..70b897d 100644
--- a/src/thermo.c
+++ b/src/thermo.c
@@ -31,7 +31,7 @@ static MD_FLOAT *tmparr;
 static MD_FLOAT *engarr;
 static MD_FLOAT *prsarr;
 static MD_FLOAT mvv2e;
-static int dof_boltz;
+static MD_FLOAT dof_boltz;
 static MD_FLOAT t_scale;
 static MD_FLOAT p_scale;
 static MD_FLOAT e_scale;
@@ -81,7 +81,6 @@ void computeThermo(int iflag, Parameter *param, Atom *atom)
 
     t = t * t_scale;
     p = (t * dof_boltz) * p_scale;
-
     int istep = iflag;
 
     if(iflag == -1){
diff --git a/src/vtk.c b/src/vtk.c
new file mode 100644
index 0000000..b3ff0e6
--- /dev/null
+++ b/src/vtk.c
@@ -0,0 +1,44 @@
+#include <stdio.h>
+#include <stdlib.h>
+
+#include <atom.h>
+
+int write_atoms_to_vtk_file(const char* filename, Atom* atom, int timestep) {
+    char timestep_filename[128];
+    snprintf(timestep_filename, sizeof timestep_filename, "%s_%d.vtk", filename, timestep);
+    FILE* fp = fopen(timestep_filename, "wb");
+
+    if(fp == NULL) {
+        fprintf(stderr, "Could not open VTK file for writing!\n");
+        return -1;
+    }
+
+    fprintf(fp, "# vtk DataFile Version 2.0\n");
+    fprintf(fp, "Particle data\n");
+    fprintf(fp, "ASCII\n");
+    fprintf(fp, "DATASET UNSTRUCTURED_GRID\n");
+    fprintf(fp, "POINTS %d double\n", atom->Nlocal);
+    for(int i = 0; i < atom->Nlocal; ++i) {
+        fprintf(fp, "%.4f %.4f %.4f\n", atom_x(i), atom_y(i), atom_z(i));
+    }
+    fprintf(fp, "\n\n");
+    fprintf(fp, "CELLS %d %d\n", atom->Nlocal, atom->Nlocal * 2);
+    for(int i = 0; i < atom->Nlocal; ++i) {
+        fprintf(fp, "1 %d\n", i);
+    }
+    fprintf(fp, "\n\n");
+    fprintf(fp, "CELL_TYPES %d\n", atom->Nlocal);
+    for(int i = 0; i < atom->Nlocal; ++i) {
+        fprintf(fp, "1\n");
+    }
+    fprintf(fp, "\n\n");
+    fprintf(fp, "POINT_DATA %d\n", atom->Nlocal);
+    fprintf(fp, "SCALARS mass double\n");
+    fprintf(fp, "LOOKUP_TABLE default\n");
+    for(int i = 0; i < atom->Nlocal; i++) {
+        fprintf(fp, "1.0\n");
+    }
+    fprintf(fp, "\n\n");
+    fclose(fp);
+    return 0;
+}