changed: ewoms/io -> opm/models/io

examples/problems/diffusionproblem.hh

@@ -30,7 +30,7 @@
 
 #include <ewoms/models/ncp/ncpproperties.hh>
 
-#include <ewoms/io/cubegridvanguard.hh>
+#include <opm/models/io/cubegridvanguard.hh>
 
 #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>

examples/problems/fingerproblem.hh

@@ -28,7 +28,7 @@
 #ifndef EWOMS_FINGER_PROBLEM_HH
 #define EWOMS_FINGER_PROBLEM_HH
 
-#include <ewoms/io/structuredgridvanguard.hh>
+#include <opm/models/io/structuredgridvanguard.hh>
 
 #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
 #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>

examples/problems/fractureproblem.hh

@@ -38,7 +38,7 @@
 #endif
 
 #include <ewoms/models/discretefracture/discretefracturemodel.hh>
-#include <ewoms/io/dgfvanguard.hh>
+#include <opm/models/io/dgfvanguard.hh>
 
 #include <opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp>
 #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>

examples/problems/lensproblem.hh

@@ -28,7 +28,7 @@
 #ifndef EWOMS_LENS_PROBLEM_HH
 #define EWOMS_LENS_PROBLEM_HH
 
-#include <ewoms/io/structuredgridvanguard.hh>
+#include <opm/models/io/structuredgridvanguard.hh>
 #include <opm/models/immiscible/immiscibleproperties.hh>
 #include <opm/models/discretization/common/fvbaseadlocallinearizer.hh>
 #include <opm/models/discretization/ecfv/ecfvdiscretization.hh>

examples/problems/powerinjectionproblem.hh

@@ -29,7 +29,7 @@
 #define EWOMS_POWER_INJECTION_PROBLEM_HH
 
 #include <opm/models/immiscible/immisciblemodel.hh>
-#include <ewoms/io/cubegridvanguard.hh>
+#include <opm/models/io/cubegridvanguard.hh>
 
 #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
 #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>

examples/tutorial1problem.hh

@@ -45,7 +45,7 @@
 
 // For the DUNE grid
 #include <dune/grid/yaspgrid.hh> /*@\label{tutorial1:include-grid-manager}@*/
-#include <ewoms/io/cubegridvanguard.hh> /*@\label{tutorial1:include-grid-manager}@*/
+#include <opm/models/io/cubegridvanguard.hh> /*@\label{tutorial1:include-grid-manager}@*/
 
 // For Dune::FieldMatrix
 #include <dune/common/fmatrix.hh>

opm/models/blackoil/blackoilenergymodules.hh

@@ -29,7 +29,7 @@
 #define EWOMS_BLACK_OIL_ENERGY_MODULE_HH
 
 #include "blackoilproperties.hh"
-#include <ewoms/io/vtkblackoilenergymodule.hh>
+#include <opm/models/io/vtkblackoilenergymodule.hh>
 #include <opm/models/common/quantitycallbacks.hh>
 
 #include <opm/material/common/Tabulated1DFunction.hpp>

opm/models/blackoil/blackoilfoammodules.hh

@@ -29,7 +29,7 @@
 #define EWOMS_BLACK_OIL_FOAM_MODULE_HH
 
 #include "blackoilproperties.hh"
-//#include <ewoms/io/vtkblackoilfoammodule.hh>
+//#include <opm/models/io/vtkblackoilfoammodule.hh>
 #include <opm/models/common/quantitycallbacks.hh>
 
 #include <opm/material/common/Tabulated1DFunction.hpp>

opm/models/blackoil/blackoilmodel.hh

@@ -47,8 +47,8 @@
 #include "blackoildarcyfluxmodule.hh"
 
 #include <opm/models/common/multiphasebasemodel.hh>
-#include <ewoms/io/vtkcompositionmodule.hh>
-#include <ewoms/io/vtkblackoilmodule.hh>
+#include <opm/models/io/vtkcompositionmodule.hh>
+#include <opm/models/io/vtkblackoilmodule.hh>
 
 #include <opm/material/fluidsystems/BlackOilFluidSystem.hpp>
 #include <opm/material/common/Unused.hpp>

opm/models/blackoil/blackoilpolymermodules.hh

@@ -29,7 +29,7 @@
 #define EWOMS_BLACK_OIL_POLYMER_MODULE_HH
 
 #include "blackoilproperties.hh"
-#include <ewoms/io/vtkblackoilpolymermodule.hh>
+#include <opm/models/io/vtkblackoilpolymermodule.hh>
 #include <opm/models/common/quantitycallbacks.hh>
 
 #include <opm/material/common/Tabulated1DFunction.hpp>

opm/models/blackoil/blackoilsolventmodules.hh

@@ -29,7 +29,7 @@
 #define EWOMS_BLACK_OIL_SOLVENT_MODULE_HH
 
 #include "blackoilproperties.hh"
-#include <ewoms/io/vtkblackoilsolventmodule.hh>
+#include <opm/models/io/vtkblackoilsolventmodule.hh>
 #include <opm/models/common/quantitycallbacks.hh>
 
 #include <opm/material/fluidsystems/blackoilpvt/SolventPvt.hpp>

opm/models/common/multiphasebaseproperties.hh

@@ -31,8 +31,8 @@
 #define EWOMS_MULTI_PHASE_BASE_PROPERTIES_HH
 
 #include <opm/models/discretization/common/fvbaseproperties.hh>
-#include <ewoms/io/vtkmultiphasemodule.hh>
-#include <ewoms/io/vtktemperaturemodule.hh>
+#include <opm/models/io/vtkmultiphasemodule.hh>
+#include <opm/models/io/vtktemperaturemodule.hh>
 
 BEGIN_PROPERTIES
 

opm/models/discretization/common/fvbasenewtonconvergencewriter.hh

@@ -28,7 +28,7 @@
 #ifndef EWOMS_FV_BASE_NEWTON_CONVERGENCE_WRITER_HH
 #define EWOMS_FV_BASE_NEWTON_CONVERGENCE_WRITER_HH
 
-#include <ewoms/io/vtkmultiwriter.hh>
+#include <opm/models/io/vtkmultiwriter.hh>
 #include <opm/models/utils/propertysystem.hh>
 
 #include <iostream>

opm/models/discretization/common/fvbaseproblem.hh

@@ -30,8 +30,8 @@
 
 #include "fvbaseproperties.hh"
 
-#include <ewoms/io/vtkmultiwriter.hh>
-#include <ewoms/io/restart.hh>
+#include <opm/models/io/vtkmultiwriter.hh>
+#include <opm/models/io/restart.hh>
 #include <opm/models/discretization/common/restrictprolong.hh>
 
 #include <opm/material/common/Unused.hpp>

opm/models/discretization/common/fvbaseproperties.hh

@@ -35,7 +35,7 @@
 #include "fvbasefdlocallinearizer.hh"
 
 #include <opm/models/utils/basicproperties.hh>
-#include <ewoms/io/vtkprimaryvarsmodule.hh>
+#include <opm/models/io/vtkprimaryvarsmodule.hh>
 #include <opm/simulators/linalg/parallelbicgstabbackend.hh>
 
 BEGIN_PROPERTIES

opm/models/discretization/ecfv/ecfvbaseoutputmodule.hh

@@ -30,7 +30,7 @@
 
 #include "ecfvproperties.hh"
 
-#include <ewoms/io/vtkmultiwriter.hh>
+#include <opm/models/io/vtkmultiwriter.hh>
 
 namespace Opm {
 /*!

opm/models/discretization/vcfv/vcfvbaseoutputmodule.hh

@@ -30,7 +30,7 @@
 
 #include "vcfvproperties.hh"
 
-#include <ewoms/io/baseoutputwriter.hh>
+#include <opm/models/io/baseoutputwriter.hh>
 
 #include <string>
 #include <vector>

opm/models/flash/flashmodel.hh

@@ -41,9 +41,9 @@
 
 #include <opm/models/common/multiphasebasemodel.hh>
 #include <opm/models/common/energymodule.hh>
-#include <ewoms/io/vtkcompositionmodule.hh>
-#include <ewoms/io/vtkenergymodule.hh>
-#include <ewoms/io/vtkdiffusionmodule.hh>
+#include <opm/models/io/vtkcompositionmodule.hh>
+#include <opm/models/io/vtkenergymodule.hh>
+#include <opm/models/io/vtkdiffusionmodule.hh>
 #include <opm/material/fluidmatrixinteractions/NullMaterial.hpp>
 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
 #include <opm/material/constraintsolvers/NcpFlash.hpp>

opm/models/flash/flashproperties.hh

@@ -31,9 +31,9 @@
 #define EWOMS_FLASH_PROPERTIES_HH
 
 #include <opm/models/common/multiphasebaseproperties.hh>
-#include <ewoms/io/vtkcompositionmodule.hh>
-#include <ewoms/io/vtkenergymodule.hh>
-#include <ewoms/io/vtkdiffusionmodule.hh>
+#include <opm/models/io/vtkcompositionmodule.hh>
+#include <opm/models/io/vtkenergymodule.hh>
+#include <opm/models/io/vtkdiffusionmodule.hh>
 
 BEGIN_PROPERTIES
 

opm/models/immiscible/immisciblemodel.hh

@@ -40,7 +40,7 @@
 
 #include <opm/models/common/multiphasebasemodel.hh>
 #include <opm/models/common/energymodule.hh>
-#include <ewoms/io/vtkenergymodule.hh>
+#include <opm/models/io/vtkenergymodule.hh>
 #include <opm/material/components/NullComponent.hpp>
 #include <opm/material/fluidsystems/GasPhase.hpp>
 #include <opm/material/fluidsystems/LiquidPhase.hpp>

opm/models/immiscible/immiscibleproperties.hh

@@ -31,7 +31,7 @@
 #define EWOMS_IMMISCIBLE_PROPERTIES_HH
 
 #include <opm/models/common/multiphasebaseproperties.hh>
-#include <ewoms/io/vtkenergymodule.hh>
+#include <opm/models/io/vtkenergymodule.hh>
 
 BEGIN_PROPERTIES
 

opm/models/io/baseoutputmodule.hh

@@ -0,0 +1,518 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::BaseOutputModule
+ */
+#ifndef EWOMS_BASE_OUTPUT_MODULE_HH
+#define EWOMS_BASE_OUTPUT_MODULE_HH
+
+#include "baseoutputwriter.hh"
+
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/utils/propertysystem.hh>
+
+#include <opm/material/common/Exceptions.hpp>
+
+#include <dune/istl/bvector.hh>
+#include <dune/common/fvector.hh>
+
+#include <vector>
+#include <sstream>
+#include <string>
+#include <array>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// forward definition of property tags
+NEW_PROP_TAG(NumPhases);
+NEW_PROP_TAG(NumComponents);
+NEW_PROP_TAG(NumEq);
+
+NEW_PROP_TAG(Model);
+NEW_PROP_TAG(Simulator);
+NEW_PROP_TAG(Scalar);
+NEW_PROP_TAG(Evaluation);
+NEW_PROP_TAG(GridView);
+NEW_PROP_TAG(ElementContext);
+NEW_PROP_TAG(FluidSystem);
+NEW_PROP_TAG(DiscBaseOutputModule);
+
+END_PROPERTIES
+
+namespace Opm {
+
+#if __GNUC__ || __clang__
+#pragma GCC diagnostic push
+#pragma GCC diagnostic ignored "-Wpragmas"
+#pragma GCC diagnostic ignored "-Wformat-nonliteral"
+#endif
+
+/*!
+ * \brief The base class for writer modules.
+ *
+ * This class also provides some convenience methods for buffer
+ * management and is the base class for all other output writer
+ * modules.
+ */
+template<class TypeTag>
+class BaseOutputModule
+{
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+    typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
+
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+    enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
+    enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
+    enum { dim = GridView::dimension };
+    enum { dimWorld = GridView::dimensionworld };
+
+    typedef BaseOutputWriter::Tensor Tensor;
+
+public:
+    typedef BaseOutputWriter::ScalarBuffer ScalarBuffer;
+    typedef BaseOutputWriter::VectorBuffer VectorBuffer;
+    typedef BaseOutputWriter::TensorBuffer TensorBuffer;
+
+    typedef std::array<ScalarBuffer, numEq> EqBuffer;
+    typedef std::array<ScalarBuffer, numPhases> PhaseBuffer;
+    typedef std::array<ScalarBuffer, numComponents> ComponentBuffer;
+    typedef std::array<std::array<ScalarBuffer, numComponents>, numPhases> PhaseComponentBuffer;
+
+    typedef std::array<VectorBuffer, numPhases> PhaseVectorBuffer;
+
+    BaseOutputModule(const Simulator& simulator)
+        : simulator_(simulator)
+    {}
+
+    virtual ~BaseOutputModule()
+    {}
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to disk.
+     *
+     * The module can dynamically cast the writer to the desired
+     * concrete class. If the writer is incompatible with the module,
+     * this method should become a no-op.
+     */
+    virtual void allocBuffers() = 0;
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quanties relevant
+     *        for an element
+     *
+     * The module can dynamically cast the writer to the desired
+     * concrete class. If the writer is incompatible with the module,
+     * this method should become a no-op.
+     */
+    virtual void processElement(const ElementContext& elemCtx) = 0;
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    virtual void commitBuffers(BaseOutputWriter& writer) = 0;
+
+    /*!
+     * \brief Returns true iff the module needs to access the extensive quantities of a
+     * context to do its job.
+     *
+     * For example, this happens if velocities or gradients should be written.
+     *
+     * Always returning true here does not do any harm from the correctness perspective,
+     * but it slows down writing the output fields. Since most output modules only write
+     * intensive quantities, this method returns 'false' by default.
+     */
+    virtual bool needExtensiveQuantities() const
+    { return false; }
+
+protected:
+    enum BufferType {
+        //! Buffer contains data associated with the degrees of freedom
+        DofBuffer,
+
+        //! Buffer contains data associated with the grid's vertices
+        VertexBuffer,
+
+        //! Buffer contains data associated with the grid's elements
+        ElementBuffer
+    };
+
+    /*!
+     * \brief Allocate the space for a buffer storing a scalar quantity
+     */
+    void resizeScalarBuffer_(ScalarBuffer& buffer,
+                             BufferType bufferType = DofBuffer)
+    {
+        size_t n;
+        if (bufferType == VertexBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(dim));
+        else if (bufferType == ElementBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(0));
+        else if (bufferType == DofBuffer)
+            n = simulator_.model().numGridDof();
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+
+        buffer.resize(n);
+        std::fill(buffer.begin(), buffer.end(), 0.0);
+    }
+
+    /*!
+     * \brief Allocate the space for a buffer storing a tensorial quantity
+     */
+    void resizeTensorBuffer_(TensorBuffer& buffer,
+                             BufferType bufferType = DofBuffer)
+    {
+        size_t n;
+        if (bufferType == VertexBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(dim));
+        else if (bufferType == ElementBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(0));
+        else if (bufferType == DofBuffer)
+            n = simulator_.model().numGridDof();
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+
+        buffer.resize(n);
+        Tensor nullMatrix(dimWorld, dimWorld, 0.0);
+        std::fill(buffer.begin(), buffer.end(), nullMatrix);
+    }
+
+    /*!
+     * \brief Allocate the space for a buffer storing a equation specific
+     *        quantity
+     */
+    void resizeEqBuffer_(EqBuffer& buffer,
+                         BufferType bufferType = DofBuffer)
+    {
+        size_t n;
+        if (bufferType == VertexBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(dim));
+        else if (bufferType == ElementBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(0));
+        else if (bufferType == DofBuffer)
+            n = simulator_.model().numGridDof();
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+
+        for (unsigned i = 0; i < numEq; ++i) {
+            buffer[i].resize(n);
+            std::fill(buffer[i].begin(), buffer[i].end(), 0.0);
+        }
+    }
+
+    /*!
+     * \brief Allocate the space for a buffer storing a phase-specific
+     *        quantity
+     */
+    void resizePhaseBuffer_(PhaseBuffer& buffer,
+                            BufferType bufferType = DofBuffer)
+    {
+        size_t n;
+        if (bufferType == VertexBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(dim));
+        else if (bufferType == ElementBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(0));
+        else if (bufferType == DofBuffer)
+            n = simulator_.model().numGridDof();
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+
+        for (unsigned i = 0; i < numPhases; ++i) {
+            buffer[i].resize(n);
+            std::fill(buffer[i].begin(), buffer[i].end(), 0.0);
+        }
+    }
+
+    /*!
+     * \brief Allocate the space for a buffer storing a component
+     *        specific quantity
+     */
+    void resizeComponentBuffer_(ComponentBuffer& buffer,
+                                BufferType bufferType = DofBuffer)
+    {
+        size_t n;
+        if (bufferType == VertexBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(dim));
+        else if (bufferType == ElementBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(0));
+        else if (bufferType == DofBuffer)
+            n = simulator_.model().numGridDof();
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+
+        for (unsigned i = 0; i < numComponents; ++i) {
+            buffer[i].resize(n);
+            std::fill(buffer[i].begin(), buffer[i].end(), 0.0);
+        }
+    }
+
+    /*!
+     * \brief Allocate the space for a buffer storing a phase and
+     *        component specific buffer
+     */
+    void resizePhaseComponentBuffer_(PhaseComponentBuffer& buffer,
+                                     BufferType bufferType = DofBuffer)
+    {
+        size_t n;
+        if (bufferType == VertexBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(dim));
+        else if (bufferType == ElementBuffer)
+            n = static_cast<size_t>(simulator_.gridView().size(0));
+        else if (bufferType == DofBuffer)
+            n = simulator_.model().numGridDof();
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+
+        for (unsigned i = 0; i < numPhases; ++i) {
+            for (unsigned j = 0; j < numComponents; ++j) {
+                buffer[i][j].resize(n);
+                std::fill(buffer[i][j].begin(), buffer[i][j].end(), 0.0);
+            }
+        }
+    }
+
+    /*!
+     * \brief Add a buffer containing scalar quantities to the result file.
+     */
+    void commitScalarBuffer_(BaseOutputWriter& baseWriter,
+                             const char *name,
+                             ScalarBuffer& buffer,
+                             BufferType bufferType = DofBuffer)
+    {
+        if (bufferType == DofBuffer)
+            DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer, name);
+        else if (bufferType == VertexBuffer)
+            attachScalarVertexData_(baseWriter, buffer, name);
+        else if (bufferType == ElementBuffer)
+            attachScalarElementData_(baseWriter, buffer, name);
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+    }
+
+    /*!
+     * \brief Add a buffer containing vectorial quantities to the result file.
+     */
+    void commitVectorBuffer_(BaseOutputWriter& baseWriter,
+                             const char *name,
+                             VectorBuffer& buffer,
+                             BufferType bufferType = DofBuffer)
+    {
+        if (bufferType == DofBuffer)
+            DiscBaseOutputModule::attachVectorDofData_(baseWriter, buffer, name);
+        else if (bufferType == VertexBuffer)
+            attachVectorVertexData_(baseWriter, buffer, name);
+        else if (bufferType == ElementBuffer)
+            attachVectorElementData_(baseWriter, buffer, name);
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+    }
+
+    /*!
+     * \brief Add a buffer containing tensorial quantities to the result file.
+     */
+    void commitTensorBuffer_(BaseOutputWriter& baseWriter,
+                             const char *name,
+                             TensorBuffer& buffer,
+                             BufferType bufferType = DofBuffer)
+    {
+        if (bufferType == DofBuffer)
+            DiscBaseOutputModule::attachTensorDofData_(baseWriter, buffer, name);
+        else if (bufferType == VertexBuffer)
+            attachTensorVertexData_(baseWriter, buffer, name);
+        else if (bufferType == ElementBuffer)
+            attachTensorElementData_(baseWriter, buffer, name);
+        else
+            throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+    }
+
+    /*!
+     * \brief Add a buffer with as many variables as PDEs to the result file.
+     */
+    void commitPriVarsBuffer_(BaseOutputWriter& baseWriter,
+                              const char *pattern,
+                              EqBuffer& buffer,
+                              BufferType bufferType = DofBuffer)
+    {
+        char name[512];
+        for (unsigned i = 0; i < numEq; ++i) {
+            std::string eqName = simulator_.model().primaryVarName(i);
+            snprintf(name, 512, pattern, eqName.c_str());
+
+            if (bufferType == DofBuffer)
+                DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer[i], name);
+            else if (bufferType == VertexBuffer)
+                attachScalarVertexData_(baseWriter, buffer[i], name);
+            else if (bufferType == ElementBuffer)
+                attachScalarElementData_(baseWriter, buffer[i], name);
+            else
+                throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+        }
+    }
+
+    /*!
+     * \brief Add a buffer with as many variables as PDEs to the result file.
+     */
+    void commitEqBuffer_(BaseOutputWriter& baseWriter,
+                         const char *pattern,
+                         EqBuffer& buffer,
+                         BufferType bufferType = DofBuffer)
+    {
+        char name[512];
+        for (unsigned i = 0; i < numEq; ++i) {
+            std::ostringstream oss;
+            oss << i;
+            snprintf(name, 512, pattern, oss.str().c_str());
+
+            if (bufferType == DofBuffer)
+                DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer[i], name);
+            else if (bufferType == VertexBuffer)
+                attachScalarVertexData_(baseWriter, buffer[i], name);
+            else if (bufferType == ElementBuffer)
+                attachScalarElementData_(baseWriter, buffer[i], name);
+            else
+                throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+        }
+    }
+
+    /*!
+     * \brief Add a phase-specific buffer to the result file.
+     */
+    void commitPhaseBuffer_(BaseOutputWriter& baseWriter,
+                            const char *pattern,
+                            PhaseBuffer& buffer,
+                            BufferType bufferType = DofBuffer)
+    {
+        char name[512];
+        for (unsigned i = 0; i < numPhases; ++i) {
+            snprintf(name, 512, pattern, FluidSystem::phaseName(i));
+
+            if (bufferType == DofBuffer)
+                DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer[i], name);
+            else if (bufferType == VertexBuffer)
+                attachScalarVertexData_(baseWriter, buffer[i], name);
+            else if (bufferType == ElementBuffer)
+                attachScalarElementData_(baseWriter, buffer[i], name);
+            else
+                throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+        }
+    }
+
+    /*!
+     * \brief Add a component-specific buffer to the result file.
+     */
+    void commitComponentBuffer_(BaseOutputWriter& baseWriter,
+                                const char *pattern,
+                                ComponentBuffer& buffer,
+                                BufferType bufferType = DofBuffer)
+    {
+        char name[512];
+        for (unsigned i = 0; i < numComponents; ++i) {
+            snprintf(name, 512, pattern, FluidSystem::componentName(i));
+
+            if (bufferType == DofBuffer)
+                DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer[i], name);
+            else if (bufferType == VertexBuffer)
+                attachScalarVertexData_(baseWriter, buffer[i], name);
+            else if (bufferType == ElementBuffer)
+                attachScalarElementData_(baseWriter, buffer[i], name);
+            else
+                throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+        }
+    }
+
+    /*!
+     * \brief Add a phase and component specific quantities to the output.
+     */
+    void commitPhaseComponentBuffer_(BaseOutputWriter& baseWriter,
+                                     const char *pattern,
+                                     PhaseComponentBuffer& buffer,
+                                     BufferType bufferType = DofBuffer)
+    {
+        char name[512];
+        for (unsigned i= 0; i < numPhases; ++i) {
+            for (unsigned j = 0; j < numComponents; ++j) {
+                snprintf(name, 512, pattern,
+                         FluidSystem::phaseName(i),
+                         FluidSystem::componentName(j));
+
+                if (bufferType == DofBuffer)
+                    DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer[i][j], name);
+                else if (bufferType == VertexBuffer)
+                    attachScalarVertexData_(baseWriter, buffer[i][j], name);
+                else if (bufferType == ElementBuffer)
+                    attachScalarElementData_(baseWriter, buffer[i][j], name);
+                else
+                    throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
+            }
+        }
+    }
+
+    void attachScalarElementData_(BaseOutputWriter& baseWriter,
+                                  ScalarBuffer& buffer,
+                                  const char *name)
+    { baseWriter.attachScalarElementData(buffer, name); }
+
+    void attachScalarVertexData_(BaseOutputWriter& baseWriter,
+                                 ScalarBuffer& buffer,
+                                 const char *name)
+    { baseWriter.attachScalarVertexData(buffer, name); }
+
+    void attachVectorElementData_(BaseOutputWriter& baseWriter,
+                                  VectorBuffer& buffer,
+                                  const char *name)
+    { baseWriter.attachVectorElementData(buffer, name); }
+
+    void attachVectorVertexData_(BaseOutputWriter& baseWriter,
+                                 VectorBuffer& buffer,
+                                 const char *name)
+    { baseWriter.attachVectorVertexData(buffer, name); }
+
+    void attachTensorElementData_(BaseOutputWriter& baseWriter,
+                                  TensorBuffer& buffer,
+                                  const char *name)
+    { baseWriter.attachTensorElementData(buffer, name); }
+
+    void attachTensorVertexData_(BaseOutputWriter& baseWriter,
+                                 TensorBuffer& buffer,
+                                 const char *name)
+    { baseWriter.attachTensorVertexData(buffer, name); }
+
+    const Simulator& simulator_;
+};
+
+#if __GNUC__ || __clang__
+#pragma GCC diagnostic pop
+#endif
+
+} // namespace Opm
+
+#endif

