diff --git a/Applications/Utils/MeshEdit/CMakeLists.txt b/Applications/Utils/MeshEdit/CMakeLists.txt
index ea9bef3f1b6f0418ea8955fbc5b259d01a380505..92fd7a7580d14568d9381bb88f381466664e675b 100644
--- a/Applications/Utils/MeshEdit/CMakeLists.txt
+++ b/Applications/Utils/MeshEdit/CMakeLists.txt
@@ -17,7 +17,8 @@ set(TOOLS
ResetPropertiesInPolygonalRegion
reviseMesh
swapNodeCoordinateAxes
- UnityPreprocessing)
+ UnityPreprocessing
+ Vtu2Grid)
foreach(TOOL ${TOOLS})
add_executable(${TOOL} ${TOOL}.cpp)
target_link_libraries(${TOOL} GitInfoLib MeshLib)
diff --git a/Applications/Utils/MeshEdit/Vtu2Grid.cpp b/Applications/Utils/MeshEdit/Vtu2Grid.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..b23adf93a73b736d4d59efca3be0799e86587746
--- /dev/null
+++ b/Applications/Utils/MeshEdit/Vtu2Grid.cpp
@@ -0,0 +1,227 @@
+/**
+ * \file
+ *
+ * @copyright
+ * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/LICENSE.txt
+ */
+
+#include <tclap/CmdLine.h>
+
+#include "Applications/ApplicationsLib/LogogSetup.h"
+#include "InfoLib/GitInfo.h"
+
+#include "MeshLib/IO/writeMeshToFile.h"
+#include "MeshLib/Mesh.h"
+#include "MeshLib/MeshEditing/RemoveMeshComponents.h"
+#include "MeshLib/MeshGenerators/MeshGenerator.h"
+#include "MeshLib/MeshSearch/ElementSearch.h"
+
+#include <vtkAbstractArray.h>
+#include <vtkCellData.h>
+#include <vtkCellLocator.h>
+#include <vtkDataArray.h>
+#include <vtkDoubleArray.h>
+#include <vtkIntArray.h>
+#include <vtkSmartPointer.h>
+#include <vtkUnstructuredGrid.h>
+#include <vtkXMLUnstructuredGridReader.h>
+
+std::string const cell_id_name = "CellIds";
+
+std::array<std::size_t, 3> getDimensions(MathLib::Point3d const& min,
+ MathLib::Point3d const& max,
+ std::array<double, 3> const& cellsize)
+{
+ return
+ {
+ static_cast<std::size_t>(std::ceil((max[0] - min[0]) / cellsize[0])),
+ static_cast<std::size_t>(std::ceil((max[1] - min[1]) / cellsize[1])),
+ static_cast<std::size_t>(std::ceil((max[2] - min[2]) / cellsize[2]))
+ };
+}
+
+std::vector<int> assignCellIds(vtkSmartPointer<vtkUnstructuredGrid> const& mesh,
+ MathLib::Point3d const& min,
+ std::array<std::size_t, 3> const& dims,
+ std::array<double, 3> const& cellsize)
+{
+ vtkSmartPointer<vtkCellLocator> locator =
+ vtkSmartPointer<vtkCellLocator>::New();
+ locator->SetDataSet(mesh);
+ locator->Update();
+
+ std::vector<int> cell_ids;
+ cell_ids.reserve(dims[0] * dims[1] * dims[2]);
+ std::array<double, 3> const grid_max = {min[0] + dims[0] * cellsize[0],
+ min[1] + dims[1] * cellsize[1],
+ min[2] + dims[2] * cellsize[2]};
+ for (double k = min[2] + (cellsize[2] / 2.0); k < grid_max[2];
+ k += cellsize[2])
+ {
+ for (double j = min[1] + (cellsize[1] / 2.0); j < grid_max[1];
+ j += cellsize[1])
+ {
+ for (double i = min[0] + (cellsize[0] / 2.0); i < grid_max[0];
+ i += cellsize[0])
+ {
+ double pnt[3] = {i, j, k};
+ cell_ids.push_back(static_cast<int>(locator->FindCell(pnt)));
+ }
+ }
+ }
+ return cell_ids;
+}
+
+bool removeUnusedGridCells(vtkSmartPointer<vtkUnstructuredGrid> const& mesh,
+ std::unique_ptr<MeshLib::Mesh>& grid)
+{
+ MeshLib::ElementSearch search(*grid);
+ std::size_t const n_elems_marked = search.searchByPropertyValueRange<int>(
+ cell_id_name, 0, static_cast<int>(mesh->GetNumberOfCells()), true);
+
+ if (n_elems_marked == grid->getNumberOfElements())
+ {
+ ERR ("No valid elements found. Aborting...");
+ return false;
+ }
+
+ if (n_elems_marked)
+ {
+ grid.reset(MeshLib::removeElements(
+ *grid, search.getSearchedElementIDs(), "trimmed_grid"));
+ }
+ return true;
+}
+
+template <typename T, typename VTK_TYPE>
+void mapArray(MeshLib::Mesh& grid, VTK_TYPE vtk_arr, std::string const& name)
+{
+ MeshLib::PropertyVector<int> const& cell_ids =
+ *grid.getProperties().getPropertyVector<int>(
+ cell_id_name, MeshLib::MeshItemType::Cell, 1);
+ std::vector<T>& arr = *grid.getProperties().createNewPropertyVector<T>(
+ name, MeshLib::MeshItemType::Cell, vtk_arr->GetNumberOfComponents());
+ std::size_t const n_elems = cell_ids.size();
+ arr.resize(n_elems);
+ for (std::size_t j = 0; j < n_elems; ++j)
+ arr[j] = vtk_arr->GetValue(cell_ids[j]);
+}
+
+void mapMeshArraysOntoGrid(vtkSmartPointer<vtkUnstructuredGrid> const& mesh,
+ std::unique_ptr<MeshLib::Mesh>& grid)
+{
+ vtkSmartPointer<vtkCellData> const cell_data = mesh->GetCellData();
+ for (int i = 0; i < cell_data->GetNumberOfArrays(); ++i)
+ {
+ std::string const& name = cell_data->GetArrayName(i);
+ vtkSmartPointer<vtkDoubleArray> const dbl_arr =
+ dynamic_cast<vtkDoubleArray*>(cell_data->GetArray(name.c_str()));
+ if (dbl_arr)
+ {
+ mapArray<double, vtkSmartPointer<vtkDoubleArray>>(*grid, dbl_arr, name);
+ continue;
+ }
+ vtkSmartPointer<vtkIntArray> const int_arr =
+ dynamic_cast<vtkIntArray*>(cell_data->GetArray(name.c_str()));
+ if (int_arr)
+ {
+ mapArray<int, vtkSmartPointer<vtkIntArray>>(*grid, int_arr, name);
+ continue;
+ }
+ WARN("Ignoring array '%s', array type not implemented...", name.c_str());
+ }
+}
+
+int main (int argc, char* argv[])
+{
+ ApplicationsLib::LogogSetup logog_setup;
+
+ TCLAP::CmdLine cmd(
+ "Reads a 3D unstructured mesh and samples it onto a structured grid of "
+ "the same extent. Cell properties are mapped onto the grid (sampled at "
+ "the centre-points of each cube), node properties are ignored. Note, "
+ "that a large cube size may result in an undersampling of the original "
+ "mesh structure.\nCube sizes are defines by x/y/z-parameters. For "
+ "equilateral cubes, only the x-parameter needs to be set.\n\n"
+ "OpenGeoSys-6 software, version " +
+ GitInfoLib::GitInfo::ogs_version +
+ ".\n"
+ "Copyright (c) 2012-2019, OpenGeoSys Community "
+ "(http://www.opengeosys.org)",
+ ' ', GitInfoLib::GitInfo::ogs_version);
+
+ TCLAP::ValueArg<double> z_arg("z", "cellsize-z",
+ "edge length of cubes in z-direction (depth)",
+ false, 1000, "floating point number");
+ cmd.add(z_arg);
+
+ TCLAP::ValueArg<double> y_arg(
+ "y", "cellsize-y", "edge length of cubes in y-direction (latitude)",
+ false, 1000, "floating point number");
+ cmd.add(y_arg);
+
+ TCLAP::ValueArg<double> x_arg(
+ "x", "cellsize-x",
+ "edge length of cubes in x-direction (longitude) or all directions, if "
+ "y and z are not set",
+ true, 1000, "floating point number");
+ cmd.add(x_arg);
+
+ TCLAP::ValueArg<std::string> output_arg(
+ "o", "output", "the output grid (*.vtu)", true, "", "output.vtu");
+ cmd.add(output_arg);
+
+ TCLAP::ValueArg<std::string> input_arg("i", "input",
+ "the 3D input mesh (*.vtu, *.msh)",
+ true, "", "input.vtu");
+ cmd.add(input_arg);
+ cmd.parse(argc, argv);
+
+ if ((y_arg.isSet() && !z_arg.isSet()) ||
+ ((!y_arg.isSet() && z_arg.isSet())))
+ {
+ ERR("For equilateral cubes, only x needs to be set. For unequal "
+ "cuboids, all three edge lengths (x/y/z) need to be specified.")
+ return -1;
+ }
+
+ double const x_size = x_arg.getValue();
+ double const y_size = (y_arg.isSet()) ? y_arg.getValue() : x_arg.getValue();
+ double const z_size = (z_arg.isSet()) ? z_arg.getValue() : x_arg.getValue();
+ std::array<double, 3> const cellsize = { x_size, y_size, z_size };
+
+ vtkSmartPointer<vtkXMLUnstructuredGridReader> reader =
+ vtkSmartPointer<vtkXMLUnstructuredGridReader>::New();
+ reader->SetFileName(input_arg.getValue().c_str());
+ reader->Update();
+ vtkSmartPointer<vtkUnstructuredGrid> mesh = reader->GetOutput();
+
+ double* const bounds = mesh->GetBounds();
+ MathLib::Point3d const min(std::array<double, 3>{bounds[0], bounds[2], bounds[4]});
+ MathLib::Point3d const max(std::array<double, 3>{bounds[1], bounds[3], bounds[5]});
+ std::array<std::size_t, 3> const dims = getDimensions(min, max, cellsize);
+ std::unique_ptr<MeshLib::Mesh> grid(
+ MeshLib::MeshGenerator::generateRegularHexMesh(
+ dims[0], dims[1], dims[2], cellsize[0], cellsize[1], cellsize[2],
+ min, "grid"));
+
+
+ std::vector<int> const tmp_ids = assignCellIds(mesh, min, dims, cellsize);
+ std::vector<int>& cell_ids =
+ *grid->getProperties().createNewPropertyVector<int>(
+ cell_id_name, MeshLib::MeshItemType::Cell, 1);
+ std::copy(tmp_ids.cbegin(), tmp_ids.cend(), std::back_inserter(cell_ids));
+
+ if (!removeUnusedGridCells(mesh, grid))
+ return EXIT_FAILURE;
+
+ mapMeshArraysOntoGrid(mesh, grid);
+
+ if (MeshLib::IO::writeMeshToFile(*grid, output_arg.getValue()) != 0)
+ return EXIT_FAILURE;
+
+ return EXIT_SUCCESS;
+}
diff --git a/Applications/Utils/Tests.cmake b/Applications/Utils/Tests.cmake
index 2c182e5519bbf7146e9d7de27d1cf3d1abd59650..6fffe0a5ed3f1d529abea405ff4131b99f866778 100644
--- a/Applications/Utils/Tests.cmake
+++ b/Applications/Utils/Tests.cmake
@@ -274,3 +274,14 @@ AddTest(
TESTER diff
DIFF_DATA output0.asc
)
+
+AddTest(
+ NAME Vtu2Grid_Test
+ PATH FileIO/
+ EXECUTABLE Vtu2Grid
+ EXECUTABLE_ARGS -i AmmerSubsurfaceCoarse.vtu -o ${Data_BINARY_DIR}/FileIO/AmmerGridOutput.vtu -x 200 -y 200 -z 20
+ REQUIREMENTS NOT OGS_USE_MPI
+ TESTER vtkdiff
+ DIFF_DATA
+ AmmerSubsurfaceGrid.vtu AmmerGridOutput.vtu MaterialIDs MaterialIDs 0 0
+)
diff --git a/Tests/Data/FileIO/AmmerSubsurfaceGrid.vtu b/Tests/Data/FileIO/AmmerSubsurfaceGrid.vtu
new file mode 100644
index 0000000000000000000000000000000000000000..81ba1fe210978559bf57a8130a60ea0158807d63
--- /dev/null
+++ b/Tests/Data/FileIO/AmmerSubsurfaceGrid.vtu
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:f9f060bd48bc8e2248ad853e29564b588427b4b2d0ac6ab10ac193102f7ef082
+size 3762510
diff --git a/web/content/docs/tools/meshing/vtu2grid/index.pandoc b/web/content/docs/tools/meshing/vtu2grid/index.pandoc
new file mode 100644
index 0000000000000000000000000000000000000000..e2794c9aecf9f1c631a1dbb29f3ceaa1a27e30c4
--- /dev/null
+++ b/web/content/docs/tools/meshing/vtu2grid/index.pandoc
@@ -0,0 +1,79 @@
++++
+date = "2019-11-27T00:00:00+01:00"
+title = "Vtu2Grid"
+author = "Karsten Rink"
+
+[menu]
+ [menu.tools]
+ parent = "meshing"
++++
+
+## Introduction
+
+This utility rasterises an existing 3D mesh at a given resolution. The result is a (quasi-)structured grid (consisting of hexahedral elements) with the same extent as the input mesh. Cell properties are mapped onto the grid (sampled at centre-points of each cube), while node properties are ignored. For large raster sizes, an undersampling of the original mesh is possible.
+
+The tool is part of the official [OpenGeoSys git repository](https://github.com/ufz/ogs).
+
+## Usage
+
+```bash
+ Vtu2Grid.exe -i <input.vtu> -o <output.vtu> -x <floating point number>
+ [-y <floating point number>] [-z <floating point number>]
+
+Where:
+
+ -i <input.vtu>, --input <input.vtu>
+ (required) the 3D input mesh (*.vtu, *.msh)
+
+ -o <output.vtu>, --output <output.vtu>
+ (required) the output grid (*.vtu)
+
+ -x <floating point number>, --cellsize-x <floating point number>
+ (required) edge length of cubes in x-direction (longitude) or all
+ directions, if y and z are not set
+
+ -y <floating point number>, --cellsize-y <floating point number>
+ edge length of cubes in y-direction (latitude)
+
+ -z <floating point number>, --cellsize-z <floating point number>
+ edge length of cubes in z-direction (depth)
+```
+
+The ```x```/```y```/```z```-parameters determine the raster size. If only ```x``` is given, all raster cells will be (equilateral) cubes. If all three parameters are specified, the programme will create rectangular cuboid cells.
+
+## Simple example
+
+{ width=66% }
+Original, unstructured grid consisting of 217,128 prism-elements. The subsurface represented by this mesh consists of a number of layers of very different thickness. The very thin layers at the top (reddish tones) and in the middle (green and cyan tones) are particularly difficult to handle during numerical simulation.
+
+**Command:**
+```
+Vtu2Grid -i input.vtu -o output.vtu -x 200
+```
+{ width=66% }
+Rasterised grid consisting of 9,240 cubes (equilateral hexahedral elements with an edge length of 200m). The result is severely undersampled and a continuous layer structure is no longer visible.
+
+**Command:**
+```
+Vtu2Grid -i input.vtu -o output.vtu -x 100
+```
+{ width=66% }
+Rasterised grid consisting of 74,048 equilateral hexahedral elements with an edge length of 100m. The result is still undersampled but layers become already visible.
+
+**Command:**
+```
+Vtu2Grid -i input.vtu -o output.vtu -x 50
+```
+{ width=66% }
+Rasterised grid consisting of 591,757 equilateral hexahedral elements with an edge length of 50m. There's still undersampling in regions containing thin layers but the overall structure is reasonably well represented.
+
+**Command:**
+```
+Vtu2Grid -i input.vtu -o output.vtu -x 50 -y 50 -z 10
+```
+{ width=66% }
+Rasterised grid consisting of 2,959,656 cuboid hexahedral elements with an edge length of 50m x 50m x 10m. The structure of the original mesh is very well represented while the number of elements has increased by an order of magnitude.
+
+## Application
+
+This utility can be used to convert complex unstructured 3D meshes with a numerically challenging geometry into simple, rasterised meshes. This will (significantly) increase the number of elements but prevents numerical issues during a subsequent simulation, thus exchanging a challenging model setup with a longer runtime. As such, the output mesh allows to quickly set-up a working model of region of interest to get first results that allow for informed decisions for a more realistic model using the unstructured input mesh.
diff --git a/web/content/docs/tools/meshing/vtu2grid/vtu2grid-100.png b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-100.png
new file mode 100644
index 0000000000000000000000000000000000000000..ebbc99a7c7d08f8344d8f67ad1ecc4f4f7d7b6a6
--- /dev/null
+++ b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-100.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:370adef3353bba3830c4036aacc2c246709cfa9446621f0edd7260030012cdac
+size 164792
diff --git a/web/content/docs/tools/meshing/vtu2grid/vtu2grid-200.png b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-200.png
new file mode 100644
index 0000000000000000000000000000000000000000..8034a65303f438e713aafcb2753593fe79191a77
--- /dev/null
+++ b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-200.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:5061070db60a27e0b0eaa8a4aeec37ff23d86a1503ecafc9821606a2a694821c
+size 106627
diff --git a/web/content/docs/tools/meshing/vtu2grid/vtu2grid-50.png b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-50.png
new file mode 100644
index 0000000000000000000000000000000000000000..3d27df2c6c25f8c159b1c817e57041f5e36c8fa5
--- /dev/null
+++ b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-50.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:8415a9b9d1879e8e59e9dc8f69da0826639ca48350c1770877ce50a20db9baa7
+size 213540
diff --git a/web/content/docs/tools/meshing/vtu2grid/vtu2grid-50x50x10.png b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-50x50x10.png
new file mode 100644
index 0000000000000000000000000000000000000000..1417cfc8b741da9319eda7141be829c9b4ef0e66
--- /dev/null
+++ b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-50x50x10.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:90afa146929b51ff8d10d98de8570ba7d137bbbd7f33b35e3d5008b0ec7bba9d
+size 259655
diff --git a/web/content/docs/tools/meshing/vtu2grid/vtu2grid-orig.png b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-orig.png
new file mode 100644
index 0000000000000000000000000000000000000000..d534624556671c79a3a50547a2e9733410a3fbdc
--- /dev/null
+++ b/web/content/docs/tools/meshing/vtu2grid/vtu2grid-orig.png
@@ -0,0 +1,3 @@
+version https://git-lfs.github.com/spec/v1
+oid sha256:a3db5df6ba9ffd6efe7740ae0ec4b38eb14869e8ec1584294b9e597ad96e9d06
+size 166845