diff --git a/test/freeflow/rans/CMakeLists.txt b/test/freeflow/rans/CMakeLists.txt
index 58bf91a55d99f54b421d38656d45a869e783f49b..0b3b205ce7aa702287aaeadd9182b8d4d3d3a1df 100644
--- a/test/freeflow/rans/CMakeLists.txt
+++ b/test/freeflow/rans/CMakeLists.txt
@@ -1,4 +1,5 @@
 add_input_file_links()
+dune_symlink_to_source_files(FILES laufer_re50000_u+y+.csv)
 
 dune_add_test(NAME test_pipe_laufer
               SOURCES test_pipe_laufer.cc
diff --git a/test/freeflow/rans/laufer_re50000_u+y+.csv b/test/freeflow/rans/laufer_re50000_u+y+.csv
new file mode 100644
index 0000000000000000000000000000000000000000..547f79ab0ced283e16dfe40aa6867440cd4b9446
--- /dev/null
+++ b/test/freeflow/rans/laufer_re50000_u+y+.csv
@@ -0,0 +1,46 @@
+###########
+# Dimensionless velocity profile for Re=50,000
+###########
+# Original source:
+# Laufer, J.
+# The structure of turbulence in fully developed pipe flow
+# NACA Report, 1954, 1174, 417-434
+###########
+# Data taken from:
+# Truckenbrodt, E.
+# Grundlagen und elementare Strömungsvorgänge dichtebeständiger Fluide
+# Fluidmechanik, Springer, 2008, XXVII, 364 , doi: 10.1007/978-3-540-79018-1
+###########
+#Y+,[-],U+,[-]
+3.12,3.243
+3.572,3.649
+4.289,4.135
+4.911,5.068
+5.327,5.27
+5.436,4.662
+5.817,5.473
+6.439,6.162
+7.373,7.095
+7.731,7.5
+8.674,7.905
+8.793,8.311
+9.345,7.5
+10.136,8.716
+10.847,9.203
+13.29,10.014
+16.284,11.351
+17.425,11.635
+26.161,13.054
+43.475,14.676
+66.608,15.486
+100,16.703
+166.184,17.838
+233.156,18.932
+295.521,19.541
+400.812,20.554
+508.022,21.568
+609.95,22.378
+712.756,22.986
+816.14,23.311
+915.724,23.716
+1034.441,23.797
diff --git a/test/freeflow/rans/test_pipe_laufer.cc b/test/freeflow/rans/test_pipe_laufer.cc
index 75de98f76841506816c4b47180a30068da3274d0..a166c766e6901d6a8b6f91abbc03850929d736ca 100644
--- a/test/freeflow/rans/test_pipe_laufer.cc
+++ b/test/freeflow/rans/test_pipe_laufer.cc
@@ -51,6 +51,7 @@
 
 #include <dumux/discretization/methods.hh>
 
+#include <dumux/io/gnuplotinterface.hh>
 #include <dumux/io/staggeredvtkoutputmodule.hh>
 
 #include <dumux/freeflow/navierstokes/staggered/fluxoverplane.hh>
@@ -203,6 +204,45 @@ int main(int argc, char** argv) try
     // finalize, print dumux message to say goodbye
     ////////////////////////////////////////////////////////////
 
+#if HAVE_PVPYTHON
+    bool shouldPlot = getParam<bool>("Output.PlotVelocityProfiles", false);
+    if (shouldPlot)
+    {
+        char fileName[255];
+        std::string fileNameFormat = "%s-%05d";
+        sprintf(fileName, fileNameFormat.c_str(), problem->name().c_str(), timeLoop->timeStepIndex());
+        std::cout << fileName << std::endl;
+        std::string vtuFileName = std::string(fileName) + ".vtu";
+        std::string script = std::string(DUMUX_SOURCE_DIR) + "/bin/postprocessing/extractlinedata.py";
+        std::string syscom;
+
+        // execute the pvpython script
+        std::string command = std::string(PVPYTHON_EXECUTABLE) + " " + script
+                              + " -f " + vtuFileName
+                              + " -v 1"
+                              + " -r 10000";
+        syscom =  command + " -p1 8.0 0.0 0.0"
+                          + " -p2 8.0 0.2469 0.0"
+                          + " -of " + std::string(fileName) + "\n";
+        system(syscom.c_str());
+
+
+        char gnuplotFileName[255];
+        Dumux::GnuplotInterface<double> gnuplot;
+        sprintf(gnuplotFileName, fileNameFormat.c_str(), "velProfiles", timeLoop->timeStepIndex());
+        gnuplot.setOpenPlotWindow(shouldPlot);
+        gnuplot.setDatafileSeparator(',');
+        gnuplot.resetPlot();
+        gnuplot.setXlabel("y^+ [-]");
+        gnuplot.setYlabel("u_+ [-]");
+        gnuplot.setOption("set log x");
+        gnuplot.setOption("set xrange [1:3000]");
+        gnuplot.addFileToPlot("laufer_re50000_u+y+.csv", "u 1:2 w p t 'Laufer1954a, Re=50000'");
+        gnuplot.addFileToPlot(std::string(fileName) + ".csv", "u 13:14 w l");
+        gnuplot.plot(std::string(gnuplotFileName));
+    }
+#endif
+
     // print dumux end message
     if (mpiHelper.rank() == 0)
     {
diff --git a/test/freeflow/rans/test_pipe_laufer.input b/test/freeflow/rans/test_pipe_laufer.input
index 0909df9c1c9ac102f5e10d40d00bef6f0df1526f..6038fce674615827e5f784a68e6312dcaaed8124 100644
--- a/test/freeflow/rans/test_pipe_laufer.input
+++ b/test/freeflow/rans/test_pipe_laufer.input
@@ -10,6 +10,9 @@ Cells0 = 25
 Cells1 = 25 25
 Grading1 = 1.2 -1.2
 
+[Output]
+PlotVelocityProfiles = true
+
 [Problem]
 Name = pipe_laufer_zeroeq
 InletVelocity = 2.5 # [m/s]