opm/models/io/baseoutputwriter.hh

@@ -0,0 +1,106 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ *
+ * \copydoc Opm::BaseOutputWriter
+ */
+#ifndef EWOMS_BASE_OUTPUT_WRITER_HH
+#define EWOMS_BASE_OUTPUT_WRITER_HH
+
+#include <dune/common/dynvector.hh>
+#include <dune/common/dynmatrix.hh>
+
+#include <vector>
+
+namespace Opm {
+/*!
+ * \brief The base class for all output writers.
+ *
+ * The sole purpose of this class is to enable RTTI
+ * (i.e. dynamic_cast) on writer objects.
+ */
+class BaseOutputWriter
+{
+public:
+    typedef double Scalar;
+    typedef Dune::DynamicVector<double> Vector;
+    typedef Dune::DynamicMatrix<double> Tensor;
+    typedef std::vector<Scalar> ScalarBuffer;
+    typedef std::vector<Vector> VectorBuffer;
+    typedef std::vector<Tensor> TensorBuffer;
+
+    BaseOutputWriter()
+    {}
+
+    virtual ~BaseOutputWriter()
+    {}
+
+    /*!
+     * \brief Called when ever a new time step or a new grid
+     *        must be written.
+     */
+    virtual void beginWrite(double t) = 0;
+
+    /*!
+     * \brief Add a scalar vertex centered vector field to the output.
+     */
+    virtual void attachScalarVertexData(ScalarBuffer& buf, std::string name) = 0;
+
+    /*!
+     * \brief Add a scalar element centered quantity to the output.
+     */
+    virtual void attachScalarElementData(ScalarBuffer& buf, std::string name) = 0;
+
+    /*!
+     * \brief Add a vectorial vertex centered vector field to the output.
+     */
+    virtual void attachVectorVertexData(VectorBuffer& buf, std::string name) = 0;
+
+    /*!
+     * \brief Add a vectorial element centered quantity to the output.
+     */
+    virtual void attachVectorElementData(VectorBuffer& buf, std::string name) = 0;
+
+    /*!
+     * \brief Add a tensorial vertex centered tensor field to the output.
+     */
+    virtual void attachTensorVertexData(TensorBuffer& buf, std::string name) = 0;
+
+    /*!
+     * \brief Add a tensorial element centered quantity to the output.
+     */
+    virtual void attachTensorElementData(TensorBuffer& buf, std::string name) = 0;
+
+    /*!
+     * \brief Finalizes the current writer.
+     *
+     * This means that everything will be written to disk, except if
+     * the onlyDiscard argument is true. In this case only all managed
+     * buffers are deleted, but no output is written.
+     */
+    virtual void endWrite(bool onlyDiscard = false) = 0;
+};
+} // namespace Opm
+
+#endif

opm/models/io/basevanguard.hh

@@ -0,0 +1,163 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::BaseVanguard
+ */
+#ifndef EWOMS_BASE_VANGUARD_HH
+#define EWOMS_BASE_VANGUARD_HH
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <dune/common/version.hh>
+
+#if HAVE_DUNE_FEM
+#include <dune/fem/space/common/dofmanager.hh>
+#endif
+
+#include <type_traits>
+#include <memory>
+
+BEGIN_PROPERTIES
+
+NEW_PROP_TAG(Grid);
+NEW_PROP_TAG(Vanguard);
+NEW_PROP_TAG(GridView);
+NEW_PROP_TAG(GridPart);
+NEW_PROP_TAG(GridViewLevel);
+NEW_PROP_TAG(GridFile);
+NEW_PROP_TAG(GridGlobalRefinements);
+NEW_PROP_TAG(Simulator);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \brief Provides the base class for most (all?) simulator vanguards.
+ */
+template <class TypeTag>
+class BaseVanguard
+{
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Vanguard) Implementation;
+
+#if HAVE_DUNE_FEM
+    typedef typename GET_PROP_TYPE(TypeTag, GridPart) GridPart;
+#endif
+
+public:
+    BaseVanguard(Simulator& simulator)
+        : simulator_(simulator)
+    {}
+
+    BaseVanguard(const BaseVanguard&) = delete;
+
+    /*!
+     * \brief Returns a reference to the grid view to be used.
+     */
+    const GridView& gridView() const
+    { return *gridView_; }
+
+#if HAVE_DUNE_FEM
+    /*!
+     * \brief Returns a reference to the grid part to be used.
+     */
+    const GridPart& gridPart() const
+    { return *gridPart_; }
+
+    /*!
+     * \brief Returns a reference to the grid part to be used.
+     */
+    GridPart& gridPart()
+    { return *gridPart_; }
+#endif
+
+    /*!
+     * \brief Returns the number of times the grid has been changed since its creation.
+     *
+     * This basically says how often the grid has been adapted in the current simulation
+     * run.
+     */
+    int gridSequenceNumber () const
+    {
+#if HAVE_DUNE_FEM
+        typedef Dune::Fem::DofManager< Grid > FemDofManager;
+        return FemDofManager::instance( asImp_().grid() ).sequence();
+#else
+        return 0; // return the same sequence number >= 0 means the grid never changes
+#endif
+    }
+
+
+    /*!
+     * \brief Distribute the grid (and attached data) over all
+     *        processes.
+     */
+    void loadBalance()
+    {
+        asImp_().grid().loadBalance();
+        updateGridView_();
+    }
+
+protected:
+    // this method should be called after the grid has been allocated
+    void finalizeInit_()
+    {
+        updateGridView_();
+    }
+
+    void updateGridView_()
+    {
+#if HAVE_DUNE_FEM
+        // first delete old grid part
+        // this is due to a bug in dune-fem (dangling reference)
+        gridPart_.reset();
+        gridPart_.reset(new GridPart(asImp_().grid()));
+        gridView_.reset(new GridView(static_cast<GridView>(*gridPart_)));
+        assert(gridView_->size(0) == asImp_().grid().leafGridView().size(0));
+#else
+        gridView_.reset(new GridView(asImp_().grid().leafGridView()));
+#endif
+    }
+
+private:
+    Implementation& asImp_()
+    { return *static_cast<Implementation*>(this); }
+
+    const Implementation& asImp_() const
+    { return *static_cast<const Implementation*>(this); }
+
+    Simulator& simulator_;
+#if HAVE_DUNE_FEM
+    std::unique_ptr<GridPart> gridPart_;
+#endif
+    std::unique_ptr<GridView> gridView_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/cubegridvanguard.hh

@@ -0,0 +1,156 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::CubeGridVanguard
+ */
+#ifndef EWOMS_CUBE_GRID_VANGUARD_HH
+#define EWOMS_CUBE_GRID_VANGUARD_HH
+
+#include <opm/models/io/basevanguard.hh>
+#include <opm/models/utils/basicproperties.hh>
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <dune/grid/utility/structuredgridfactory.hh>
+
+#include <dune/common/fvector.hh>
+
+#include <memory>
+
+BEGIN_PROPERTIES
+
+NEW_PROP_TAG(Scalar);
+NEW_PROP_TAG(Grid);
+
+NEW_PROP_TAG(DomainSizeX);
+NEW_PROP_TAG(DomainSizeY);
+NEW_PROP_TAG(DomainSizeZ);
+
+NEW_PROP_TAG(CellsX);
+NEW_PROP_TAG(CellsY);
+NEW_PROP_TAG(CellsZ);
+
+NEW_PROP_TAG(GridGlobalRefinements);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \brief Provides a simulator vanguad which creates a regular grid made of
+ *        quadrilaterals.
+ *
+ * A quadrilateral is a line segment in 1D, a rectangle in 2D and a
+ * cube in 3D.
+ */
+template <class TypeTag>
+class CubeGridVanguard : public BaseVanguard<TypeTag>
+{
+    typedef BaseVanguard<TypeTag> ParentType;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
+
+    typedef Dune::shared_ptr<Grid> GridPointer;
+    typedef typename Grid::ctype CoordScalar;
+    enum { dimWorld = Grid::dimensionworld };
+    typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
+
+public:
+    /*!
+     * \brief Register all run-time parameters for the simulator vanguad.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, GridGlobalRefinements,
+                             "The number of global refinements of the grid "
+                             "executed after it was loaded");
+        EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeX,
+                             "The size of the domain in x direction");
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsX,
+                             "The number of intervalls in x direction");
+        if (dimWorld > 1) {
+            EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeY,
+                                 "The size of the domain in y direction");
+            EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsY,
+                                 "The number of intervalls in y direction");
+        }
+        if (dimWorld > 2) {
+            EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeZ,
+                                 "The size of the domain in z direction");
+            EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsZ,
+                                 "The number of intervalls in z direction");
+        }
+    }
+
+    /*!
+     * \brief Create the grid
+     */
+    CubeGridVanguard(Simulator& simulator)
+        : ParentType(simulator)
+    {
+        std::array<unsigned int, dimWorld> cellRes;
+        GlobalPosition upperRight(0.0);
+        GlobalPosition lowerLeft(0.0);
+
+        for (unsigned i = 0; i < dimWorld; ++i)
+            cellRes[i] = 0;
+
+        upperRight[0] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeX);
+        cellRes[0] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsX);
+        if (dimWorld > 1) {
+            upperRight[1] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeY);
+            cellRes[1] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsY);
+        }
+        if (dimWorld > 2) {
+            upperRight[2] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeZ);
+            cellRes[2] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsZ);
+        }
+
+        unsigned numRefinements = EWOMS_GET_PARAM(TypeTag, unsigned, GridGlobalRefinements);
+        cubeGrid_ = Dune::StructuredGridFactory<Grid>::createCubeGrid(lowerLeft, upperRight, cellRes);
+        cubeGrid_->globalRefine(static_cast<int>(numRefinements));
+
+        this->finalizeInit_();
+    }
+
+    /*!
+     * \brief Returns a reference to the grid.
+     */
+    Grid& grid()
+    { return *cubeGrid_; }
+
+    /*!
+     * \brief Returns a reference to the grid.
+     */
+    const Grid& grid() const
+    { return *cubeGrid_; }
+
+protected:
+    GridPointer cubeGrid_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/dgfvanguard.hh

@@ -0,0 +1,208 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::DgfVanguard
+ */
+#ifndef EWOMS_DGF_GRID_VANGUARD_HH
+#define EWOMS_DGF_GRID_VANGUARD_HH
+
+#include <dune/grid/io/file/dgfparser/dgfparser.hh>
+#include <dune/grid/common/mcmgmapper.hh>
+#include <opm/models/discretefracture/fracturemapper.hh>
+
+#include <opm/models/io/basevanguard.hh>
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+
+#include <type_traits>
+#include <string>
+
+BEGIN_PROPERTIES
+
+NEW_PROP_TAG(Grid);
+NEW_PROP_TAG(GridFile);
+NEW_PROP_TAG(Vanguard);
+NEW_PROP_TAG(GridGlobalRefinements);
+NEW_PROP_TAG(Scalar);
+NEW_PROP_TAG(Simulator);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \brief Provides a simulator vanguard which creates a grid by parsing a Dune Grid
+ *        Format (DGF) file.
+ */
+template <class TypeTag>
+class DgfVanguard : public BaseVanguard<TypeTag>
+{
+    typedef BaseVanguard<TypeTag> ParentType;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
+    typedef Opm::FractureMapper<TypeTag> FractureMapper;
+
+    typedef std::unique_ptr< Grid > GridPointer;
+
+public:
+    /*!
+     * \brief Register all run-time parameters for the DGF simulator vanguard.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, std::string, GridFile,
+                             "The file name of the DGF file to load");
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, GridGlobalRefinements,
+                             "The number of global refinements of the grid "
+                             "executed after it was loaded");
+    }
+
+    /*!
+     * \brief Load the grid from the file.
+     */
+    DgfVanguard(Simulator& simulator)
+        : ParentType(simulator)
+    {
+        const std::string dgfFileName = EWOMS_GET_PARAM(TypeTag, std::string, GridFile);
+        unsigned numRefinments = EWOMS_GET_PARAM(TypeTag, unsigned, GridGlobalRefinements);
+
+        {
+            // create DGF GridPtr from a dgf file
+            Dune::GridPtr< Grid > dgfPointer( dgfFileName );
+
+            // this is only implemented for 2d currently
+            addFractures_( dgfPointer );
+
+            // store pointer to dune grid
+            gridPtr_.reset( dgfPointer.release() );
+        }
+
+        if (numRefinments > 0)
+            gridPtr_->globalRefine(static_cast<int>(numRefinments));
+
+        this->finalizeInit_();
+    }
+
+    /*!
+     * \brief Returns a reference to the grid.
+     */
+    Grid& grid()
+    { return *gridPtr_; }
+
+    /*!
+     * \brief Returns a reference to the grid.
+     */
+    const Grid& grid() const
+    { return *gridPtr_; }
+
+    /*!
+     * \brief Distributes the grid on all processes of a parallel
+     *        computation.
+     *
+     * This grid manager plays nice and also distributes the data of
+     * the DGF...
+     */
+    void loadBalance()
+    { gridPtr_->loadBalance(); }
+
+    /*!
+     * \brief Returns the fracture mapper
+     *
+     * The fracture mapper determines the topology of the fractures.
+     */
+    FractureMapper& fractureMapper()
+    { return fractureMapper_; }
+
+    /*!
+     * \brief Returns the fracture mapper
+     *
+     * The fracture mapper determines the topology of the fractures.
+     */
+    const FractureMapper& fractureMapper() const
+    { return fractureMapper_; }
+
+protected:
+    void addFractures_(Dune::GridPtr<Grid>& dgfPointer)
+    {
+        typedef typename Grid::LevelGridView LevelGridView;
+
+        // check if fractures are available (only 2d currently)
+        if (dgfPointer.nofParameters(static_cast<int>(Grid::dimension)) == 0)
+            return;
+
+        LevelGridView gridView = dgfPointer->levelGridView(/*level=*/0);
+        const unsigned edgeCodim = Grid::dimension - 1;
+
+#if DUNE_VERSION_NEWER(DUNE_GRID, 2,6)
+        typedef Dune::MultipleCodimMultipleGeomTypeMapper<LevelGridView> VertexMapper;
+        VertexMapper vertexMapper(gridView, Dune::mcmgVertexLayout());
+#else
+        typedef Dune::MultipleCodimMultipleGeomTypeMapper<LevelGridView, Dune::MCMGVertexLayout> VertexMapper;
+        VertexMapper vertexMapper(gridView);
+#endif
+
+        // first create a map of the dune to ART vertex indices
+        auto eIt = gridView.template begin</*codim=*/0>();
+        const auto eEndIt = gridView.template end</*codim=*/0>();
+        for (; eIt != eEndIt; ++eIt) {
+            const auto& element = *eIt;
+            const auto& refElem =
+                Dune::ReferenceElements<Scalar, Grid::dimension>::general(element.type());
+
+            const int edges = refElem.size( edgeCodim );
+            for (int edge = 0; edge < edges; ++edge) {
+                const int vertices = refElem.size(edge, edgeCodim, Grid::dimension);
+                std::vector<unsigned> vertexIndices;
+                vertexIndices.reserve(Grid::dimension);
+                for (int vx = 0; vx < vertices; ++vx) {
+                    // get local vertex number from edge
+                    const int localVx = refElem.subEntity(edge, edgeCodim, vx, Grid::dimension);
+
+                    // get vertex
+                    const auto vertex = element.template subEntity<Grid::dimension>(localVx);
+
+                    // if vertex has parameter 1 insert as a fracture vertex
+                    if (dgfPointer.parameters( vertex )[ 0 ] > 0)
+                        vertexIndices.push_back(
+                            static_cast<unsigned>(vertexMapper.subIndex(element,
+                                                                        static_cast<int>(localVx),
+                                                                        Grid::dimension)));
+                }
+                // if 2 vertices have been found with flag 1 insert a fracture edge
+                if (static_cast<int>(vertexIndices.size()) == Grid::dimension)
+                    fractureMapper_.addFractureEdge(vertexIndices[0], vertexIndices[1]);
+            }
+        }
+    }
+
+private:
+    GridPointer    gridPtr_;
+    FractureMapper fractureMapper_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/restart.hh

@@ -0,0 +1,275 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::Restart
+ */
+#ifndef EWOMS_RESTART_HH
+#define EWOMS_RESTART_HH
+
+#include <string>
+#include <fstream>
+#include <iostream>
+#include <sstream>
+
+namespace Opm {
+
+/*!
+ * \brief Load or save a state of a problem to/from the harddisk.
+ */
+class Restart
+{
+    /*!
+     * \brief Create a magic cookie for restart files, so that it is
+     *        unlikely to load a restart file for an incorrectly.
+     */
+    template <class GridView>
+    static const std::string magicRestartCookie_(const GridView& gridView)
+    {
+        static const std::string gridName = "blubb"; // gridView.grid().name();
+        static const int dim = GridView::dimension;
+
+        int numVertices = gridView.size(dim);
+        int numElements = gridView.size(0);
+        int numEdges = gridView.size(dim - 1);
+        int numCPUs = gridView.comm().size();
+        int rank = gridView.comm().rank();
+
+        std::ostringstream oss;
+        oss << "eWoms restart file: "
+            << "gridName='" << gridName << "' "
+            << "numCPUs=" << numCPUs << " "
+            << "myRank=" << rank << " "
+            << "numElements=" << numElements << " "
+            << "numEdges=" << numEdges << " "
+            << "numVertices=" << numVertices;
+        return oss.str();
+    }
+
+    /*!
+     * \brief Return the restart file name.
+     */
+    template <class GridView, class Scalar>
+    static const std::string restartFileName_(const GridView& gridView,
+                                              const std::string& outputDir,
+                                              const std::string& simName,
+                                              Scalar t)
+    {
+        std::string dir = outputDir;
+        if (dir == ".")
+            dir = "";
+        else if (!dir.empty() && dir.back() != '/')
+            dir += "/";
+
+        int rank = gridView.comm().rank();
+        std::ostringstream oss;
+        oss << dir << simName << "_time=" << t << "_rank=" << rank << ".ers";
+        return oss.str();
+    }
+
+public:
+    /*!
+     * \brief Returns the name of the file which is (de-)serialized.
+     */
+    const std::string& fileName() const
+    { return fileName_; }
+
+    /*!
+     * \brief Write the current state of the model to disk.
+     */
+    template <class Simulator>
+    void serializeBegin(Simulator& simulator)
+    {
+        const std::string magicCookie = magicRestartCookie_(simulator.gridView());
+        fileName_ = restartFileName_(simulator.gridView(),
+                                     simulator.problem().outputDir(),
+                                     simulator.problem().name(),
+                                     simulator.time());
+
+        // open output file and write magic cookie
+        outStream_.open(fileName_.c_str());
+        outStream_.precision(20);
+
+        serializeSectionBegin(magicCookie);
+        serializeSectionEnd();
+    }
+
+    /*!
+     * \brief The output stream to write the serialized data.
+     */
+    std::ostream& serializeStream()
+    { return outStream_; }
+
+    /*!
+     * \brief Start a new section in the serialized output.
+     */
+    void serializeSectionBegin(const std::string& cookie)
+    { outStream_ << cookie << "\n"; }
+
+    /*!
+     * \brief End of a section in the serialized output.
+     */
+    void serializeSectionEnd()
+    { outStream_ << "\n"; }
+
+    /*!
+     * \brief Serialize all leaf entities of a codim in a gridView.
+     *
+     * The actual work is done by Serializer::serialize(Entity)
+     */
+    template <int codim, class Serializer, class GridView>
+    void serializeEntities(Serializer& serializer, const GridView& gridView)
+    {
+        std::ostringstream oss;
+        oss << "Entities: Codim " << codim;
+        std::string cookie = oss.str();
+        serializeSectionBegin(cookie);
+
+        // write element data
+        typedef typename GridView::template Codim<codim>::Iterator Iterator;
+
+        Iterator it = gridView.template begin<codim>();
+        const Iterator& endIt = gridView.template end<codim>();
+        for (; it != endIt; ++it) {
+            serializer.serializeEntity(outStream_, *it);
+            outStream_ << "\n";
+        }
+
+        serializeSectionEnd();
+    }
+
+    /*!
+     * \brief Finish the restart file.
+     */
+    void serializeEnd()
+    { outStream_.close(); }
+
+    /*!
+     * \brief Start reading a restart file at a certain simulated
+     *        time.
+     */
+    template <class Simulator, class Scalar>
+    void deserializeBegin(Simulator& simulator, Scalar t)
+    {
+        fileName_ = restartFileName_(simulator.gridView(), simulator.problem().outputDir(), simulator.problem().name(), t);
+
+        // open input file and read magic cookie
+        inStream_.open(fileName_.c_str());
+        if (!inStream_.good()) {
+            throw std::runtime_error("Restart file '"+fileName_+"' could not be opened properly");
+        }
+
+        // make sure that we don't open an empty file
+        inStream_.seekg(0, std::ios::end);
+        auto pos = inStream_.tellg();
+        if (pos == 0) {
+            throw std::runtime_error("Restart file '"+fileName_+"' is empty");
+        }
+        inStream_.seekg(0, std::ios::beg);
+
+        const std::string magicCookie = magicRestartCookie_(simulator.gridView());
+
+        deserializeSectionBegin(magicCookie);
+        deserializeSectionEnd();
+    }
+
+    /*!
+     * \brief The input stream to read the data which ought to be
+     *        deserialized.
+     */
+    std::istream& deserializeStream()
+    { return inStream_; }
+
+    /*!
+     * \brief Start reading a new section of the restart file.
+     */
+    void deserializeSectionBegin(const std::string& cookie)
+    {
+        if (!inStream_.good())
+            throw std::runtime_error("Encountered unexpected EOF in restart file.");
+        std::string buf;
+        std::getline(inStream_, buf);
+        if (buf != cookie)
+            throw std::runtime_error("Could not start section '"+cookie+"'");
+    }
+
+    /*!
+     * \brief End of a section in the serialized output.
+     */
+    void deserializeSectionEnd()
+    {
+        std::string dummy;
+        std::getline(inStream_, dummy);
+        for (unsigned i = 0; i < dummy.length(); ++i) {
+            if (!std::isspace(dummy[i])) {
+                throw std::logic_error("Encountered unread values while deserializing");
+            }
+        }
+    }
+
+    /*!
+     * \brief Deserialize all leaf entities of a codim in a grid.
+     *
+     * The actual work is done by Deserializer::deserialize(Entity)
+     */
+    template <int codim, class Deserializer, class GridView>
+    void deserializeEntities(Deserializer& deserializer, const GridView& gridView)
+    {
+        std::ostringstream oss;
+        oss << "Entities: Codim " << codim;
+        std::string cookie = oss.str();
+        deserializeSectionBegin(cookie);
+
+        std::string curLine;
+
+        // read entity data
+        typedef typename GridView::template Codim<codim>::Iterator Iterator;
+        Iterator it = gridView.template begin<codim>();
+        const Iterator& endIt = gridView.template end<codim>();
+        for (; it != endIt; ++it) {
+            if (!inStream_.good()) {
+                throw std::runtime_error("Restart file is corrupted");
+            }
+
+            std::getline(inStream_, curLine);
+            std::istringstream curLineStream(curLine);
+            deserializer.deserializeEntity(curLineStream, *it);
+        }
+
+        deserializeSectionEnd();
+    }
+
+    /*!
+     * \brief Stop reading the restart file.
+     */
+    void deserializeEnd()
+    { inStream_.close(); }
+
+private:
+    std::string fileName_;
+    std::ifstream inStream_;
+    std::ofstream outStream_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/simplexvanguard.hh

@@ -0,0 +1,151 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::SimplexGridVanguard
+ */
+#ifndef EWOMS_SIMPLEX_GRID_VANGUARD_HH
+#define EWOMS_SIMPLEX_GRID_VANGUARD_HH
+
+#include <opm/models/utils/basicproperties.hh>
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <dune/grid/utility/structuredgridfactory.hh>
+#include <dune/common/fvector.hh>
+
+#include <memory>
+
+BEGIN_PROPERTIES
+
+NEW_PROP_TAG(Scalar);
+NEW_PROP_TAG(Grid);
+
+NEW_PROP_TAG(DomainSizeX);
+NEW_PROP_TAG(DomainSizeY);
+NEW_PROP_TAG(DomainSizeZ);
+
+NEW_PROP_TAG(CellsX);
+NEW_PROP_TAG(CellsY);
+NEW_PROP_TAG(CellsZ);
+
+NEW_PROP_TAG(GridGlobalRefinements);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \brief Provides a simulator vanguard which a creates regular grid made of simplices.
+ */
+template <class TypeTag>
+class SimplexGridVanguard
+{
+    typedef BaseVanguard<TypeTag> ParentType;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
+
+    typedef Dune::shared_ptr<Grid> GridPointer;
+    typedef typename Grid::ctype CoordScalar;
+    enum { dimWorld = Grid::dimensionworld };
+    typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
+
+public:
+    /*!
+     * \brief Register all run-time parameters for the grid manager.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, GridGlobalRefinements,
+                             "The number of global refinements of the grid "
+                             "executed after it was loaded");
+        EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeX,
+                             "The size of the domain in x direction");
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsX,
+                             "The number of intervalls in x direction");
+        if (dimWorld > 1) {
+            EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeY,
+                                 "The size of the domain in y direction");
+            EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsY,
+                                 "The number of intervalls in y direction");
+        }
+        if (dimWorld > 2) {
+            EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeZ,
+                                 "The size of the domain in z direction");
+            EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsZ,
+                                 "The number of intervalls in z direction");
+        }
+    }
+
+    /*!
+     * \brief Create the Grid
+     */
+    SimplexGridVanguard(Simulator& simulator)
+        : ParentType(simulator)
+    {
+        Dune::array<unsigned, dimWorld> cellRes;
+        GlobalPosition upperRight;
+        GlobalPosition lowerLeft;
+
+        lowerLeft[0] = 0.0;
+        upperRight[0] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeX);
+        cellRes[0] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsX);
+        if (dimWorld > 1) {
+            lowerLeft[1] = 0.0;
+            upperRight[1] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeY);
+            cellRes[1] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsY);
+        }
+        if (dimWorld > 2) {
+            lowerLeft[2] = 0.0;
+            upperRight[2] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeZ);
+            cellRes[2] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsZ);
+        }
+
+        simplexGrid_ = Dune::StructuredGridFactory<Grid>::createSimplexGrid(lowerLeft,
+                                                                            upperRight,
+                                                                            cellRes);
+
+        unsigned numRefinments = EWOMS_GET_PARAM(TypeTag, unsigned, GridGlobalRefinements);
+        simplexGrid_->globalRefine(numRefinments);
+
+        this->finalizeInit_();
+    }
+
+    /*!
+     * \brief Returns a reference to the grid.
+     */
+    Grid& grid()
+    { return simplexGrid_; }
+
+    /*!
+     * \brief Returns a reference to the grid.
+     */
+    const Grid& grid() const
+    { return *simplexGrid_; }
+
+private:
+    GridPointer simplexGrid_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/structuredgridvanguard.hh

@@ -0,0 +1,199 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::StructuredGridVanguard
+ */
+#ifndef EWOMS_STRUCTURED_GRID_VANGUARD_HH
+#define EWOMS_STRUCTURED_GRID_VANGUARD_HH
+
+#include <opm/models/io/basevanguard.hh>
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <dune/grid/yaspgrid.hh>
+#include <dune/grid/io/file/dgfparser/dgfyasp.hh>
+
+#if HAVE_DUNE_ALUGRID
+#include <dune/alugrid/grid.hh>
+#include <dune/alugrid/dgf.hh>
+#endif
+
+#include <dune/common/fvector.hh>
+#include <dune/common/version.hh>
+
+#include <vector>
+#include <memory>
+
+namespace Opm {
+
+template <class TypeTag>
+class StructuredGridVanguard;
+
+} // namespace Opm
+
+BEGIN_PROPERTIES
+
+NEW_TYPE_TAG(StructuredGridVanguard);
+
+// declare the properties required by the for the structured grid simulator vanguard
+NEW_PROP_TAG(Grid);
+NEW_PROP_TAG(Scalar);
+
+NEW_PROP_TAG(DomainSizeX);
+NEW_PROP_TAG(DomainSizeY);
+NEW_PROP_TAG(DomainSizeZ);
+
+NEW_PROP_TAG(CellsX);
+NEW_PROP_TAG(CellsY);
+NEW_PROP_TAG(CellsZ);
+
+NEW_PROP_TAG(GridGlobalRefinements);
+
+// GRIDDIM is only set by the finger problem
+#ifndef GRIDDIM
+static const int dim = 2;
+#else
+static const int dim = GRIDDIM;
+#endif
+
+// set the Grid and Vanguard properties
+#if HAVE_DUNE_ALUGRID
+SET_TYPE_PROP(StructuredGridVanguard, Grid, Dune::ALUGrid< dim, dim, Dune::cube, Dune::nonconforming >);
+#else
+SET_TYPE_PROP(StructuredGridVanguard, Grid, Dune::YaspGrid< dim >);
+#endif
+
+SET_TYPE_PROP(StructuredGridVanguard, Vanguard, Opm::StructuredGridVanguard<TypeTag>);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \ingroup TestProblems
+ *
+ * \brief Helper class for grid instantiation of the lens problem.
+ */
+template <class TypeTag>
+class StructuredGridVanguard : public BaseVanguard<TypeTag>
+{
+    typedef BaseVanguard<TypeTag> ParentType;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
+
+    typedef std::unique_ptr<Grid> GridPointer;
+
+    static const int dim = Grid::dimension;
+
+public:
+    /*!
+     * \brief Register all run-time parameters for the structured grid simulator vanguard.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, GridGlobalRefinements,
+                             "The number of global refinements of the grid "
+                             "executed after it was loaded");
+        EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeX,
+                             "The size of the domain in x direction");
+        EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsX,
+                             "The number of intervalls in x direction");
+        if (dim > 1) {
+            EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeY,
+                                 "The size of the domain in y direction");
+            EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsY,
+                                 "The number of intervalls in y direction");
+        }
+        if (dim > 2) {
+            EWOMS_REGISTER_PARAM(TypeTag, Scalar, DomainSizeZ,
+                                 "The size of the domain in z direction");
+            EWOMS_REGISTER_PARAM(TypeTag, unsigned, CellsZ,
+                                 "The number of intervalls in z direction");
+        }
+    }
+
+    /*!
+     * \brief Create the grid for the lens problem
+     */
+    StructuredGridVanguard(Simulator& simulator)
+        : ParentType(simulator)
+    {
+        Dune::FieldVector<int, dim> cellRes;
+
+        typedef double GridScalar;
+        Dune::FieldVector<GridScalar, dim> upperRight;
+        Dune::FieldVector<GridScalar, dim> lowerLeft( 0 );
+
+        upperRight[0] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeX);
+        upperRight[1] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeY);
+
+        cellRes[0] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsX);
+        cellRes[1] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsY);
+        if (dim == 3) {
+            upperRight[2] = EWOMS_GET_PARAM(TypeTag, Scalar, DomainSizeZ);
+            cellRes[2] = EWOMS_GET_PARAM(TypeTag, unsigned, CellsZ);
+        }
+
+        std::stringstream dgffile;
+        dgffile << "DGF" << std::endl;
+        dgffile << "INTERVAL" << std::endl;
+        dgffile << lowerLeft  << std::endl;
+        dgffile << upperRight << std::endl;
+        dgffile << cellRes    << std::endl;
+        dgffile << "#" << std::endl;
+        dgffile << "GridParameter" << std::endl;
+        dgffile << "overlap 1" << std::endl;
+        dgffile << "#" << std::endl;
+        dgffile << "Simplex" << std::endl;
+        dgffile << "#" << std::endl;
+
+        // use DGF parser to create a grid from interval block
+        gridPtr_.reset( Dune::GridPtr< Grid >( dgffile ).release() );
+
+        unsigned numRefinements = EWOMS_GET_PARAM(TypeTag, unsigned, GridGlobalRefinements);
+        gridPtr_->globalRefine(static_cast<int>(numRefinements));
+
+        this->finalizeInit_();
+    }
+
+    /*!
+     * \brief Return a reference to the grid object.
+     */
+    Grid& grid()
+    { return *gridPtr_; }
+
+    /*!
+     * \brief Return a constant reference to the grid object.
+     */
+    const Grid& grid() const
+    { return *gridPtr_; }
+
+private:
+    GridPointer gridPtr_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkblackoilenergymodule.hh

@@ -0,0 +1,233 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkBlackOilEnergyModule
+ */
+#ifndef EWOMS_VTK_BLACK_OIL_ENERGY_MODULE_HH
+#define EWOMS_VTK_BLACK_OIL_ENERGY_MODULE_HH
+
+#include <opm/material/densead/Math.hpp>
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/blackoil/blackoilproperties.hh>
+
+#include <dune/common/fvector.hh>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK multi-phase output
+NEW_TYPE_TAG(VtkBlackOilEnergy);
+
+// create the property tags needed for the energy module
+NEW_PROP_TAG(VtkWriteRockInternalEnergy);
+NEW_PROP_TAG(VtkWriteTotalThermalConductivity);
+NEW_PROP_TAG(VtkWriteFluidInternalEnergies);
+NEW_PROP_TAG(VtkWriteFluidEnthalpies);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkBlackOilEnergy, VtkWriteRockInternalEnergy, true);
+SET_BOOL_PROP(VtkBlackOilEnergy, VtkWriteTotalThermalConductivity, true);
+SET_BOOL_PROP(VtkBlackOilEnergy, VtkWriteFluidInternalEnergies, true);
+SET_BOOL_PROP(VtkBlackOilEnergy, VtkWriteFluidEnthalpies, true);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the black oil model's energy related quantities.
+ */
+template <class TypeTag>
+class VtkBlackOilEnergyModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy) };
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+    typedef typename ParentType::PhaseBuffer PhaseBuffer;
+
+public:
+    VtkBlackOilEnergyModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the multi-phase VTK output
+     * module.
+     */
+    static void registerParameters()
+    {
+        if (!enableEnergy)
+            return;
+
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteRockInternalEnergy,
+                             "Include the volumetric internal energy of rock "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteTotalThermalConductivity,
+                             "Include the total thermal conductivity of the medium and the fluids "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFluidInternalEnergies,
+                             "Include the internal energies of the fluids "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFluidEnthalpies,
+                             "Include the enthalpies of the fluids "
+                             "in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        if (!enableEnergy)
+            return;
+
+        if (rockInternalEnergyOutput_())
+            this->resizeScalarBuffer_(rockInternalEnergy_);
+        if (totalThermalConductivityOutput_())
+            this->resizeScalarBuffer_(totalThermalConductivity_);
+        if (fluidInternalEnergiesOutput_())
+            this->resizePhaseBuffer_(fluidInternalEnergies_);
+        if (fluidEnthalpiesOutput_())
+            this->resizePhaseBuffer_(fluidEnthalpies_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant for
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        if (!enableEnergy)
+            return;
+
+        for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
+            const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
+            unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
+
+            if (rockInternalEnergyOutput_())
+                rockInternalEnergy_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.rockInternalEnergy());
+
+            if (totalThermalConductivityOutput_())
+                totalThermalConductivity_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.totalThermalConductivity());
+
+            for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
+                if (fluidInternalEnergiesOutput_())
+                    fluidInternalEnergies_[phaseIdx][globalDofIdx] =
+                        Opm::scalarValue(intQuants.fluidState().internalEnergy(phaseIdx));
+
+                if (fluidEnthalpiesOutput_())
+                    fluidEnthalpies_[phaseIdx][globalDofIdx] =
+                        Opm::scalarValue(intQuants.fluidState().enthalpy(phaseIdx));
+            }
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter)
+            return;
+
+        if (!enableEnergy)
+            return;
+
+        if (rockInternalEnergyOutput_())
+            this->commitScalarBuffer_(baseWriter, "volumetric internal energy rock", rockInternalEnergy_);
+
+        if (totalThermalConductivityOutput_())
+            this->commitScalarBuffer_(baseWriter, "total thermal conductivity", totalThermalConductivity_);
+
+        if (fluidInternalEnergiesOutput_())
+            this->commitPhaseBuffer_(baseWriter, "internal energy_%s", fluidInternalEnergies_);
+
+        if (fluidEnthalpiesOutput_())
+            this->commitPhaseBuffer_(baseWriter, "enthalpy_%s", fluidEnthalpies_);
+    }
+
+private:
+    static bool rockInternalEnergyOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteRockInternalEnergy);
+        return val;
+    }
+
+    static bool totalThermalConductivityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteTotalThermalConductivity);
+        return val;
+    }
+
+    static bool fluidInternalEnergiesOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFluidInternalEnergies);
+        return val;
+    }
+
+    static bool fluidEnthalpiesOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFluidEnthalpies);
+        return val;
+    }
+
+    ScalarBuffer rockInternalEnergy_;
+    ScalarBuffer totalThermalConductivity_;
+    PhaseBuffer fluidInternalEnergies_;
+    PhaseBuffer fluidEnthalpies_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/vtkblackoilmodule.hh

@@ -0,0 +1,390 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkBlackOilModule
+ */
+#ifndef EWOMS_VTK_BLACK_OIL_MODULE_HH
+#define EWOMS_VTK_BLACK_OIL_MODULE_HH
+
+#include <opm/material/densead/Math.hpp>
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/blackoil/blackoilproperties.hh>
+
+#include <dune/common/fvector.hh>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK multi-phase output
+NEW_TYPE_TAG(VtkBlackOil);
+
+// create the property tags needed for the multi phase module
+NEW_PROP_TAG(EnableVtkOutput);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(VtkWriteGasDissolutionFactor);
+NEW_PROP_TAG(VtkWriteOilVaporizationFactor);
+NEW_PROP_TAG(VtkWriteOilFormationVolumeFactor);
+NEW_PROP_TAG(VtkWriteGasFormationVolumeFactor);
+NEW_PROP_TAG(VtkWriteWaterFormationVolumeFactor);
+NEW_PROP_TAG(VtkWriteOilSaturationPressure);
+NEW_PROP_TAG(VtkWriteGasSaturationPressure);
+NEW_PROP_TAG(VtkWriteSaturationRatios);
+NEW_PROP_TAG(VtkWriteSaturatedOilGasDissolutionFactor);
+NEW_PROP_TAG(VtkWriteSaturatedGasOilVaporizationFactor);
+NEW_PROP_TAG(VtkWritePrimaryVarsMeaning);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkBlackOil, VtkWriteGasDissolutionFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteOilVaporizationFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteOilFormationVolumeFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteGasFormationVolumeFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteWaterFormationVolumeFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteOilSaturationPressure, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteGasSaturationPressure, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteSaturationRatios, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteSaturatedOilGasDissolutionFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWriteSaturatedGasOilVaporizationFactor, false);
+SET_BOOL_PROP(VtkBlackOil, VtkWritePrimaryVarsMeaning, false);
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the black oil model's parameters.
+ */
+template <class TypeTag>
+class VtkBlackOilModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
+    enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
+    enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
+
+    enum { gasCompIdx = FluidSystem::gasCompIdx };
+    enum { oilCompIdx = FluidSystem::oilCompIdx };
+    enum { waterCompIdx = FluidSystem::waterCompIdx };
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+
+public:
+    VtkBlackOilModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the multi-phase VTK output
+     * module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteGasDissolutionFactor,
+                             "Include the gas dissolution factor (R_s) of the observed oil "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteOilVaporizationFactor,
+                             "Include the oil vaporization factor (R_v) of the observed gas "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteOilFormationVolumeFactor,
+                             "Include the oil formation volume factor (B_o) in the "
+                             "VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteGasFormationVolumeFactor,
+                             "Include the gas formation volume factor (B_g) in the "
+                             "VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteWaterFormationVolumeFactor,
+                             "Include the water formation volume factor (B_w) in the "
+                             "VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteOilSaturationPressure,
+                             "Include the saturation pressure of oil (p_o,sat) in the "
+                             "VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteGasSaturationPressure,
+                             "Include the saturation pressure of gas (p_g,sat) in the "
+                             "VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSaturatedOilGasDissolutionFactor,
+                             "Include the gas dissolution factor (R_s,sat) of gas saturated "
+                             "oil in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSaturatedGasOilVaporizationFactor,
+                             "Include the oil vaporization factor (R_v,sat) of oil saturated "
+                             "gas in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSaturationRatios,
+                             "Write the ratio of the actually and maximum dissolved component of "
+                             "the mixtures");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePrimaryVarsMeaning,
+                             "Include how the primary variables should be interpreted to the "
+                             "VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (gasDissolutionFactorOutput_())
+            this->resizeScalarBuffer_(gasDissolutionFactor_);
+        if (oilVaporizationFactorOutput_())
+            this->resizeScalarBuffer_(oilVaporizationFactor_);
+        if (oilFormationVolumeFactorOutput_())
+            this->resizeScalarBuffer_(oilFormationVolumeFactor_);
+        if (gasFormationVolumeFactorOutput_())
+            this->resizeScalarBuffer_(gasFormationVolumeFactor_);
+        if (waterFormationVolumeFactorOutput_())
+            this->resizeScalarBuffer_(waterFormationVolumeFactor_);
+        if (oilSaturationPressureOutput_())
+            this->resizeScalarBuffer_(oilSaturationPressure_);
+        if (gasSaturationPressureOutput_())
+            this->resizeScalarBuffer_(gasSaturationPressure_);
+        if (saturatedOilGasDissolutionFactorOutput_())
+            this->resizeScalarBuffer_(saturatedOilGasDissolutionFactor_);
+        if (saturatedGasOilVaporizationFactorOutput_())
+            this->resizeScalarBuffer_(saturatedGasOilVaporizationFactor_);
+        if (saturationRatiosOutput_()) {
+            this->resizeScalarBuffer_(oilSaturationRatio_);
+            this->resizeScalarBuffer_(gasSaturationRatio_);
+        }
+        if (primaryVarsMeaningOutput_())
+            this->resizeScalarBuffer_(primaryVarsMeaning_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant for
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
+            const auto& fs = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0).fluidState();
+            typedef typename std::remove_const<typename std::remove_reference<decltype(fs)>::type>::type FluidState;
+            unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
+
+            const auto& primaryVars = elemCtx.primaryVars(dofIdx, /*timeIdx=*/0);
+
+            unsigned pvtRegionIdx = elemCtx.primaryVars(dofIdx, /*timeIdx=*/0).pvtRegionIndex();
+            Scalar SoMax = std::max(Opm::getValue(fs.saturation(oilPhaseIdx)),
+                                    elemCtx.problem().maxOilSaturation(globalDofIdx));
+            Scalar x_oG = Opm::getValue(fs.moleFraction(oilPhaseIdx, gasCompIdx));
+            Scalar x_gO = Opm::getValue(fs.moleFraction(gasPhaseIdx, oilCompIdx));
+            Scalar X_oG = Opm::getValue(fs.massFraction(oilPhaseIdx, gasCompIdx));
+            Scalar X_gO = Opm::getValue(fs.massFraction(gasPhaseIdx, oilCompIdx));
+            Scalar Rs = FluidSystem::convertXoGToRs(X_oG, pvtRegionIdx);
+            Scalar Rv = FluidSystem::convertXgOToRv(X_gO, pvtRegionIdx);
+
+            Scalar RsSat =
+                FluidSystem::template saturatedDissolutionFactor<FluidState, Scalar>(fs,
+                                                                                     oilPhaseIdx,
+                                                                                     pvtRegionIdx,
+                                                                                     SoMax);
+            Scalar X_oG_sat = FluidSystem::convertRsToXoG(RsSat, pvtRegionIdx);
+            Scalar x_oG_sat = FluidSystem::convertXoGToxoG(X_oG_sat, pvtRegionIdx);
+
+            Scalar RvSat =
+                FluidSystem::template saturatedDissolutionFactor<FluidState, Scalar>(fs,
+                                                                                     gasPhaseIdx,
+                                                                                     pvtRegionIdx,
+                                                                                     SoMax);
+            Scalar X_gO_sat = FluidSystem::convertRvToXgO(RvSat, pvtRegionIdx);
+            Scalar x_gO_sat = FluidSystem::convertXgOToxgO(X_gO_sat, pvtRegionIdx);
+
+            if (gasDissolutionFactorOutput_())
+                gasDissolutionFactor_[globalDofIdx] = Rs;
+            if (oilVaporizationFactorOutput_())
+                oilVaporizationFactor_[globalDofIdx] = Rv;
+            if (oilFormationVolumeFactorOutput_())
+                oilFormationVolumeFactor_[globalDofIdx] =
+                    1.0/FluidSystem::template inverseFormationVolumeFactor<FluidState, Scalar>(fs, oilPhaseIdx, pvtRegionIdx);
+            if (gasFormationVolumeFactorOutput_())
+                gasFormationVolumeFactor_[globalDofIdx] =
+                    1.0/FluidSystem::template inverseFormationVolumeFactor<FluidState, Scalar>(fs, gasPhaseIdx, pvtRegionIdx);
+            if (waterFormationVolumeFactorOutput_())
+                waterFormationVolumeFactor_[globalDofIdx] =
+                    1.0/FluidSystem::template inverseFormationVolumeFactor<FluidState, Scalar>(fs, waterPhaseIdx, pvtRegionIdx);
+            if (oilSaturationPressureOutput_())
+                oilSaturationPressure_[globalDofIdx] =
+                    FluidSystem::template saturationPressure<FluidState, Scalar>(fs, oilPhaseIdx, pvtRegionIdx);
+            if (gasSaturationPressureOutput_())
+                gasSaturationPressure_[globalDofIdx] =
+                    FluidSystem::template saturationPressure<FluidState, Scalar>(fs, gasPhaseIdx, pvtRegionIdx);
+            if (saturatedOilGasDissolutionFactorOutput_())
+                saturatedOilGasDissolutionFactor_[globalDofIdx] = RsSat;
+            if (saturatedGasOilVaporizationFactorOutput_())
+                saturatedGasOilVaporizationFactor_[globalDofIdx] = RvSat;
+            if (saturationRatiosOutput_()) {
+                if (x_oG_sat <= 0.0)
+                    oilSaturationRatio_[globalDofIdx] = 1.0;
+                else
+                    oilSaturationRatio_[globalDofIdx] = x_oG / x_oG_sat;
+
+                if (x_gO_sat <= 0.0)
+                    gasSaturationRatio_[globalDofIdx] = 1.0;
+                else
+                    gasSaturationRatio_[globalDofIdx] = x_gO / x_gO_sat;
+            }
+
+            if (primaryVarsMeaningOutput_())
+                primaryVarsMeaning_[globalDofIdx] =
+                    primaryVars.primaryVarsMeaning();
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter)
+            return;
+
+        if (gasDissolutionFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "R_s", gasDissolutionFactor_);
+        if (oilVaporizationFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "R_v", oilVaporizationFactor_);
+        if (oilFormationVolumeFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "B_o", oilFormationVolumeFactor_);
+        if (gasFormationVolumeFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "B_g", gasFormationVolumeFactor_);
+        if (waterFormationVolumeFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "B_w", waterFormationVolumeFactor_);
+        if (oilSaturationPressureOutput_())
+            this->commitScalarBuffer_(baseWriter, "p_o,sat", oilSaturationPressure_);
+        if (gasSaturationPressureOutput_())
+            this->commitScalarBuffer_(baseWriter, "p_g,sat", gasSaturationPressure_);
+        if (saturatedOilGasDissolutionFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "R_s,sat", saturatedOilGasDissolutionFactor_);
+        if (saturatedGasOilVaporizationFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "R_v,sat", saturatedGasOilVaporizationFactor_);
+        if (saturationRatiosOutput_()) {
+            this->commitScalarBuffer_(baseWriter, "saturation ratio_oil", oilSaturationRatio_);
+            this->commitScalarBuffer_(baseWriter, "saturation ratio_gas", gasSaturationRatio_);
+        }
+
+        if (primaryVarsMeaningOutput_())
+            this->commitScalarBuffer_(baseWriter, "primary vars meaning", primaryVarsMeaning_);
+    }
+
+private:
+    static bool gasDissolutionFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteGasDissolutionFactor);
+        return val;
+    }
+
+    static bool oilVaporizationFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteOilVaporizationFactor);
+        return val;
+    }
+
+    static bool oilFormationVolumeFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteOilFormationVolumeFactor);
+        return val;
+    }
+
+    static bool gasFormationVolumeFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteGasFormationVolumeFactor);
+        return val;
+    }
+
+    static bool waterFormationVolumeFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteWaterFormationVolumeFactor);
+        return val;
+    }
+
+    static bool oilSaturationPressureOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteOilSaturationPressure);
+        return val;
+    }
+
+    static bool gasSaturationPressureOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteGasSaturationPressure);
+        return val;
+    }
+
+    static bool saturatedOilGasDissolutionFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSaturatedOilGasDissolutionFactor);
+        return val;
+    }
+
+    static bool saturatedGasOilVaporizationFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSaturatedGasOilVaporizationFactor);
+        return val;
+    }
+
+    static bool saturationRatiosOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSaturationRatios);
+        return val;
+    }
+
+    static bool primaryVarsMeaningOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePrimaryVarsMeaning);
+        return val;
+    }
+
+    ScalarBuffer gasDissolutionFactor_;
+    ScalarBuffer oilVaporizationFactor_;
+    ScalarBuffer oilFormationVolumeFactor_;
+    ScalarBuffer gasFormationVolumeFactor_;
+    ScalarBuffer waterFormationVolumeFactor_;
+    ScalarBuffer oilSaturationPressure_;
+    ScalarBuffer gasSaturationPressure_;
+
+    ScalarBuffer saturatedOilGasDissolutionFactor_;
+    ScalarBuffer saturatedGasOilVaporizationFactor_;
+    ScalarBuffer oilSaturationRatio_;
+    ScalarBuffer gasSaturationRatio_;
+
+    ScalarBuffer primaryVarsMeaning_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/vtkblackoilpolymermodule.hh

@@ -0,0 +1,271 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkBlackOilPolymerModule
+ */
+#ifndef EWOMS_VTK_BLACK_OIL_POLYMER_MODULE_HH
+#define EWOMS_VTK_BLACK_OIL_POLYMER_MODULE_HH
+
+#include <opm/material/densead/Math.hpp>
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/blackoil/blackoilproperties.hh>
+
+#include <dune/common/fvector.hh>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK multi-phase output
+NEW_TYPE_TAG(VtkBlackOilPolymer);
+
+// create the property tags needed for the polymer output module
+NEW_PROP_TAG(EnablePolymer);
+NEW_PROP_TAG(EnableVtkOutput);
+NEW_PROP_TAG(VtkWritePolymerConcentration);
+NEW_PROP_TAG(VtkWritePolymerDeadPoreVolume);
+NEW_PROP_TAG(VtkWritePolymerAdsorption);
+NEW_PROP_TAG(VtkWritePolymerRockDensity);
+NEW_PROP_TAG(VtkWritePolymerViscosityCorrection);
+NEW_PROP_TAG(VtkWriteWaterViscosityCorrection);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerConcentration, true);
+SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerDeadPoreVolume, true);
+SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerViscosityCorrection, true);
+SET_BOOL_PROP(VtkBlackOilPolymer, VtkWriteWaterViscosityCorrection, true);
+SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerRockDensity, true);
+SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerAdsorption, true);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the black oil model's polymer related quantities.
+ */
+template <class TypeTag>
+class VtkBlackOilPolymerModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { enablePolymer = GET_PROP_VALUE(TypeTag, EnablePolymer) };
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+
+public:
+    VtkBlackOilPolymerModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the multi-phase VTK output
+     * module.
+     */
+    static void registerParameters()
+    {
+        if (!enablePolymer)
+            return;
+
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerConcentration,
+                             "Include the concentration of the polymer component in the water phase "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerDeadPoreVolume,
+                             "Include the fraction of the \"dead\" pore volume "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerRockDensity,
+                             "Include the amount of already adsorbed polymer component"
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerAdsorption,
+                             "Include the adsorption rate of the polymer component"
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerViscosityCorrection,
+                             "Include the viscosity correction of the polymer component "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteWaterViscosityCorrection,
+                             "Include the viscosity correction of the water component "
+                             "due to polymers in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        if (!enablePolymer)
+            return;
+
+        if (polymerConcentrationOutput_())
+            this->resizeScalarBuffer_(polymerConcentration_);
+        if (polymerDeadPoreVolumeOutput_())
+            this->resizeScalarBuffer_(polymerDeadPoreVolume_);
+        if (polymerRockDensityOutput_())
+            this->resizeScalarBuffer_(polymerRockDensity_);
+        if (polymerAdsorptionOutput_())
+            this->resizeScalarBuffer_(polymerAdsorption_);
+        if (polymerViscosityCorrectionOutput_())
+            this->resizeScalarBuffer_(polymerViscosityCorrection_);
+        if (waterViscosityCorrectionOutput_())
+            this->resizeScalarBuffer_(waterViscosityCorrection_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant for
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        if (!enablePolymer)
+            return;
+
+        for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
+            const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
+            unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
+
+            if (polymerConcentrationOutput_())
+                polymerConcentration_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.polymerConcentration());
+
+            if (polymerDeadPoreVolumeOutput_())
+                polymerDeadPoreVolume_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.polymerDeadPoreVolume());
+
+            if (polymerRockDensityOutput_())
+                polymerRockDensity_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.polymerRockDensity());
+
+            if (polymerAdsorptionOutput_())
+                polymerAdsorption_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.polymerAdsorption());
+
+            if (polymerViscosityCorrectionOutput_())
+                polymerViscosityCorrection_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.polymerViscosityCorrection());
+
+            if (waterViscosityCorrectionOutput_())
+                waterViscosityCorrection_[globalDofIdx] =
+                    Opm::scalarValue(intQuants.waterViscosityCorrection());
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter)
+            return;
+
+        if (!enablePolymer)
+            return;
+
+        if (polymerConcentrationOutput_())
+            this->commitScalarBuffer_(baseWriter, "polymer concentration", polymerConcentration_);
+
+        if (polymerDeadPoreVolumeOutput_())
+            this->commitScalarBuffer_(baseWriter, "dead pore volume fraction", polymerDeadPoreVolume_);
+
+        if (polymerRockDensityOutput_())
+            this->commitScalarBuffer_(baseWriter, "polymer rock density", polymerRockDensity_);
+
+        if (polymerAdsorptionOutput_())
+            this->commitScalarBuffer_(baseWriter, "polymer adsorption", polymerAdsorption_);
+
+        if (polymerViscosityCorrectionOutput_())
+            this->commitScalarBuffer_(baseWriter, "polymer viscosity correction", polymerViscosityCorrection_);
+
+        if (waterViscosityCorrectionOutput_())
+            this->commitScalarBuffer_(baseWriter, "water viscosity correction", waterViscosityCorrection_);
+    }
+
+private:
+    static bool polymerConcentrationOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerConcentration);
+        return val;
+    }
+
+    static bool polymerDeadPoreVolumeOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerDeadPoreVolume);
+        return val;
+    }
+
+    static bool polymerRockDensityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerRockDensity);
+        return val;
+    }
+
+    static bool polymerAdsorptionOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerAdsorption);
+        return val;
+    }
+
+    static bool polymerViscosityCorrectionOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerViscosityCorrection);
+        return val;
+    }
+
+    static bool waterViscosityCorrectionOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerViscosityCorrection);
+        return val;
+    }
+
+    ScalarBuffer polymerConcentration_;
+    ScalarBuffer polymerDeadPoreVolume_;
+    ScalarBuffer polymerRockDensity_;
+    ScalarBuffer polymerAdsorption_;
+    ScalarBuffer polymerViscosityCorrection_;
+    ScalarBuffer waterViscosityCorrection_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/vtkblackoilsolventmodule.hh

@@ -0,0 +1,230 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkBlackOilSolventModule
+ */
+#ifndef EWOMS_VTK_BLACK_OIL_SOLVENT_MODULE_HH
+#define EWOMS_VTK_BLACK_OIL_SOLVENT_MODULE_HH
+
+#include <opm/material/densead/Math.hpp>
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/blackoil/blackoilproperties.hh>
+
+#include <dune/common/fvector.hh>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK multi-phase output
+NEW_TYPE_TAG(VtkBlackOilSolvent);
+
+// create the property tags needed for the solvent output module
+NEW_PROP_TAG(EnableSolvent);
+NEW_PROP_TAG(EnableVtkOutput);
+NEW_PROP_TAG(VtkWriteSolventSaturation);
+NEW_PROP_TAG(VtkWriteSolventDensity);
+NEW_PROP_TAG(VtkWriteSolventViscosity);
+NEW_PROP_TAG(VtkWriteSolventMobility);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkBlackOilSolvent, VtkWriteSolventSaturation, true);
+SET_BOOL_PROP(VtkBlackOilSolvent, VtkWriteSolventDensity, true);
+SET_BOOL_PROP(VtkBlackOilSolvent, VtkWriteSolventViscosity, true);
+SET_BOOL_PROP(VtkBlackOilSolvent, VtkWriteSolventMobility, true);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the black oil model's solvent related quantities.
+ */
+template <class TypeTag>
+class VtkBlackOilSolventModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { enableSolvent = GET_PROP_VALUE(TypeTag, EnableSolvent) };
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+
+public:
+    VtkBlackOilSolventModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the multi-phase VTK output
+     * module.
+     */
+    static void registerParameters()
+    {
+        if (!enableSolvent)
+            return;
+
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSolventSaturation,
+                             "Include the \"saturation\" of the solvent component "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSolventDensity,
+                             "Include the \"density\" of the solvent component "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSolventViscosity,
+                             "Include the \"viscosity\" of the solvent component "
+                             "in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSolventMobility,
+                             "Include the \"mobility\" of the solvent component "
+                             "in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        if (!enableSolvent)
+            return;
+
+        if (solventSaturationOutput_())
+            this->resizeScalarBuffer_(solventSaturation_);
+        if (solventDensityOutput_())
+            this->resizeScalarBuffer_(solventDensity_);
+        if (solventViscosityOutput_())
+            this->resizeScalarBuffer_(solventViscosity_);
+        if (solventMobilityOutput_())
+            this->resizeScalarBuffer_(solventMobility_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant for
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        if (!enableSolvent)
+            return;
+
+        typedef Opm::MathToolbox<Evaluation> Toolbox;
+        for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
+            const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
+            unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
+
+            if (solventSaturationOutput_())
+                solventSaturation_[globalDofIdx] =
+                    Toolbox::scalarValue(intQuants.solventSaturation());
+
+            if (solventDensityOutput_())
+                solventDensity_[globalDofIdx] =
+                    Toolbox::scalarValue(intQuants.solventDensity());
+
+            if (solventViscosityOutput_())
+                solventViscosity_[globalDofIdx] =
+                    Toolbox::scalarValue(intQuants.solventViscosity());
+
+            if (solventMobilityOutput_())
+                solventMobility_[globalDofIdx] =
+                    Toolbox::scalarValue(intQuants.solventMobility());
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter)
+            return;
+
+        if (!enableSolvent)
+            return;
+
+        if (solventSaturationOutput_())
+            this->commitScalarBuffer_(baseWriter, "saturation_solvent", solventSaturation_);
+
+        if (solventDensityOutput_())
+            this->commitScalarBuffer_(baseWriter, "density_solvent", solventDensity_);
+
+        if (solventViscosityOutput_())
+            this->commitScalarBuffer_(baseWriter, "viscosity_solvent", solventViscosity_);
+
+        if (solventMobilityOutput_())
+            this->commitScalarBuffer_(baseWriter, "mobility_solvent", solventMobility_);
+    }
+
+private:
+    static bool solventSaturationOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSolventSaturation);
+        return val;
+    }
+
+    static bool solventDensityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSolventDensity);
+        return val;
+    }
+
+    static bool solventViscosityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSolventViscosity);
+        return val;
+    }
+
+    static bool solventMobilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSolventMobility);
+        return val;
+    }
+
+    ScalarBuffer solventSaturation_;
+    ScalarBuffer solventDensity_;
+    ScalarBuffer solventViscosity_;
+    ScalarBuffer solventMobility_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/vtkcompositionmodule.hh

@@ -0,0 +1,295 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkCompositionModule
+ */
+#ifndef EWOMS_VTK_COMPOSITION_MODULE_HH
+#define EWOMS_VTK_COMPOSITION_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <opm/material/common/MathToolbox.hpp>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK composition output
+NEW_TYPE_TAG(VtkComposition);
+
+// create the property tags needed for the composition module
+NEW_PROP_TAG(VtkWriteMassFractions);
+NEW_PROP_TAG(VtkWriteMoleFractions);
+NEW_PROP_TAG(VtkWriteTotalMassFractions);
+NEW_PROP_TAG(VtkWriteTotalMoleFractions);
+NEW_PROP_TAG(VtkWriteMolarities);
+NEW_PROP_TAG(VtkWriteFugacities);
+NEW_PROP_TAG(VtkWriteFugacityCoeffs);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkComposition, VtkWriteMassFractions, false);
+SET_BOOL_PROP(VtkComposition, VtkWriteMoleFractions, true);
+SET_BOOL_PROP(VtkComposition, VtkWriteTotalMassFractions, false);
+SET_BOOL_PROP(VtkComposition, VtkWriteTotalMoleFractions, false);
+SET_BOOL_PROP(VtkComposition, VtkWriteMolarities, false);
+SET_BOOL_PROP(VtkComposition, VtkWriteFugacities, false);
+SET_BOOL_PROP(VtkComposition, VtkWriteFugacityCoeffs, false);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the fluid composition
+ *
+ * This module deals with the following quantities:
+ * - Mole fraction of a component in a fluid phase
+ * - Mass fraction of a component in a fluid phase
+ * - Molarity (i.e. molar concentration) of a component in a fluid phase
+ * - Fugacity of all components
+ * - FugacityCoefficient of all components in all phases
+ */
+template <class TypeTag>
+class VtkCompositionModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+    enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    typedef typename ParentType::ComponentBuffer ComponentBuffer;
+    typedef typename ParentType::PhaseComponentBuffer PhaseComponentBuffer;
+
+public:
+    VtkCompositionModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the Vtk output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteMassFractions,
+                             "Include mass fractions in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteMoleFractions,
+                             "Include mole fractions in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteTotalMassFractions,
+                             "Include total mass fractions in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteTotalMoleFractions,
+                             "Include total mole fractions in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteMolarities,
+                             "Include component molarities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFugacities,
+                             "Include component fugacities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFugacityCoeffs,
+                             "Include component fugacity coefficients in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (moleFracOutput_())
+            this->resizePhaseComponentBuffer_(moleFrac_);
+        if (massFracOutput_())
+            this->resizePhaseComponentBuffer_(massFrac_);
+        if (totalMassFracOutput_())
+            this->resizeComponentBuffer_(totalMassFrac_);
+        if (totalMoleFracOutput_())
+            this->resizeComponentBuffer_(totalMoleFrac_);
+        if (molarityOutput_())
+            this->resizePhaseComponentBuffer_(molarity_);
+
+        if (fugacityOutput_())
+            this->resizeComponentBuffer_(fugacity_);
+        if (fugacityCoeffOutput_())
+            this->resizePhaseComponentBuffer_(fugacityCoeff_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant
+     *        for an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        typedef Opm::MathToolbox<Evaluation> Toolbox;
+
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
+            const auto& fs = intQuants.fluidState();
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
+                    if (moleFracOutput_())
+                        moleFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
+                    if (massFracOutput_())
+                        massFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.massFraction(phaseIdx, compIdx));
+                    if (molarityOutput_())
+                        molarity_[phaseIdx][compIdx][I] = Toolbox::value(fs.molarity(phaseIdx, compIdx));
+
+                    if (fugacityCoeffOutput_())
+                        fugacityCoeff_[phaseIdx][compIdx][I] =
+                            Toolbox::value(fs.fugacityCoefficient(phaseIdx, compIdx));
+                }
+            }
+
+            for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
+                if (totalMassFracOutput_()) {
+                    Scalar compMass = 0;
+                    Scalar totalMass = 0;
+                    for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                        totalMass += Toolbox::value(fs.density(phaseIdx)) * Toolbox::value(fs.saturation(phaseIdx));
+                        compMass +=
+                            Toolbox::value(fs.density(phaseIdx))
+                            *Toolbox::value(fs.saturation(phaseIdx))
+                            *Toolbox::value(fs.massFraction(phaseIdx, compIdx));
+                    }
+                    totalMassFrac_[compIdx][I] = compMass / totalMass;
+                }
+                if (totalMoleFracOutput_()) {
+                    Scalar compMoles = 0;
+                    Scalar totalMoles = 0;
+                    for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                        totalMoles +=
+                            Toolbox::value(fs.molarDensity(phaseIdx))
+                            *Toolbox::value(fs.saturation(phaseIdx));
+                        compMoles +=
+                            Toolbox::value(fs.molarDensity(phaseIdx))
+                            *Toolbox::value(fs.saturation(phaseIdx))
+                            *Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
+                    }
+                    totalMoleFrac_[compIdx][I] = compMoles / totalMoles;
+                }
+                if (fugacityOutput_())
+                    fugacity_[compIdx][I] = Toolbox::value(intQuants.fluidState().fugacity(/*phaseIdx=*/0, compIdx));
+            }
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (moleFracOutput_())
+            this->commitPhaseComponentBuffer_(baseWriter, "moleFrac_%s^%s", moleFrac_);
+        if (massFracOutput_())
+            this->commitPhaseComponentBuffer_(baseWriter, "massFrac_%s^%s", massFrac_);
+        if (molarityOutput_())
+            this->commitPhaseComponentBuffer_(baseWriter, "molarity_%s^%s", molarity_);
+        if (totalMassFracOutput_())
+            this->commitComponentBuffer_(baseWriter, "totalMassFrac^%s", totalMassFrac_);
+        if (totalMoleFracOutput_())
+            this->commitComponentBuffer_(baseWriter, "totalMoleFrac^%s", totalMoleFrac_);
+
+        if (fugacityOutput_())
+            this->commitComponentBuffer_(baseWriter, "fugacity^%s", fugacity_);
+        if (fugacityCoeffOutput_())
+            this->commitPhaseComponentBuffer_(baseWriter, "fugacityCoeff_%s^%s", fugacityCoeff_);
+    }
+
+private:
+    static bool massFracOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteMassFractions);
+        return val;
+    }
+
+    static bool moleFracOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteMoleFractions);
+        return val;
+    }
+
+    static bool totalMassFracOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteTotalMassFractions);
+        return val;
+    }
+
+    static bool totalMoleFracOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteTotalMoleFractions);
+        return val;
+    }
+
+    static bool molarityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteMolarities);
+        return val;
+    }
+
+    static bool fugacityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFugacities);
+        return val;
+    }
+
+    static bool fugacityCoeffOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFugacityCoeffs);
+        return val;
+    }
+
+    PhaseComponentBuffer moleFrac_;
+    PhaseComponentBuffer massFrac_;
+    PhaseComponentBuffer molarity_;
+    ComponentBuffer totalMassFrac_;
+    ComponentBuffer totalMoleFrac_;
+
+    ComponentBuffer fugacity_;
+    PhaseComponentBuffer fugacityCoeff_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkdiffusionmodule.hh

@@ -0,0 +1,204 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ *
+ * \copydoc Opm::VtkDiffusionModule
+ */
+#ifndef EWOMS_VTK_DIFFUSION_MODULE_HH
+#define EWOMS_VTK_DIFFUSION_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <opm/material/densead/Evaluation.hpp>
+#include <opm/material/densead/Math.hpp>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK output of the quantities for molecular
+// diffusion
+NEW_TYPE_TAG(VtkDiffusion);
+
+// create the property tags needed for the diffusion module
+NEW_PROP_TAG(VtkWriteTortuosities);
+NEW_PROP_TAG(VtkWriteDiffusionCoefficients);
+NEW_PROP_TAG(VtkWriteEffectiveDiffusionCoefficients);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkDiffusion, VtkWriteTortuosities, false);
+SET_BOOL_PROP(VtkDiffusion, VtkWriteDiffusionCoefficients, false);
+SET_BOOL_PROP(VtkDiffusion, VtkWriteEffectiveDiffusionCoefficients, false);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for quantities which make sense for models which
+ *        incorperate molecular diffusion.
+ *
+ * This module deals with the following quantities:
+ * - Molecular diffusion coefficients of all components in all fluid phases
+ * - Effective molecular diffusion coefficients of the porous medium of all
+ *   components in all fluid phases
+ */
+template <class TypeTag>
+class VtkDiffusionModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+
+    typedef Opm::MathToolbox<Evaluation> Toolbox;
+
+    typedef typename ParentType::PhaseComponentBuffer PhaseComponentBuffer;
+    typedef typename ParentType::PhaseBuffer PhaseBuffer;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+    enum { numComponents = GET_PROP_VALUE(TypeTag, NumComponents) };
+
+public:
+    VtkDiffusionModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the Vtk output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteTortuosities,
+                             "Include the tortuosity for each phase in the VTK "
+                             "output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteDiffusionCoefficients,
+                             "Include the molecular diffusion coefficients in "
+                             "the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool,
+                             VtkWriteEffectiveDiffusionCoefficients,
+                             "Include the effective molecular diffusion "
+                             "coefficients the medium in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (tortuosityOutput_())
+            this->resizePhaseBuffer_(tortuosity_);
+        if (diffusionCoefficientOutput_())
+            this->resizePhaseComponentBuffer_(diffusionCoefficient_);
+        if (effectiveDiffusionCoefficientOutput_())
+            this->resizePhaseComponentBuffer_(effectiveDiffusionCoefficient_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quanties relevant
+     *        for an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                if (tortuosityOutput_())
+                    tortuosity_[phaseIdx][I] = Toolbox::value(intQuants.tortuosity(phaseIdx));
+                for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
+                    if (diffusionCoefficientOutput_())
+                        diffusionCoefficient_[phaseIdx][compIdx][I] =
+                            Toolbox::value(intQuants.diffusionCoefficient(phaseIdx, compIdx));
+                    if (effectiveDiffusionCoefficientOutput_())
+                        effectiveDiffusionCoefficient_[phaseIdx][compIdx][I] =
+                            Toolbox::value(intQuants.effectiveDiffusionCoefficient(phaseIdx, compIdx));
+                }
+            }
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (tortuosityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "tortuosity", tortuosity_);
+        if (diffusionCoefficientOutput_())
+            this->commitPhaseComponentBuffer_(baseWriter, "diffusionCoefficient",
+                                              diffusionCoefficient_);
+        if (effectiveDiffusionCoefficientOutput_())
+            this->commitPhaseComponentBuffer_(baseWriter,
+                                              "effectiveDiffusionCoefficient",
+                                              effectiveDiffusionCoefficient_);
+    }
+
+private:
+    static bool tortuosityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteTortuosities);
+        return val;
+    }
+
+    static bool diffusionCoefficientOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteDiffusionCoefficients);
+        return val;
+    }
+
+    static bool effectiveDiffusionCoefficientOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteEffectiveDiffusionCoefficients);
+        return val;
+    }
+
+    PhaseBuffer tortuosity_;
+    PhaseComponentBuffer diffusionCoefficient_;
+    PhaseComponentBuffer effectiveDiffusionCoefficient_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkdiscretefracturemodule.hh

@@ -0,0 +1,364 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkDiscreteFractureModule
+ */
+#ifndef EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
+#define EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <opm/material/common/Valgrind.hpp>
+
+#include <dune/common/fvector.hh>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK multi-phase output
+NEW_TYPE_TAG(VtkDiscreteFracture);
+
+// create the property tags needed for the multi phase module
+NEW_PROP_TAG(Vanguard);
+NEW_PROP_TAG(VtkWriteFractureSaturations);
+NEW_PROP_TAG(VtkWriteFractureMobilities);
+NEW_PROP_TAG(VtkWriteFractureRelativePermeabilities);
+NEW_PROP_TAG(VtkWriteFracturePorosity);
+NEW_PROP_TAG(VtkWriteFractureIntrinsicPermeabilities);
+NEW_PROP_TAG(VtkWriteFractureFilterVelocities);
+NEW_PROP_TAG(VtkWriteFractureVolumeFraction);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+NEW_PROP_TAG(DiscBaseOutputModule);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureSaturations, true);
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureMobilities, false);
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureRelativePermeabilities, true);
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFracturePorosity, true);
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureIntrinsicPermeabilities, false);
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureFilterVelocities, false);
+SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureVolumeFraction, true);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for quantities which make sense for all
+ *        models which deal with discrete fractures in porous media.
+ *
+ * This module deals with the following quantities:
+ * - Saturations of all fluid phases in the fracture
+ * - Mobilities of all fluid phases in the fracture
+ * - Relative permeabilities of all fluid phases in the fracture
+ * - Porosity of the medium in the fracture
+ * - Norm of the intrinsic permeability of the medium in the fracture
+ */
+template <class TypeTag>
+class VtkDiscreteFractureModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Vanguard) Vanguard;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+
+    typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { dim = GridView::dimension };
+    enum { dimWorld = GridView::dimensionworld };
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+    typedef typename ParentType::PhaseBuffer PhaseBuffer;
+    typedef typename ParentType::PhaseVectorBuffer PhaseVectorBuffer;
+
+public:
+    VtkDiscreteFractureModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the multi-phase VTK output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureSaturations,
+                             "Include the phase saturations in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureMobilities,
+                             "Include the phase mobilities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureRelativePermeabilities,
+                             "Include the phase relative permeabilities in the "
+                             "VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFracturePorosity,
+                             "Include the porosity in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureIntrinsicPermeabilities,
+                             "Include the intrinsic permeability in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureFilterVelocities,
+                             "Include in the filter velocities of the phases "
+                             "the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureVolumeFraction,
+                             "Add the fraction of the total volume which is "
+                             "occupied by fractures in the VTK output");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (saturationOutput_())
+            this->resizePhaseBuffer_(fractureSaturation_);
+        if (mobilityOutput_())
+            this->resizePhaseBuffer_(fractureMobility_);
+        if (relativePermeabilityOutput_())
+            this->resizePhaseBuffer_(fractureRelativePermeability_);
+
+        if (porosityOutput_())
+            this->resizeScalarBuffer_(fracturePorosity_);
+        if (intrinsicPermeabilityOutput_())
+            this->resizeScalarBuffer_(fractureIntrinsicPermeability_);
+        if (volumeFractionOutput_())
+            this->resizeScalarBuffer_(fractureVolumeFraction_);
+
+        if (velocityOutput_()) {
+            size_t nDof = this->simulator_.model().numGridDof();
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                fractureVelocity_[phaseIdx].resize(nDof);
+                for (unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
+                    fractureVelocity_[phaseIdx][dofIdx].resize(dimWorld);
+                    fractureVelocity_[phaseIdx][dofIdx] = 0.0;
+                }
+            }
+            this->resizePhaseBuffer_(fractureVelocityWeight_);
+        }
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant for
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        const auto& fractureMapper = elemCtx.simulator().vanguard().fractureMapper();
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            if (!fractureMapper.isFractureVertex(I))
+                continue;
+
+            const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
+            const auto& fs = intQuants.fractureFluidState();
+
+            if (porosityOutput_()) {
+                Opm::Valgrind::CheckDefined(intQuants.fracturePorosity());
+                fracturePorosity_[I] = intQuants.fracturePorosity();
+            }
+            if (intrinsicPermeabilityOutput_()) {
+                const auto& K = intQuants.fractureIntrinsicPermeability();
+                fractureIntrinsicPermeability_[I] = K[0][0];
+            }
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                if (saturationOutput_()) {
+                    Opm::Valgrind::CheckDefined(fs.saturation(phaseIdx));
+                    fractureSaturation_[phaseIdx][I] = fs.saturation(phaseIdx);
+                }
+                if (mobilityOutput_()) {
+                    Opm::Valgrind::CheckDefined(intQuants.fractureMobility(phaseIdx));
+                    fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
+                }
+                if (relativePermeabilityOutput_()) {
+                    Opm::Valgrind::CheckDefined(intQuants.fractureRelativePermeability(phaseIdx));
+                    fractureRelativePermeability_[phaseIdx][I] =
+                        intQuants.fractureRelativePermeability(phaseIdx);
+                }
+                if (volumeFractionOutput_()) {
+                    Opm::Valgrind::CheckDefined(intQuants.fractureVolume());
+                    fractureVolumeFraction_[I] += intQuants.fractureVolume();
+                }
+            }
+        }
+
+        if (velocityOutput_()) {
+            // calculate velocities if requested by the simulator
+            for (unsigned scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++ scvfIdx) {
+                const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, /*timeIdx=*/0);
+
+                unsigned i = extQuants.interiorIndex();
+                unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+
+                unsigned j = extQuants.exteriorIndex();
+                unsigned J = elemCtx.globalSpaceIndex(j, /*timeIdx=*/0);
+
+                if (!fractureMapper.isFractureEdge(I, J))
+                    continue;
+
+                for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                    Scalar weight =
+                        std::max<Scalar>(1e-16, std::abs(extQuants.fractureVolumeFlux(phaseIdx)));
+                    Opm::Valgrind::CheckDefined(extQuants.extrusionFactor());
+                    assert(extQuants.extrusionFactor() > 0);
+                    weight *= extQuants.extrusionFactor();
+
+                    Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
+                    v *= weight;
+
+                    for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
+                        fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
+                        fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
+                    }
+
+                    fractureVelocityWeight_[phaseIdx][I] += weight;
+                    fractureVelocityWeight_[phaseIdx][J] += weight;
+                }
+            }
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (saturationOutput_())
+            this->commitPhaseBuffer_(baseWriter, "fractureSaturation_%s", fractureSaturation_);
+        if (mobilityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "fractureMobility_%s", fractureMobility_);
+        if (relativePermeabilityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "fractureRelativePerm_%s", fractureRelativePermeability_);
+
+        if (porosityOutput_())
+            this->commitScalarBuffer_(baseWriter, "fracturePorosity", fracturePorosity_);
+        if (intrinsicPermeabilityOutput_())
+            this->commitScalarBuffer_(baseWriter, "fractureIntrinsicPerm", fractureIntrinsicPermeability_);
+        if (volumeFractionOutput_()) {
+            // divide the fracture volume by the total volume of the finite volumes
+            for (unsigned I = 0; I < fractureVolumeFraction_.size(); ++I)
+                fractureVolumeFraction_[I] /= this->simulator_.model().dofTotalVolume(I);
+            this->commitScalarBuffer_(baseWriter, "fractureVolumeFraction", fractureVolumeFraction_);
+        }
+
+        if (velocityOutput_()) {
+            size_t nDof = this->simulator_.model().numGridDof();
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                // first, divide the velocity field by the
+                // respective finite volume's surface area
+                for (unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx)
+                    fractureVelocity_[phaseIdx][dofIdx] /=
+                        std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
+                // commit the phase velocity
+                char name[512];
+                snprintf(name, 512, "fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx));
+
+                DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
+            }
+        }
+    }
+
+private:
+    static bool saturationOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureSaturations);
+        return val;
+    }
+
+    static bool mobilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureMobilities);
+        return val;
+    }
+
+    static bool relativePermeabilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureRelativePermeabilities);
+        return val;
+    }
+
+    static bool porosityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFracturePorosity);
+        return val;
+    }
+
+    static bool intrinsicPermeabilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureIntrinsicPermeabilities);
+        return val;
+    }
+
+    static bool volumeFractionOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureVolumeFraction);
+        return val;
+    }
+
+    static bool velocityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureFilterVelocities);
+        return val;
+    }
+
+    PhaseBuffer fractureSaturation_;
+    PhaseBuffer fractureMobility_;
+    PhaseBuffer fractureRelativePermeability_;
+
+    ScalarBuffer fracturePorosity_;
+    ScalarBuffer fractureVolumeFraction_;
+    ScalarBuffer fractureIntrinsicPermeability_;
+
+    PhaseVectorBuffer fractureVelocity_;
+    PhaseBuffer fractureVelocityWeight_;
+
+    PhaseVectorBuffer potentialGradient_;
+    PhaseBuffer potentialWeight_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkenergymodule.hh

@@ -0,0 +1,218 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkEnergyModule
+ */
+#ifndef EWOMS_VTK_ENERGY_MODULE_HH
+#define EWOMS_VTK_ENERGY_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <opm/material/common/MathToolbox.hpp>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK energy output
+NEW_TYPE_TAG(VtkEnergy);
+
+// create the property tags needed for the energy module
+NEW_PROP_TAG(VtkWriteSolidInternalEnergy);
+NEW_PROP_TAG(VtkWriteThermalConductivity);
+NEW_PROP_TAG(VtkWriteInternalEnergies);
+NEW_PROP_TAG(VtkWriteEnthalpies);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkEnergy, VtkWriteSolidInternalEnergy, false);
+SET_BOOL_PROP(VtkEnergy, VtkWriteThermalConductivity, false);
+SET_BOOL_PROP(VtkEnergy, VtkWriteInternalEnergies, false);
+SET_BOOL_PROP(VtkEnergy, VtkWriteEnthalpies, false);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for quantities which make sense for models which
+ *        assume thermal equilibrium.
+ *
+ * This module deals with the following quantities:
+ * - Specific enthalpy of all fluid phases
+ * - Specific internal energy of all fluid phases
+ * - Volumetric internal energy of the solid phase
+ * - Total thermal conductivity, i.e. the conductivity of the solid and all fluid phases
+ *   combined
+ */
+template <class TypeTag>
+class VtkEnergyModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+    typedef typename ParentType::PhaseBuffer PhaseBuffer;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+
+    typedef typename Opm::MathToolbox<Evaluation> Toolbox;
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+public:
+    VtkEnergyModule(const Simulator& simulator)
+        : ParentType(simulator)
+    {
+    }
+
+    /*!
+     * \brief Register all run-time parameters for the Vtk output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSolidInternalEnergy,
+                             "Include the volumetric internal energy of solid"
+                             "matrix in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteThermalConductivity,
+                             "Include the total thermal conductivity of the"
+                             "medium in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteEnthalpies,
+                             "Include the specific enthalpy of the phases in "
+                             "the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteInternalEnergies,
+                             "Include the specific internal energy of the "
+                             "phases in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (enthalpyOutput_())
+            this->resizePhaseBuffer_(enthalpy_);
+        if (internalEnergyOutput_())
+            this->resizePhaseBuffer_(internalEnergy_);
+
+        if (solidInternalEnergyOutput_())
+            this->resizeScalarBuffer_(solidInternalEnergy_);
+        if (thermalConductivityOutput_())
+            this->resizeScalarBuffer_(thermalConductivity_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quanties relevant
+     *        for an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
+            const auto& fs = intQuants.fluidState();
+
+            if (solidInternalEnergyOutput_())
+                solidInternalEnergy_[I] = Toolbox::value(intQuants.solidInternalEnergy());
+            if (thermalConductivityOutput_())
+                thermalConductivity_[I] = Toolbox::value(intQuants.thermalConductivity());
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                if (enthalpyOutput_())
+                    enthalpy_[phaseIdx][I] = Toolbox::value(fs.enthalpy(phaseIdx));
+                if (internalEnergyOutput_())
+                    internalEnergy_[phaseIdx][I] = Toolbox::value(fs.internalEnergy(phaseIdx));
+            }
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (solidInternalEnergyOutput_())
+            this->commitScalarBuffer_(baseWriter, "internalEnergySolid", solidInternalEnergy_);
+        if (thermalConductivityOutput_())
+            this->commitScalarBuffer_(baseWriter, "thermalConductivity", thermalConductivity_);
+
+        if (enthalpyOutput_())
+            this->commitPhaseBuffer_(baseWriter, "enthalpy_%s", enthalpy_);
+        if (internalEnergyOutput_())
+            this->commitPhaseBuffer_(baseWriter, "internalEnergy_%s", internalEnergy_);
+    }
+
+private:
+    static bool solidInternalEnergyOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSolidInternalEnergy);
+        return val;
+    }
+
+    static bool thermalConductivityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteThermalConductivity);
+        return val;
+    }
+
+    static bool enthalpyOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteEnthalpies);
+        return val;
+    }
+
+    static bool internalEnergyOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteInternalEnergies);
+        return val;
+    }
+
+    PhaseBuffer enthalpy_;
+    PhaseBuffer internalEnergy_;
+
+    ScalarBuffer thermalConductivity_;
+    ScalarBuffer solidInternalEnergy_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkmultiphasemodule.hh

@@ -0,0 +1,482 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkMultiPhaseModule
+ */
+#ifndef EWOMS_VTK_MULTI_PHASE_MODULE_HH
+#define EWOMS_VTK_MULTI_PHASE_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/propertysystem.hh>
+#include <opm/models/utils/parametersystem.hh>
+
+#include <opm/material/common/MathToolbox.hpp>
+#include <opm/material/common/Valgrind.hpp>
+
+#include <dune/common/fvector.hh>
+
+#include <cstdio>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK multi-phase output
+NEW_TYPE_TAG(VtkMultiPhase);
+
+// create the property tags needed for the multi phase module
+NEW_PROP_TAG(VtkWriteExtrusionFactor);
+NEW_PROP_TAG(VtkWritePressures);
+NEW_PROP_TAG(VtkWriteDensities);
+NEW_PROP_TAG(VtkWriteSaturations);
+NEW_PROP_TAG(VtkWriteMobilities);
+NEW_PROP_TAG(VtkWriteRelativePermeabilities);
+NEW_PROP_TAG(VtkWriteViscosities);
+NEW_PROP_TAG(VtkWriteAverageMolarMasses);
+NEW_PROP_TAG(VtkWritePorosity);
+NEW_PROP_TAG(VtkWriteIntrinsicPermeabilities);
+NEW_PROP_TAG(VtkWritePotentialGradients);
+NEW_PROP_TAG(VtkWriteFilterVelocities);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteExtrusionFactor, false);
+SET_BOOL_PROP(VtkMultiPhase, VtkWritePressures, true);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteDensities, true);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteSaturations, true);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteMobilities, false);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteRelativePermeabilities, true);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteViscosities, false);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteAverageMolarMasses, false);
+SET_BOOL_PROP(VtkMultiPhase, VtkWritePorosity, true);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteIntrinsicPermeabilities, false);
+SET_BOOL_PROP(VtkMultiPhase, VtkWritePotentialGradients, false);
+SET_BOOL_PROP(VtkMultiPhase, VtkWriteFilterVelocities, false);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for quantities which make sense for all
+ *        models which deal with multiple fluid phases in porous media
+ *        that don't use flashy concepts like interfacial area.
+ *
+ * This module deals with the following quantities:
+ * - Pressures of all fluid phases
+ * - Densities of all fluid phases
+ * - Saturations of all fluid phases
+ * - Mobilities of all fluid phases
+ * - Relative permeabilities of all fluid phases
+ * - Viscosities of all fluid phases
+ * - Average molar masses of all fluid phases
+ * - Porosity of the medium
+ * - Norm of the intrinsic permeability of the medium
+ */
+template<class TypeTag>
+class VtkMultiPhaseModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+    typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    enum { dimWorld = GridView::dimensionworld };
+    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+    typedef typename ParentType::VectorBuffer VectorBuffer;
+    typedef typename ParentType::TensorBuffer TensorBuffer;
+    typedef typename ParentType::PhaseBuffer PhaseBuffer;
+
+    typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
+
+    typedef std::array<VectorBuffer, numPhases> PhaseVectorBuffer;
+
+public:
+    VtkMultiPhaseModule(const Simulator& simulator)
+        : ParentType(simulator)
+    {}
+
+    /*!
+     * \brief Register all run-time parameters for the multi-phase VTK output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteExtrusionFactor,
+                             "Include the extrusion factor of the degrees of freedom into the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePressures,
+                             "Include the phase pressures in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteDensities,
+                             "Include the phase densities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteSaturations,
+                             "Include the phase saturations in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteMobilities,
+                             "Include the phase mobilities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteRelativePermeabilities,
+                             "Include the phase relative permeabilities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteViscosities,
+                             "Include component phase viscosities in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteAverageMolarMasses,
+                             "Include the average phase mass in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePorosity,
+                             "Include the porosity in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteIntrinsicPermeabilities,
+                             "Include the intrinsic permeability in the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFilterVelocities,
+                             "Include in the filter velocities of the phases the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePotentialGradients,
+                             "Include the phase pressure potential gradients in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (extrusionFactorOutput_()) this->resizeScalarBuffer_(extrusionFactor_);
+        if (pressureOutput_()) this->resizePhaseBuffer_(pressure_);
+        if (densityOutput_()) this->resizePhaseBuffer_(density_);
+        if (saturationOutput_()) this->resizePhaseBuffer_(saturation_);
+        if (mobilityOutput_()) this->resizePhaseBuffer_(mobility_);
+        if (relativePermeabilityOutput_()) this->resizePhaseBuffer_(relativePermeability_);
+        if (viscosityOutput_()) this->resizePhaseBuffer_(viscosity_);
+        if (averageMolarMassOutput_()) this->resizePhaseBuffer_(averageMolarMass_);
+
+        if (porosityOutput_()) this->resizeScalarBuffer_(porosity_);
+        if (intrinsicPermeabilityOutput_()) this->resizeTensorBuffer_(intrinsicPermeability_);
+
+        if (velocityOutput_()) {
+            size_t nDof = this->simulator_.model().numGridDof();
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
+                velocity_[phaseIdx].resize(nDof);
+                for (unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
+                    velocity_[phaseIdx][dofIdx].resize(dimWorld);
+                    velocity_[phaseIdx][dofIdx] = 0.0;
+                }
+            }
+            this->resizePhaseBuffer_(velocityWeight_);
+        }
+
+        if (potentialGradientOutput_()) {
+            size_t nDof = this->simulator_.model().numGridDof();
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
+                potentialGradient_[phaseIdx].resize(nDof);
+                for (unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
+                    potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
+                    potentialGradient_[phaseIdx][dofIdx] = 0.0;
+                }
+            }
+
+            this->resizePhaseBuffer_(potentialWeight_);
+        }
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities seen on
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        const auto& problem = elemCtx.problem();
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
+            const auto& fs = intQuants.fluidState();
+
+            if (extrusionFactorOutput_()) extrusionFactor_[I] = intQuants.extrusionFactor();
+            if (porosityOutput_()) porosity_[I] = Opm::getValue(intQuants.porosity());
+
+            if (intrinsicPermeabilityOutput_()) {
+                const auto& K = problem.intrinsicPermeability(elemCtx, i, /*timeIdx=*/0);
+                for (unsigned rowIdx = 0; rowIdx < K.rows; ++rowIdx)
+                    for (unsigned colIdx = 0; colIdx < K.cols; ++colIdx)
+                        intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
+            }
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                if (pressureOutput_())
+                    pressure_[phaseIdx][I] = Opm::getValue(fs.pressure(phaseIdx));
+                if (densityOutput_())
+                    density_[phaseIdx][I] = Opm::getValue(fs.density(phaseIdx));
+                if (saturationOutput_())
+                    saturation_[phaseIdx][I] = Opm::getValue(fs.saturation(phaseIdx));
+                if (mobilityOutput_())
+                    mobility_[phaseIdx][I] = Opm::getValue(intQuants.mobility(phaseIdx));
+                if (relativePermeabilityOutput_())
+                    relativePermeability_[phaseIdx][I] = Opm::getValue(intQuants.relativePermeability(phaseIdx));
+                if (viscosityOutput_())
+                    viscosity_[phaseIdx][I] = Opm::getValue(fs.viscosity(phaseIdx));
+                if (averageMolarMassOutput_())
+                    averageMolarMass_[phaseIdx][I] = Opm::getValue(fs.averageMolarMass(phaseIdx));
+            }
+        }
+
+        if (potentialGradientOutput_()) {
+            // calculate velocities if requested
+            for (unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++ faceIdx) {
+                const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, /*timeIdx=*/0);
+
+                unsigned i = extQuants.interiorIndex();
+                unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+
+                for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                    Scalar weight = extQuants.extrusionFactor();
+
+                    potentialWeight_[phaseIdx][I] += weight;
+
+                    const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
+                    DimVector pGrad;
+                    for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
+                        pGrad[dimIdx] = Opm::getValue(inputPGrad[dimIdx])*weight;
+                    potentialGradient_[phaseIdx][I] += pGrad;
+                } // end for all phases
+            } // end for all faces
+        }
+
+        if (velocityOutput_()) {
+            // calculate velocities if requested
+            for (unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++ faceIdx) {
+                const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, /*timeIdx=*/0);
+
+                unsigned i = extQuants.interiorIndex();
+                unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+
+                unsigned j = extQuants.exteriorIndex();
+                unsigned J = elemCtx.globalSpaceIndex(j, /*timeIdx=*/0);
+
+                for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                    Scalar weight = std::max<Scalar>(1e-16,
+                                                     std::abs(Opm::getValue(extQuants.volumeFlux(phaseIdx))));
+                    Opm::Valgrind::CheckDefined(extQuants.extrusionFactor());
+                    assert(extQuants.extrusionFactor() > 0);
+                    weight *= extQuants.extrusionFactor();
+
+                    const auto& inputV = extQuants.filterVelocity(phaseIdx);
+                    DimVector v;
+                    for (unsigned k = 0; k < dimWorld; ++k)
+                        v[k] = Opm::getValue(inputV[k]);
+                    if (v.two_norm() > 1e-20)
+                        weight /= v.two_norm();
+                    v *= weight;
+
+                    velocity_[phaseIdx][I] += v;
+                    velocity_[phaseIdx][J] += v;
+
+                    velocityWeight_[phaseIdx][I] += weight;
+                    velocityWeight_[phaseIdx][J] += weight;
+                } // end for all phases
+            } // end for all faces
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter)
+            return;
+
+        if (extrusionFactorOutput_())
+            this->commitScalarBuffer_(baseWriter, "extrusionFactor", extrusionFactor_);
+        if (pressureOutput_())
+            this->commitPhaseBuffer_(baseWriter, "pressure_%s", pressure_);
+        if (densityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "density_%s", density_);
+        if (saturationOutput_())
+            this->commitPhaseBuffer_(baseWriter, "saturation_%s", saturation_);
+        if (mobilityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "mobility_%s", mobility_);
+        if (relativePermeabilityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "relativePerm_%s", relativePermeability_);
+        if (viscosityOutput_())
+            this->commitPhaseBuffer_(baseWriter, "viscosity_%s", viscosity_);
+        if (averageMolarMassOutput_())
+            this->commitPhaseBuffer_(baseWriter, "averageMolarMass_%s", averageMolarMass_);
+
+        if (porosityOutput_())
+            this->commitScalarBuffer_(baseWriter, "porosity", porosity_);
+        if (intrinsicPermeabilityOutput_())
+            this->commitTensorBuffer_(baseWriter, "intrinsicPerm", intrinsicPermeability_);
+
+        if (velocityOutput_()) {
+            size_t numDof = this->simulator_.model().numGridDof();
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                // first, divide the velocity field by the
+                // respective finite volume's surface area
+                for (unsigned i = 0; i < numDof; ++i)
+                    velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
+                // commit the phase velocity
+                char name[512];
+                snprintf(name, 512, "filterVelocity_%s", FluidSystem::phaseName(phaseIdx));
+
+                DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
+            }
+        }
+
+        if (potentialGradientOutput_()) {
+            size_t numDof = this->simulator_.model().numGridDof();
+
+            for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+                // first, divide the velocity field by the
+                // respective finite volume's surface area
+                for (unsigned i = 0; i < numDof; ++i)
+                    potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
+                // commit the phase velocity
+                char name[512];
+                snprintf(name, 512, "gradP_%s", FluidSystem::phaseName(phaseIdx));
+
+                DiscBaseOutputModule::attachVectorDofData_(baseWriter,
+                                                           potentialGradient_[phaseIdx],
+                                                           name);
+            }
+        }
+    }
+
+    /*!
+     * \brief Returns true iff the module needs to access the extensive quantities of a
+     * context to do its job.
+     *
+     * For example, this happens if velocities or gradients should be written. Always
+     * returning true here does not do any harm from the correctness perspective, but it
+     * slows down writing the output fields.
+     */
+    virtual bool needExtensiveQuantities() const final
+    {
+        return velocityOutput_() || potentialGradientOutput_();
+    }
+
+private:
+    static bool extrusionFactorOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteExtrusionFactor);
+        return val;
+    }
+
+    static bool pressureOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePressures);
+        return val;
+    }
+
+    static bool densityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteDensities);
+        return val;
+    }
+
+    static bool saturationOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteSaturations);
+        return val;
+    }
+
+    static bool mobilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteMobilities);
+        return val;
+    }
+
+    static bool relativePermeabilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteRelativePermeabilities);
+        return val;
+    }
+
+    static bool viscosityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteViscosities);
+        return val;
+    }
+
+    static bool averageMolarMassOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteAverageMolarMasses);
+        return val;
+    }
+
+    static bool porosityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePorosity);
+        return val;
+    }
+
+    static bool intrinsicPermeabilityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteIntrinsicPermeabilities);
+        return val;
+    }
+
+    static bool velocityOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFilterVelocities);
+        return val;
+    }
+
+    static bool potentialGradientOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePotentialGradients);
+        return val;
+    }
+
+    ScalarBuffer extrusionFactor_;
+    PhaseBuffer pressure_;
+    PhaseBuffer density_;
+    PhaseBuffer saturation_;
+    PhaseBuffer mobility_;
+    PhaseBuffer relativePermeability_;
+    PhaseBuffer viscosity_;
+    PhaseBuffer averageMolarMass_;
+
+    ScalarBuffer porosity_;
+    TensorBuffer intrinsicPermeability_;
+
+    PhaseVectorBuffer velocity_;
+    PhaseBuffer velocityWeight_;
+
+    PhaseVectorBuffer potentialGradient_;
+    PhaseBuffer potentialWeight_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkmultiwriter.hh

@@ -0,0 +1,575 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ *
+ * \copydoc Opm::VtkMultiWriter
+ */
+#ifndef EWOMS_VTK_MULTI_WRITER_HH
+#define EWOMS_VTK_MULTI_WRITER_HH
+
+#include "vtkscalarfunction.hh"
+#include "vtkvectorfunction.hh"
+#include "vtktensorfunction.hh"
+
+#include <opm/models/io/baseoutputwriter.hh>
+#include <opm/models/parallel/tasklets.hh>
+
+#include <opm/material/common/Valgrind.hpp>
+#include <opm/material/common/Unused.hpp>
+
+#include <dune/common/fvector.hh>
+#include <dune/istl/bvector.hh>
+#include <dune/grid/io/file/vtk/vtkwriter.hh>
+
+#if HAVE_MPI
+#include <mpi.h>
+#endif
+
+#include <list>
+#include <string>
+#include <limits>
+#include <sstream>
+#include <fstream>
+
+namespace Opm {
+/*!
+ * \brief Simplifies writing multi-file VTK datasets.
+ *
+ * This class automatically keeps the meta file up to date and
+ * simplifies writing datasets consisting of multiple files. (i.e.
+ * multiple time steps or grid refinements within a time step.)
+ */
+template <class GridView, int vtkFormat>
+class VtkMultiWriter : public BaseOutputWriter
+{
+    class WriteDataTasklet : public TaskletInterface
+    {
+    public:
+        WriteDataTasklet(VtkMultiWriter& multiWriter)
+            : multiWriter_(multiWriter)
+        { }
+
+        void run() final
+        {
+            std::string fileName;
+            // write the actual data as vtu or vtp (plus the pieces file in the parallel case)
+            if (multiWriter_.commSize_ > 1)
+                fileName = multiWriter_.curWriter_->pwrite(/*name=*/multiWriter_.curOutFileName_,
+                                                           /*path=*/multiWriter_.outputDir_,
+                                                           /*extendPath=*/"",
+                                                           static_cast<Dune::VTK::OutputType>(vtkFormat));
+            else
+                fileName = multiWriter_.curWriter_->write(/*name=*/multiWriter_.outputDir_ + "/" + multiWriter_.curOutFileName_,
+                                                          static_cast<Dune::VTK::OutputType>(vtkFormat));
+
+            // determine name to write into the multi-file for the
+            // current time step
+            multiWriter_.multiFile_.precision(16);
+            multiWriter_.multiFile_ << "   <DataSet timestep=\"" << multiWriter_.curTime_ << "\" file=\""
+                                    << fileName << "\"/>\n";
+        }
+
+    private:
+        VtkMultiWriter& multiWriter_;
+    };
+
+    enum { dim = GridView::dimension };
+
+#if DUNE_VERSION_NEWER(DUNE_GRID, 2,6)
+    typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView> VertexMapper;
+    typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView> ElementMapper;
+#else
+    typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGElementLayout> ElementMapper;
+    typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGVertexLayout> VertexMapper;
+#endif
+
+public:
+    typedef BaseOutputWriter::Scalar Scalar;
+    typedef BaseOutputWriter::Vector Vector;
+    typedef BaseOutputWriter::Tensor Tensor;
+    typedef BaseOutputWriter::ScalarBuffer ScalarBuffer;
+    typedef BaseOutputWriter::VectorBuffer VectorBuffer;
+    typedef BaseOutputWriter::TensorBuffer TensorBuffer;
+
+    typedef Dune::VTKWriter<GridView> VtkWriter;
+#if DUNE_VERSION_NEWER(DUNE_GRID, 2,5)
+    typedef std::shared_ptr< Dune::VTKFunction< GridView > > FunctionPtr;
+#else
+    typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
+#endif
+
+    VtkMultiWriter(bool asyncWriting,
+                   const GridView& gridView,
+                   const std::string& outputDir,
+                   const std::string& simName = "",
+                   std::string multiFileName = "")
+        : gridView_(gridView)
+#if DUNE_VERSION_NEWER(DUNE_GRID, 2,6)
+        , elementMapper_(gridView, Dune::mcmgElementLayout())
+        , vertexMapper_(gridView, Dune::mcmgVertexLayout())
+#else
+        , elementMapper_(gridView)
+        , vertexMapper_(gridView)
+#endif
+        , curWriter_(nullptr)
+        , curWriterNum_(0)
+        , taskletRunner_(/*numThreads=*/asyncWriting?1:0)
+    {
+        outputDir_ = outputDir;
+        if (outputDir == "")
+            outputDir_ = ".";
+
+        simName_ = (simName.empty()) ? "sim" : simName;
+        multiFileName_ = multiFileName;
+        if (multiFileName_.empty())
+            multiFileName_ = outputDir_+"/"+simName_+".pvd";
+
+        commRank_ = gridView.comm().rank();
+        commSize_ = gridView.comm().size();
+    }
+
+    ~VtkMultiWriter()
+    {
+        taskletRunner_.barrier();
+        releaseBuffers_();
+        finishMultiFile_();
+
+        if (commRank_ == 0)
+            multiFile_.close();
+    }
+
+    /*!
+     * \brief Returns the number of the current VTK file.
+     */
+    int curWriterNum() const
+    { return curWriterNum_; }
+
+    /*!
+     * \brief Updates the internal data structures after mesh
+     *        refinement.
+     *
+     * If the grid changes between two calls of beginWrite(), this
+     * method _must_ be called before the second beginWrite()!
+     */
+    void gridChanged()
+    {
+        elementMapper_.update();
+        vertexMapper_.update();
+    }
+
+    /*!
+     * \brief Called whenever a new time step must be written.
+     */
+    void beginWrite(double t)
+    {
+        if (!multiFile_.is_open()) {
+            startMultiFile_(multiFileName_);
+        }
+
+        // make sure that all previous output has been written and no other thread
+        // accesses the memory used as the target for the extracted quantities
+        taskletRunner_.barrier();
+        releaseBuffers_();
+
+        curTime_ = t;
+        curOutFileName_ = fileName_();
+
+        curWriter_ = new VtkWriter(gridView_, Dune::VTK::conforming);
+        ++curWriterNum_;
+    }
+
+    /*!
+     * \brief Allocate a managed buffer for a scalar field
+     *
+     * The buffer will be deleted automatically after the data has
+     * been written by to disk.
+     */
+    ScalarBuffer *allocateManagedScalarBuffer(size_t numEntities)
+    {
+        ScalarBuffer *buf = new ScalarBuffer(numEntities);
+        managedScalarBuffers_.push_back(buf);
+        return buf;
+    }
+
+    /*!
+     * \brief Allocate a managed buffer for a vector field
+     *
+     * The buffer will be deleted automatically after the data has
+     * been written by to disk.
+     */
+    VectorBuffer *allocateManagedVectorBuffer(size_t numOuter, size_t numInner)
+    {
+        VectorBuffer *buf = new VectorBuffer(numOuter);
+        for (size_t i = 0; i < numOuter; ++ i)
+            (*buf)[i].resize(numInner);
+
+        managedVectorBuffers_.push_back(buf);
+        return buf;
+    }
+
+    /*!
+     * \brief Add a finished vertex centered vector field to the
+     *        output.
+     *
+     * If the buffer is managed by the VtkMultiWriter, it must have
+     * been created using allocateManagedBuffer() and may not be used
+     * anywhere after calling this method. After the data is written
+     * to disk, it will be deleted automatically.
+     *
+     * If the buffer is not managed by the MultiWriter, the buffer
+     * must exist at least until the call to endWrite()
+     * finishes.
+     *
+     * In both cases, modifying the buffer between the call to this
+     * method and endWrite() results in _undefined behavior_.
+     */
+    void attachScalarVertexData(ScalarBuffer& buf, std::string name)
+    {
+        sanitizeScalarBuffer_(buf);
+
+        typedef Opm::VtkScalarFunction<GridView, VertexMapper> VtkFn;
+        FunctionPtr fnPtr(new VtkFn(name,
+                                    gridView_,
+                                    vertexMapper_,
+                                    buf,
+                                    /*codim=*/dim));
+        curWriter_->addVertexData(fnPtr);
+    }
+
+    /*!
+     * \brief Add a element centered quantity to the output.
+     *
+     * If the buffer is managed by the VtkMultiWriter, it must have
+     * been created using createField() and may not be used by
+     * anywhere after calling this method. After the data is written
+     * to disk, it will be deleted automatically.
+     *
+     * If the buffer is not managed by the MultiWriter, the buffer
+     * must exist at least until the call to endWrite()
+     * finishes.
+     *
+     * In both cases, modifying the buffer between the call to this
+     * method and endWrite() results in _undefined behaviour_.
+     */
+    void attachScalarElementData(ScalarBuffer& buf, std::string name)
+    {
+        sanitizeScalarBuffer_(buf);
+
+        typedef Opm::VtkScalarFunction<GridView, ElementMapper> VtkFn;
+        FunctionPtr fnPtr(new VtkFn(name,
+                                    gridView_,
+                                    elementMapper_,
+                                    buf,
+                                    /*codim=*/0));
+        curWriter_->addCellData(fnPtr);
+    }
+
+    /*!
+     * \brief Add a finished vertex centered vector field to the
+     *        output.
+     *
+     * If the buffer is managed by the VtkMultiWriter, it must have
+     * been created using allocateManagedBuffer() and may not be used
+     * anywhere after calling this method. After the data is written
+     * to disk, it will be deleted automatically.
+     *
+     * If the buffer is not managed by the MultiWriter, the buffer
+     * must exist at least until the call to endWrite()
+     * finishes.
+     *
+     * In both cases, modifying the buffer between the call to this
+     * method and endWrite() results in _undefined behavior_.
+     */
+    void attachVectorVertexData(VectorBuffer& buf, std::string name)
+    {
+        sanitizeVectorBuffer_(buf);
+
+        typedef Opm::VtkVectorFunction<GridView, VertexMapper> VtkFn;
+        FunctionPtr fnPtr(new VtkFn(name,
+                                    gridView_,
+                                    vertexMapper_,
+                                    buf,
+                                    /*codim=*/dim));
+        curWriter_->addVertexData(fnPtr);
+    }
+
+    /*!
+     * \brief Add a finished vertex-centered tensor field to the output.
+     */
+    void attachTensorVertexData(TensorBuffer& buf, std::string name)
+    {
+        typedef Opm::VtkTensorFunction<GridView, VertexMapper> VtkFn;
+
+        for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
+            std::ostringstream oss;
+            oss << name <<  "[" << colIdx << "]";
+
+            FunctionPtr fnPtr(new VtkFn(oss.str(),
+                                        gridView_,
+                                        vertexMapper_,
+                                        buf,
+                                        /*codim=*/dim,
+                                        colIdx));
+            curWriter_->addVertexData(fnPtr);
+        }
+    }
+
+    /*!
+     * \brief Add a element centered quantity to the output.
+     *
+     * If the buffer is managed by the VtkMultiWriter, it must have
+     * been created using createField() and may not be used by
+     * anywhere after calling this method. After the data is written
+     * to disk, it will be deleted automatically.
+     *
+     * If the buffer is not managed by the MultiWriter, the buffer
+     * must exist at least until the call to endWrite()
+     * finishes.
+     *
+     * In both cases, modifying the buffer between the call to this
+     * method and endWrite() results in _undefined behaviour_.
+     */
+    void attachVectorElementData(VectorBuffer& buf, std::string name)
+    {
+        sanitizeVectorBuffer_(buf);
+
+        typedef Opm::VtkVectorFunction<GridView, ElementMapper> VtkFn;
+        FunctionPtr fnPtr(new VtkFn(name,
+                                    gridView_,
+                                    elementMapper_,
+                                    buf,
+                                    /*codim=*/0));
+        curWriter_->addCellData(fnPtr);
+    }
+
+    /*!
+     * \brief Add a finished element-centered tensor field to the output.
+     */
+    void attachTensorElementData(TensorBuffer& buf, std::string name)
+    {
+        typedef Opm::VtkTensorFunction<GridView, ElementMapper> VtkFn;
+
+        for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
+            std::ostringstream oss;
+            oss << name <<  "[" << colIdx << "]";
+
+            FunctionPtr fnPtr(new VtkFn(oss.str(),
+                                        gridView_,
+                                        elementMapper_,
+                                        buf,
+                                        /*codim=*/0,
+                                        colIdx));
+            curWriter_->addCellData(fnPtr);
+        }
+    }
+
+    /*!
+     * \brief Finalizes the current writer.
+     *
+     * This means that everything will be written to disk, except if
+     * the onlyDiscard argument is true. In this case only all managed
+     * buffers are deleted, but no output is written.
+     */
+    void endWrite(bool onlyDiscard = false)
+    {
+        if (!onlyDiscard) {
+            auto tasklet = std::make_shared<WriteDataTasklet>(*this);
+            taskletRunner_.dispatch(tasklet);
+        }
+        else
+            --curWriterNum_;
+
+        // temporarily write the closing XML mumbo-jumbo to the mashup
+        // file so that the data set can be loaded even if the
+        // simulation is aborted (or not yet finished)
+        finishMultiFile_();
+    }
+
+    /*!
+     * \brief Write the multi-writer's state to a restart file.
+     */
+    template <class Restarter>
+    void serialize(Restarter& res)
+    {
+        res.serializeSectionBegin("VTKMultiWriter");
+        res.serializeStream() << curWriterNum_ << "\n";
+
+        if (commRank_ == 0) {
+            std::streamsize fileLen = 0;
+            std::streamoff filePos = 0;
+            if (multiFile_.is_open()) {
+                // write the meta file into the restart file
+                filePos = multiFile_.tellp();
+                multiFile_.seekp(0, std::ios::end);
+                fileLen = multiFile_.tellp();
+                multiFile_.seekp(filePos);
+            }
+
+            res.serializeStream() << fileLen << "  " << filePos << "\n";
+
+            if (fileLen > 0) {
+                std::ifstream multiFileIn(multiFileName_.c_str());
+                char *tmp = new char[fileLen];
+                multiFileIn.read(tmp, static_cast<long>(fileLen));
+                res.serializeStream().write(tmp, fileLen);
+                delete[] tmp;
+            }
+        }
+
+        res.serializeSectionEnd();
+    }
+
+    /*!
+     * \brief Read the multi-writer's state from a restart file.
+     */
+    template <class Restarter>
+    void deserialize(Restarter& res)
+    {
+        res.deserializeSectionBegin("VTKMultiWriter");
+        res.deserializeStream() >> curWriterNum_;
+
+        if (commRank_ == 0) {
+            std::string dummy;
+            std::getline(res.deserializeStream(), dummy);
+
+            // recreate the meta file from the restart file
+            std::streamoff filePos;
+            std::streamsize fileLen;
+            res.deserializeStream() >> fileLen >> filePos;
+            std::getline(res.deserializeStream(), dummy);
+            if (multiFile_.is_open())
+                multiFile_.close();
+
+            if (fileLen > 0) {
+                multiFile_.open(multiFileName_.c_str());
+
+                char *tmp = new char[fileLen];
+                res.deserializeStream().read(tmp, fileLen);
+                multiFile_.write(tmp, fileLen);
+                delete[] tmp;
+            }
+
+            multiFile_.seekp(filePos);
+        }
+        else {
+            std::string tmp;
+            std::getline(res.deserializeStream(), tmp);
+        }
+        res.deserializeSectionEnd();
+    }
+
+private:
+    std::string fileName_()
+    {
+        // use a new file name for each time step
+        std::ostringstream oss;
+        oss << simName_ << "-"
+            << std::setw(5) << std::setfill('0') << curWriterNum_;
+        return oss.str();
+    }
+
+    std::string fileSuffix_()
+    { return (GridView::dimension == 1) ? "vtp" : "vtu"; }
+
+    void startMultiFile_(const std::string& multiFileName)
+    {
+        // only the first process writes to the multi-file
+        if (commRank_ == 0) {
+            // generate one meta vtk-file holding the individual time steps
+            multiFile_.open(multiFileName.c_str());
+            multiFile_ << "<?xml version=\"1.0\"?>\n"
+                          "<VTKFile type=\"Collection\"\n"
+                          "         version=\"0.1\"\n"
+                          "         byte_order=\"LittleEndian\"\n"
+                          "         compressor=\"vtkZLibDataCompressor\">\n"
+                          " <Collection>\n";
+        }
+    }
+
+    void finishMultiFile_()
+    {
+        // only the first process writes to the multi-file
+        if (commRank_ == 0) {
+            // make sure that we always have a working meta file
+            std::ofstream::pos_type pos = multiFile_.tellp();
+            multiFile_ << " </Collection>\n"
+                          "</VTKFile>\n";
+            multiFile_.seekp(pos);
+            multiFile_.flush();
+        }
+    }
+
+    // make sure the field is well defined if running under valgrind
+    // and make sure that all values can be displayed by paraview
+    void sanitizeScalarBuffer_(ScalarBuffer& b OPM_UNUSED)
+    {
+        // nothing to do: this is done by VtkScalarFunction
+    }
+
+    void sanitizeVectorBuffer_(VectorBuffer& b OPM_UNUSED)
+    {
+        // nothing to do: this is done by VtkVectorFunction
+    }
+
+    // release the memory occupied by all buffer objects managed by the multi-writer
+    void releaseBuffers_()
+    {
+        // discard managed objects and the current VTK writer
+        delete curWriter_;
+        curWriter_ = nullptr;
+        while (managedScalarBuffers_.begin() != managedScalarBuffers_.end()) {
+            delete managedScalarBuffers_.front();
+            managedScalarBuffers_.pop_front();
+        }
+        while (managedVectorBuffers_.begin() != managedVectorBuffers_.end()) {
+            delete managedVectorBuffers_.front();
+            managedVectorBuffers_.pop_front();
+        }
+    }
+
+    const GridView gridView_;
+    ElementMapper elementMapper_;
+    VertexMapper vertexMapper_;
+
+    std::string outputDir_;
+    std::string simName_;
+    std::ofstream multiFile_;
+    std::string multiFileName_;
+
+    int commSize_; // number of processes in the communicator
+    int commRank_; // rank of the current process in the communicator
+
+    VtkWriter *curWriter_;
+    double curTime_;
+    std::string curOutFileName_;
+    int curWriterNum_;
+
+    std::list<ScalarBuffer *> managedScalarBuffers_;
+    std::list<VectorBuffer *> managedVectorBuffers_;
+
+    TaskletRunner taskletRunner_;
+};
+} // namespace Opm
+
+#endif

opm/models/io/vtkphasepresencemodule.hh

@@ -0,0 +1,141 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkPhasePresenceModule
+ */
+#ifndef EWOMS_VTK_PHASE_PRESENCE_MODULE_HH
+#define EWOMS_VTK_PHASE_PRESENCE_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/utils/propertysystem.hh>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK primary variables output
+NEW_TYPE_TAG(VtkPhasePresence);
+
+// create the property tags needed for the primary variables module
+NEW_PROP_TAG(VtkWritePhasePresence);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+SET_BOOL_PROP(VtkPhasePresence, VtkWritePhasePresence, false);
+
+END_PROPERTIES
+
+namespace Opm {
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the fluid composition
+ */
+template<class TypeTag>
+class VtkPhasePresenceModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+
+
+public:
+    VtkPhasePresenceModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the Vtk output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePhasePresence,
+                             "Include the phase presence pseudo primary "
+                             "variable in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (phasePresenceOutput_()) this->resizeScalarBuffer_(phasePresence_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quanties relevant
+     *        for an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            // calculate the phase presence
+            int phasePresence = elemCtx.primaryVars(i, /*timeIdx=*/0).phasePresence();
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+
+            if (phasePresenceOutput_())
+                phasePresence_[I] = phasePresence;
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (phasePresenceOutput_())
+            this->commitScalarBuffer_(baseWriter, "phase presence", phasePresence_);
+    }
+
+private:
+    static bool phasePresenceOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePhasePresence);
+        return val;
+    }
+
+    ScalarBuffer phasePresence_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkprimaryvarsmodule.hh

@@ -0,0 +1,185 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkPrimaryVarsModule
+ */
+#ifndef EWOMS_VTK_PRIMARY_VARS_MODULE_HH
+#define EWOMS_VTK_PRIMARY_VARS_MODULE_HH
+
+#include <opm/models/io/baseoutputmodule.hh>
+#include <opm/models/io/vtkmultiwriter.hh>
+
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/utils/propertysystem.hh>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK primary variables output
+NEW_PROP_TAG(EnableVtkOutput);
+NEW_TYPE_TAG(VtkPrimaryVars);
+
+// create the property tags needed for the primary variables module
+NEW_PROP_TAG(VtkWritePrimaryVars);
+NEW_PROP_TAG(VtkWriteProcessRank);
+NEW_PROP_TAG(VtkWriteDofIndex);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+SET_BOOL_PROP(VtkPrimaryVars, VtkWritePrimaryVars, false);
+SET_BOOL_PROP(VtkPrimaryVars, VtkWriteProcessRank, false);
+SET_BOOL_PROP(VtkPrimaryVars, VtkWriteDofIndex, false);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the fluid composition
+ */
+template<class TypeTag>
+class VtkPrimaryVarsModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+    typedef typename ParentType::EqBuffer EqBuffer;
+
+    enum { numEq = GET_PROP_VALUE(TypeTag, NumEq) };
+
+public:
+    VtkPrimaryVarsModule(const Simulator& simulator)
+        : ParentType(simulator)
+    { }
+
+    /*!
+     * \brief Register all run-time parameters for the Vtk output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePrimaryVars,
+                             "Include the primary variables into the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteProcessRank,
+                             "Include the MPI process rank into the VTK output files");
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteDofIndex,
+                             "Include the index of the degrees of freedom into the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (primaryVarsOutput_())
+            this->resizeEqBuffer_(primaryVars_);
+        if (processRankOutput_())
+            this->resizeScalarBuffer_(processRank_,
+                                      /*bufferType=*/ParentType::ElementBuffer);
+        if (dofIndexOutput_())
+            this->resizeScalarBuffer_(dofIndex_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant for
+     *        an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        const auto& elementMapper = elemCtx.model().elementMapper();
+        unsigned elemIdx = static_cast<unsigned>(elementMapper.index(elemCtx.element()));
+        if (processRankOutput_() && !processRank_.empty())
+            processRank_[elemIdx] = static_cast<unsigned>(this->simulator_.gridView().comm().rank());
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            const auto& priVars = elemCtx.primaryVars(i, /*timeIdx=*/0);
+
+            if (dofIndexOutput_())
+                dofIndex_[I] = I;
+
+            for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
+                if (primaryVarsOutput_() && !primaryVars_[eqIdx].empty())
+                    primaryVars_[eqIdx][I] = priVars[eqIdx];
+            }
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (primaryVarsOutput_())
+            this->commitPriVarsBuffer_(baseWriter, "PV_%s", primaryVars_);
+        if (processRankOutput_())
+            this->commitScalarBuffer_(baseWriter,
+                                      "process rank",
+                                      processRank_,
+                                      /*bufferType=*/ParentType::ElementBuffer);
+        if (dofIndexOutput_())
+            this->commitScalarBuffer_(baseWriter, "DOF index", dofIndex_);
+    }
+
+private:
+    static bool primaryVarsOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWritePrimaryVars);
+        return val;
+    }
+    static bool processRankOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteProcessRank);
+        return val;
+    }
+    static bool dofIndexOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteDofIndex);
+        return val;
+    }
+
+    EqBuffer primaryVars_;
+    ScalarBuffer processRank_;
+    ScalarBuffer dofIndex_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkscalarfunction.hh

@@ -0,0 +1,126 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkScalarFunction
+ */
+#ifndef VTK_SCALAR_FUNCTION_HH
+#define VTK_SCALAR_FUNCTION_HH
+
+#include <opm/models/io/baseoutputwriter.hh>
+
+#include <dune/grid/io/file/vtk/function.hh>
+#include <dune/istl/bvector.hh>
+#include <dune/common/fvector.hh>
+#include <dune/common/version.hh>
+
+#include <opm/material/common/Unused.hpp>
+
+#include <opm/material/common/Exceptions.hpp>
+
+#include <string>
+#include <limits>
+#include <vector>
+
+namespace Opm {
+
+/*!
+ * \brief Provides a vector-valued function using Dune::FieldVectors
+ *        as elements.
+ */
+template <class GridView, class Mapper>
+class VtkScalarFunction : public Dune::VTKFunction<GridView>
+{
+    enum { dim = GridView::dimension };
+    typedef typename GridView::ctype ctype;
+    typedef typename GridView::template Codim<0>::Entity Element;
+
+    typedef BaseOutputWriter::ScalarBuffer ScalarBuffer;
+
+public:
+    VtkScalarFunction(std::string name,
+                      const GridView& gridView,
+                      const Mapper& mapper,
+                      const ScalarBuffer& buf,
+                      unsigned codim)
+        : name_(name)
+        , gridView_(gridView)
+        , mapper_(mapper)
+        , buf_(buf)
+        , codim_(codim)
+    { assert(int(buf_.size()) == int(mapper_.size())); }
+
+    virtual std::string name() const
+    { return name_; }
+
+    virtual int ncomps() const
+    { return 1; }
+
+    virtual double evaluate(int mycomp OPM_UNUSED,
+                            const Element& e,
+                            const Dune::FieldVector<ctype, dim>& xi) const
+    {
+        unsigned idx;
+        if (codim_ == 0) {
+            // cells. map element to the index
+            idx = static_cast<unsigned>(mapper_.index(e));
+        }
+        else if (codim_ == dim) {
+            // find vertex which is closest to xi in local
+            // coordinates. This code is based on Dune::P1VTKFunction
+            double min = 1e100;
+            int imin = -1;
+            Dune::GeometryType gt = e.type();
+            int n = static_cast<int>(e.subEntities(dim));
+            for (int i = 0; i < n; ++i) {
+                Dune::FieldVector<ctype, dim> local =
+                    Dune::ReferenceElements<ctype, dim>::general(gt).position(i, dim);
+
+                local -= xi;
+                if (local.infinity_norm() < min) {
+                    min = local.infinity_norm();
+                    imin = static_cast<int>(i);
+                }
+            }
+
+            // map vertex to an index
+            idx = static_cast<unsigned>(mapper_.subIndex(e, imin, codim_));
+        }
+        else
+            throw std::logic_error("Only element and vertex based vector fields are"
+                                   " supported so far.");
+
+        return static_cast<double>(static_cast<float>(buf_[idx]));
+    }
+
+private:
+    const std::string name_;
+    const GridView gridView_;
+    const Mapper& mapper_;
+    const ScalarBuffer& buf_;
+    unsigned codim_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtktemperaturemodule.hh

@@ -0,0 +1,147 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkTemperatureModule
+ */
+#ifndef EWOMS_VTK_TEMPERATURE_MODULE_HH
+#define EWOMS_VTK_TEMPERATURE_MODULE_HH
+
+#include "vtkmultiwriter.hh"
+#include "baseoutputmodule.hh"
+
+#include <opm/models/utils/parametersystem.hh>
+#include <opm/models/utils/propertysystem.hh>
+
+#include <opm/material/common/MathToolbox.hpp>
+
+BEGIN_PROPERTIES
+
+// create new type tag for the VTK temperature output
+NEW_TYPE_TAG(VtkTemperature);
+
+// create the property tags needed for the temperature module
+NEW_PROP_TAG(VtkWriteTemperature);
+NEW_PROP_TAG(VtkOutputFormat);
+NEW_PROP_TAG(EnableVtkOutput);
+
+// set default values for what quantities to output
+SET_BOOL_PROP(VtkTemperature, VtkWriteTemperature, true);
+
+END_PROPERTIES
+
+namespace Opm {
+
+/*!
+ * \ingroup Vtk
+ *
+ * \brief VTK output module for the temperature in which assume
+ *        thermal equilibrium
+ */
+template<class TypeTag>
+class VtkTemperatureModule : public BaseOutputModule<TypeTag>
+{
+    typedef BaseOutputModule<TypeTag> ParentType;
+
+    typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
+
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
+
+    typedef typename ParentType::ScalarBuffer ScalarBuffer;
+
+    static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
+    typedef Opm::VtkMultiWriter<GridView, vtkFormat> VtkMultiWriter;
+
+public:
+    VtkTemperatureModule(const Simulator& simulator)
+        : ParentType(simulator)
+    {}
+
+    /*!
+     * \brief Register all run-time parameters for the Vtk output module.
+     */
+    static void registerParameters()
+    {
+        EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteTemperature,
+                             "Include the temperature in the VTK output files");
+    }
+
+    /*!
+     * \brief Allocate memory for the scalar fields we would like to
+     *        write to the VTK file.
+     */
+    void allocBuffers()
+    {
+        if (temperatureOutput_()) this->resizeScalarBuffer_(temperature_);
+    }
+
+    /*!
+     * \brief Modify the internal buffers according to the intensive quantities relevant
+     *        for an element
+     */
+    void processElement(const ElementContext& elemCtx)
+    {
+        typedef Opm::MathToolbox<Evaluation> Toolbox;
+
+        if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
+            return;
+
+        for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
+            unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
+            const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
+            const auto& fs = intQuants.fluidState();
+
+            if (temperatureOutput_())
+                temperature_[I] = Toolbox::value(fs.temperature(/*phaseIdx=*/0));
+        }
+    }
+
+    /*!
+     * \brief Add all buffers to the VTK output writer.
+     */
+    void commitBuffers(BaseOutputWriter& baseWriter)
+    {
+        VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
+        if (!vtkWriter) {
+            return;
+        }
+
+        if (temperatureOutput_())
+            this->commitScalarBuffer_(baseWriter, "temperature", temperature_);
+    }
+
+private:
+    static bool temperatureOutput_()
+    {
+        static bool val = EWOMS_GET_PARAM(TypeTag, bool, VtkWriteTemperature);
+        return val;
+    }
+
+    ScalarBuffer temperature_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtktensorfunction.hh

@@ -0,0 +1,127 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkTensorFunction
+ */
+#ifndef VTK_TENSOR_FUNCTION_HH
+#define VTK_TENSOR_FUNCTION_HH
+
+#include <opm/models/io/baseoutputwriter.hh>
+
+#include <dune/grid/io/file/vtk/function.hh>
+#include <dune/common/fvector.hh>
+#include <dune/common/version.hh>
+
+#include <opm/material/common/Exceptions.hpp>
+
+#include <string>
+#include <limits>
+#include <vector>
+
+namespace Opm {
+
+/*!
+ * \brief Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
+ */
+template <class GridView, class Mapper>
+class VtkTensorFunction : public Dune::VTKFunction<GridView>
+{
+    enum { dim = GridView::dimension };
+    typedef typename GridView::ctype ctype;
+    typedef typename GridView::template Codim<0>::Entity Element;
+
+    typedef BaseOutputWriter::TensorBuffer TensorBuffer;
+
+public:
+    VtkTensorFunction(std::string name,
+                      const GridView& gridView,
+                      const Mapper& mapper,
+                      const TensorBuffer& buf,
+                      unsigned codim,
+                      unsigned matrixColumnIdx)
+        : name_(name)
+        , gridView_(gridView)
+        , mapper_(mapper)
+        , buf_(buf)
+        , codim_(codim)
+        , matrixColumnIdx_(matrixColumnIdx)
+    { assert(int(buf_.size()) == int(mapper_.size())); }
+
+    virtual std::string name() const
+    { return name_; }
+
+    virtual int ncomps() const
+    { return static_cast<int>(buf_[0].M()); }
+
+    virtual double evaluate(int mycomp,
+                            const Element& e,
+                            const Dune::FieldVector<ctype, dim>& xi) const
+    {
+        size_t idx;
+        if (codim_ == 0) {
+            // cells. map element to the index
+            idx = static_cast<size_t>(mapper_.index(e));
+        }
+        else if (codim_ == dim) {
+            // find vertex which is closest to xi in local
+            // coordinates. This code is based on Dune::P1VTKFunction
+            double min = 1e100;
+            int imin = -1;
+            Dune::GeometryType gt = e.type();
+            int n = static_cast<int>(e.subEntities(dim));
+            for (int i = 0; i < n; ++i) {
+                Dune::FieldVector<ctype, dim> local =
+                    Dune::ReferenceElements<ctype, dim>::general(gt).position(i, dim);
+
+                local -= xi;
+                if (local.infinity_norm() < min) {
+                    min = local.infinity_norm();
+                    imin = i;
+                }
+            }
+
+            // map vertex to an index
+            idx = static_cast<size_t>(mapper_.subIndex(e, imin, codim_));
+        }
+        else
+            throw std::logic_error("Only element and vertex based tensor fields are supported so far.");
+
+        unsigned i = static_cast<unsigned>(mycomp);
+        unsigned j = static_cast<unsigned>(matrixColumnIdx_);
+
+        return static_cast<double>(static_cast<float>(buf_[idx][i][j]));
+    }
+
+private:
+    const std::string name_;
+    const GridView gridView_;
+    const Mapper& mapper_;
+    const TensorBuffer& buf_;
+    unsigned codim_;
+    unsigned matrixColumnIdx_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/io/vtkvectorfunction.hh

@@ -0,0 +1,124 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM 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 General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ * \copydoc Opm::VtkVectorFunction
+ */
+#ifndef VTK_VECTOR_FUNCTION_HH
+#define VTK_VECTOR_FUNCTION_HH
+
+#include <opm/models/io/baseoutputwriter.hh>
+
+#include <dune/grid/io/file/vtk/function.hh>
+#include <dune/istl/bvector.hh>
+#include <dune/common/fvector.hh>
+#include <dune/common/version.hh>
+
+#include <opm/material/common/Exceptions.hpp>
+
+#include <string>
+#include <limits>
+#include <vector>
+
+namespace Opm {
+
+/*!
+ * \brief Provides a vector-valued function using Dune::FieldVectors
+ *        as elements.
+ */
+template <class GridView, class Mapper>
+class VtkVectorFunction : public Dune::VTKFunction<GridView>
+{
+    enum { dim = GridView::dimension };
+    typedef typename GridView::ctype ctype;
+    typedef typename GridView::template Codim<0>::Entity Element;
+
+    typedef BaseOutputWriter::VectorBuffer VectorBuffer;
+
+public:
+    VtkVectorFunction(std::string name,
+                      const GridView& gridView,
+                      const Mapper& mapper,
+                      const VectorBuffer& buf,
+                      unsigned codim)
+        : name_(name)
+        , gridView_(gridView)
+        , mapper_(mapper)
+        , buf_(buf)
+        , codim_(codim)
+    { assert(int(buf_.size()) == int(mapper_.size())); }
+
+    virtual std::string name() const
+    { return name_; }
+
+    virtual int ncomps() const
+    { return static_cast<int>(buf_[0].size()); }
+
+    virtual double evaluate(int mycomp,
+                            const Element& e,
+                            const Dune::FieldVector<ctype, dim>& xi) const
+    {
+        unsigned idx;
+        if (codim_ == 0) {
+            // cells. map element to the index
+            idx = static_cast<unsigned>(mapper_.index(e));
+        }
+        else if (codim_ == dim) {
+            // find vertex which is closest to xi in local
+            // coordinates. This code is based on Dune::P1VTKFunction
+            double min = 1e100;
+            int imin = -1;
+            Dune::GeometryType gt = e.type();
+            int n = static_cast<int>(e.subEntities(dim));
+            for (int i = 0; i < n; ++i) {
+                Dune::FieldVector<ctype, dim> local =
+                    Dune::ReferenceElements<ctype, dim>::general(gt).position(i, dim);
+
+                local -= xi;
+                if (local.infinity_norm() < min) {
+                    min = local.infinity_norm();
+                    imin = i;
+                }
+            }
+
+            // map vertex to an index
+            idx = static_cast<unsigned>(mapper_.subIndex(e, imin, codim_));
+        }
+        else
+            throw std::logic_error("Only element and vertex based vector fields are "
+                                   "supported so far.");
+
+        return static_cast<double>(static_cast<float>(buf_[idx][static_cast<unsigned>(mycomp)]));
+    }
+
+private:
+    const std::string name_;
+    const GridView gridView_;
+    const Mapper& mapper_;
+    const VectorBuffer& buf_;
+    unsigned codim_;
+};
+
+} // namespace Opm
+
+#endif

opm/models/utils/basicproperties.hh

@@ -32,7 +32,7 @@
 
 #include <opm/models/utils/propertysystem.hh>
 #include <opm/models/utils/parametersystem.hh>
-#include <ewoms/io/dgfvanguard.hh>
+#include <opm/models/io/dgfvanguard.hh>
 
 #if HAVE_DUNE_FEM
 #include <dune/fem/gridpart/adaptiveleafgridpart.hh>

opm/models/utils/simulator.hh

@@ -28,7 +28,7 @@
 #ifndef EWOMS_SIMULATOR_HH
 #define EWOMS_SIMULATOR_HH
 
-#include <ewoms/io/restart.hh>
+#include <opm/models/io/restart.hh>
 #include <opm/models/utils/parametersystem.hh>
 
 #include <opm/models/utils/propertysystem.hh>