[test][multidomain][facet] Add three domain (3d-2d-1d) fracture test

test/multidomain/facet/1p_1p/CMakeLists.txt

+add_subdirectory(threedomain)

test/multidomain/facet/1p_1p/threedomain/CMakeLists.txt

+dune_symlink_to_source_files(FILES "grids" "test_facetcoupling_1p1p_threedomain.input")
+
+dune_add_test(NAME test_facetcoupling_tpfa_1p1p_threedomain
+              SOURCES test_facetcoupling_tpfa_1p1p_threedomain.cc
+              CMAKE_GUARD dune-foamgrid_FOUND
+              CMAKE_GUARD dune-alugrid_FOUND
+              COMPILE_DEFINITIONS BULKTYPETAG=TissueCCTypeTag LOWDIMTYPETAG=BloodFlowCCTypeTag COUPLINGMODE=EmbeddedCouplingMode::average
+              COMMAND ${CMAKE_SOURCE_DIR}/bin/testing/runtest.py
+              CMD_ARGS  --script fuzzy
+                        --files ${CMAKE_SOURCE_DIR}/test/references/facetcoupling_1p1p_threedomain_bulk.vtu
+                                ${CMAKE_CURRENT_BINARY_DIR}/test_1p_1p_bulk-00001.vtu
+                                ${CMAKE_SOURCE_DIR}/test/references/facetcoupling_1p1p_threedomain_facet.vtu
+                                ${CMAKE_CURRENT_BINARY_DIR}/test_1p_1p_facet-00001.vtu
+                                ${CMAKE_SOURCE_DIR}/test/references/facetcoupling_1p1p_threedomain_edge.vtp
+                                ${CMAKE_CURRENT_BINARY_DIR}/test_1p_1p_edge-00001.vtp
+                       --command "${CMAKE_CURRENT_BINARY_DIR}/test_facetcoupling_tpfa_1p1p_threedomain test_facetcoupling_1p1p_threedomain.input")
+
+set(CMAKE_BUILD_TYPE Release)
+
+#install sources
+install(FILES
+test_facetcoupling_tpfa_1p1p_threedomain.cc
+bulkProblem.hh
+facetproblem.hh
+edgeproblem.hh
+spatialparams.hh
+DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/test/multidomain/facet/1p_1p/threedomain)

test/multidomain/facet/1p_1p/threedomain/bulkproblem.hh

+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program 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.                                     *
+ *                                                                           *
+ *   This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ * \ingroup MixedDimension
+ * \ingroup MixedDimensionFacet
+ * \brief The problem for the bulk domain in the single-phase facet coupling test
+ */
+#ifndef DUMUX_TEST_FACETCOUPLING_THREEDOMAIN_ONEP_BULKPROBLEM_HH
+#define DUMUX_TEST_FACETCOUPLING_THREEDOMAIN_ONEP_BULKPROBLEM_HH
+
+#include <dune/alugrid/grid.hh>
+
+#include <dumux/material/components/constant.hh>
+#include <dumux/material/fluidsystems/1pliquid.hh>
+
+#include <dumux/multidomain/facet/cellcentered/tpfa/properties.hh>
+#include <dumux/porousmediumflow/problem.hh>
+#include <dumux/porousmediumflow/1p/model.hh>
+
+#include "spatialparams.hh"
+
+namespace Dumux {
+// forward declarations
+template<class TypeTag> class OnePBulkProblem;
+
+namespace Properties {
+// create the type tag nodes
+NEW_TYPE_TAG(OnePBulk, INHERITS_FROM(OneP));
+NEW_TYPE_TAG(OnePBulkTpfa, INHERITS_FROM(OnePBulk, CCTpfaFacetCouplingModel));
+
+// Set the grid type
+SET_TYPE_PROP(OnePBulk, Grid, Dune::ALUGrid<3, 3, Dune::simplex, Dune::nonconforming>);
+// Set the problem type
+SET_TYPE_PROP(OnePBulk, Problem, OnePBulkProblem<TypeTag>);
+// set the spatial params
+SET_TYPE_PROP(OnePBulk, SpatialParams, OnePSpatialParams<TypeTag>);
+
+// the fluid system
+SET_PROP(OnePBulk, FluidSystem)
+{
+private:
+    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+public:
+    using type = FluidSystems::OnePLiquid< Scalar, Components::Constant<1, Scalar> >;
+};
+
+} // end namespace Properties
+
+/*!
+ * \ingroup OnePTests
+ * \brief Test problem for the incompressible one-phase model
+ *        with coupling across the bulk grid facets
+ */
+template<class TypeTag>
+class OnePBulkProblem : public PorousMediumFlowProblem<TypeTag>
+{
+    using ParentType = PorousMediumFlowProblem<TypeTag>;
+
+    using NumEqVector = typename GET_PROP_TYPE(TypeTag, NumEqVector);
+    using GridVariables = typename GET_PROP_TYPE(TypeTag, GridVariables);
+    using PrimaryVariables = typename GridVariables::PrimaryVariables;
+    using Scalar = typename GridVariables::Scalar;
+
+    using FVGridGeometry = typename GridVariables::GridGeometry;
+    using SubControlVolumeFace = typename FVGridGeometry::SubControlVolumeFace;
+    using GridView = typename FVGridGeometry::GridView;
+    using Element = typename GridView::template Codim<0>::Entity;
+    using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
+
+    using BoundaryTypes = typename GET_PROP_TYPE(TypeTag, BoundaryTypes);
+    using CouplingManager = typename GET_PROP_TYPE(TypeTag, CouplingManager);
+
+public:
+    //! The constructor
+    OnePBulkProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
+                    std::shared_ptr<typename ParentType::SpatialParams> spatialParams,
+                    const std::string& paramGroup = "")
+    : ParentType(fvGridGeometry, spatialParams, paramGroup)
+    {}
+
+    //! Specifies the kind of boundary condition on a given boundary position.
+    BoundaryTypes boundaryTypesAtPos(const GlobalPosition& globalPos) const
+    {
+        BoundaryTypes values;
+        values.setAllNeumann();
+        if (globalPos[1] < this->fvGridGeometry().bBoxMin()[1] + 1e-6 || globalPos[1] > this->fvGridGeometry().bBoxMax()[1] - 1e-6)
+            values.setAllDirichlet();
+        return values;
+    }
+
+    /*!
+     * \brief Specifies which kind of interior boundary condition should be
+     *        used for which equation on a given sub-control volume face
+     *        that couples to a facet element.
+     *
+     * \param element The finite element the scvf is embedded in
+     * \param scvf The sub-control volume face
+     */
+    BoundaryTypes interiorBoundaryTypes(const Element& element, const SubControlVolumeFace& scvf) const
+    {
+        BoundaryTypes values;
+        values.setAllNeumann();
+        return values;
+    }
+
+    //! Evaluate the Dirichlet boundary conditions at a given position
+    PrimaryVariables dirichletAtPos(const GlobalPosition& globalPos) const
+    {
+        const auto y = globalPos[1];
+        const auto yMin = this->fvGridGeometry().bBoxMin()[1];
+        const auto yMax = this->fvGridGeometry().bBoxMax()[1];
+
+        return PrimaryVariables( {2.0 - (y-yMin)/(yMax-yMin)} );
+    }
+
+    //! Evaluat the initial conditions
+    PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
+    { return PrimaryVariables(1.0); }
+
+    //! Return the temperature in \f$\mathrm{[K]}\f$ in the domain
+    Scalar temperature() const
+    { return 283.15; /*10°C*/ }
+
+    //! Return const reference to the coupling manager.
+    const CouplingManager& couplingManager() const
+    { return *couplingManagerPtr_; }
+
+    //! sets the pointer to the coupling manager.
+    void setCouplingManager(std::shared_ptr<CouplingManager> cm)
+    { couplingManagerPtr_ = cm; }
+
+private:
+    std::shared_ptr<CouplingManager> couplingManagerPtr_;
+};
+
+} // end namespace Dumux
+
+#endif

test/multidomain/facet/1p_1p/threedomain/edgeproblem.hh

+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program 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.                                     *
+ *                                                                           *
+ *   This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ * \ingroup MixedDimension
+ * \ingroup MixedDimensionFacet
+ * \brief The problem for the (d-2)-dimensional edge domain in the single-phase
+ *        facet coupling test involving three domains.
+ */
+#ifndef DUMUX_TEST_FACETCOUPLING_THREEDOMAIN_ONEP_EDGEPROBLEM_HH
+#define DUMUX_TEST_FACETCOUPLING_THREEDOMAIN_ONEP_EDGEPROBLEM_HH
+
+#include <dune/foamgrid/foamgrid.hh>
+
+#include <dumux/material/components/constant.hh>
+#include <dumux/material/fluidsystems/1pliquid.hh>
+
+#include <dumux/discretization/cellcentered/tpfa/properties.hh>
+
+#include <dumux/porousmediumflow/problem.hh>
+#include <dumux/porousmediumflow/1p/model.hh>
+
+#include "spatialparams.hh"
+
+namespace Dumux {
+// forward declarations
+template<class TypeTag> class OnePEdgeProblem;
+
+namespace Properties {
+// create the type tag nodes
+NEW_TYPE_TAG(OnePEdge, INHERITS_FROM(OneP));
+NEW_TYPE_TAG(OnePEdgeTpfa, INHERITS_FROM(CCTpfaModel, OnePEdge));
+
+// Set the grid type
+SET_TYPE_PROP(OnePEdge, Grid, Dune::FoamGrid<1, 3>);
+// Set the problem type
+SET_TYPE_PROP(OnePEdge, Problem, OnePEdgeProblem<TypeTag>);
+// set the spatial params
+SET_TYPE_PROP(OnePEdge, SpatialParams, OnePSpatialParams<TypeTag>);
+
+// the fluid system
+SET_PROP(OnePEdge, FluidSystem)
+{
+private:
+    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+public:
+    using type = FluidSystems::OnePLiquid< Scalar, Components::Constant<1, Scalar> >;
+};
+
+} // end namespace Properties
+/*!
+ * \ingroup OnePTests
+ * \brief The (d-2)-dimensional test problem for the incompressible
+ *        one-phase model with coupling across the bulk grid facets
+ */
+template<class TypeTag>
+class OnePEdgeProblem : public PorousMediumFlowProblem<TypeTag>
+{
+    using ParentType = PorousMediumFlowProblem<TypeTag>;
+
+    using NumEqVector = typename GET_PROP_TYPE(TypeTag, NumEqVector);
+    using GridVariables = typename GET_PROP_TYPE(TypeTag, GridVariables);
+    using ElementVolumeVariables = typename GridVariables::GridVolumeVariables::LocalView;
+    using PrimaryVariables = typename GridVariables::PrimaryVariables;
+    using Scalar = typename GridVariables::Scalar;
+
+    using FVGridGeometry = typename GridVariables::GridGeometry;
+    using FVElementGeometry = typename FVGridGeometry::LocalView;
+    using SubControlVolume = typename FVGridGeometry::SubControlVolume;
+    using SubControlVolumeFace = typename FVGridGeometry::SubControlVolumeFace;
+    using GridView = typename FVGridGeometry::GridView;
+    using Element = typename GridView::template Codim<0>::Entity;
+    using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
+
+    using BoundaryTypes = typename GET_PROP_TYPE(TypeTag, BoundaryTypes);
+    using CouplingManager = typename GET_PROP_TYPE(TypeTag, CouplingManager);
+
+public:
+    //! The constructor
+    OnePEdgeProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
+                    std::shared_ptr<typename ParentType::SpatialParams> spatialParams,
+                    const std::string& paramGroup = "")
+    : ParentType(fvGridGeometry, spatialParams, paramGroup)
+    {
+        const auto a = getParam<Scalar>("Extrusion.Aperture");
+        exFactor_ = a*a;
+    }
+
+    //!Specifies the type of boundary condition on a boundary position
+    BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
+    {
+        BoundaryTypes values;
+        values.setAllNeumann();
+        return values;
+    }
+
+    /*!
+     * \brief Evaluate the source term for all phases within a given
+     *        sub-control-volume.
+     */
+    NumEqVector source(const Element& element,
+                       const FVElementGeometry& fvGeometry,
+                       const ElementVolumeVariables& elemVolVars,
+                       const SubControlVolume& scv) const
+    {
+        // forward to solution independent, fully-implicit specific interface
+        auto source = couplingManagerPtr_->evalSourcesFromBulk(element, fvGeometry, elemVolVars, scv);
+        source /= scv.volume()*elemVolVars[scv].extrusionFactor();
+        return source;
+    }
+
+    //! Set the aperture squared as extrusion factor.
+    Scalar extrusionFactorAtPos(const GlobalPosition& globalPos) const
+    { return exFactor_; }
+
+    //! Evaluate the initial conditions
+    PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
+    { return PrimaryVariables(1.0); }
+
+    //! Returns the temperature \f$\mathrm{[K]}\f$ for an isothermal problem.
+    Scalar temperature() const
+    { return 283.15; /*10°*/ }
+
+    //! Return const reference to the coupling manager.
+    const CouplingManager& couplingManager() const
+    { return *couplingManagerPtr_; }
+
+    //! sets the pointer to the coupling manager.
+    void setCouplingManager(std::shared_ptr<CouplingManager> cm)
+    { couplingManagerPtr_ = cm; }
+
+private:
+    Scalar exFactor_;
+    std::shared_ptr<CouplingManager> couplingManagerPtr_;
+};
+
+} // end namespace Dumux
+
+#endif

test/multidomain/facet/1p_1p/threedomain/facetproblem.hh

+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program 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.                                     *
+ *                                                                           *
+ *   This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ * \ingroup MixedDimension
+ * \ingroup MixedDimensionFacet
+ * \brief The problem for the (d-1)-dimensional facet domain in the single-phase
+ *        facet coupling test involving three domains.
+ */
+#ifndef DUMUX_TEST_FACETCOUPLING_THREEDOMAIN_ONEP_FACETPROBLEM_HH
+#define DUMUX_TEST_FACETCOUPLING_THREEDOMAIN_ONEP_FACETPROBLEM_HH
+
+#include <dune/foamgrid/foamgrid.hh>
+
+#include <dumux/material/components/constant.hh>
+#include <dumux/material/fluidsystems/1pliquid.hh>
+
+#include <dumux/multidomain/facet/cellcentered/tpfa/properties.hh>
+#include <dumux/porousmediumflow/problem.hh>
+#include <dumux/porousmediumflow/1p/model.hh>
+
+#include "spatialparams.hh"
+
+namespace Dumux {
+// forward declarations
+template<class TypeTag> class OnePFacetProblem;
+
+namespace Properties {
+// create the type tag nodes
+NEW_TYPE_TAG(OnePFacet, INHERITS_FROM(OneP));
+NEW_TYPE_TAG(OnePFacetTpfa, INHERITS_FROM(OnePFacet, CCTpfaFacetCouplingModel));
+
+// Set the grid type
+SET_TYPE_PROP(OnePFacet, Grid, Dune::FoamGrid<2, 3>);
+// Set the problem type
+SET_TYPE_PROP(OnePFacet, Problem, OnePFacetProblem<TypeTag>);
+// set the spatial params
+SET_TYPE_PROP(OnePFacet, SpatialParams, OnePSpatialParams<TypeTag>);
+
+// the fluid system
+SET_PROP(OnePFacet, FluidSystem)
+{
+private:
+    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+public:
+    using type = FluidSystems::OnePLiquid< Scalar, Components::Constant<1, Scalar> >;
+};
+
+} // end namespace Properties
+
+/*!
+ * \ingroup OnePTests
+ * \brief The (d-1)-dimensional test problem for the incompressible
+ *        one-phase model with coupling across the bulk grid facets
+ */
+template<class TypeTag>
+class OnePFacetProblem : public PorousMediumFlowProblem<TypeTag>
+{
+    using ParentType = PorousMediumFlowProblem<TypeTag>;
+
+    using NumEqVector = typename GET_PROP_TYPE(TypeTag, NumEqVector);
+    using GridVariables = typename GET_PROP_TYPE(TypeTag, GridVariables);
+    using ElementVolumeVariables = typename GridVariables::GridVolumeVariables::LocalView;
+    using PrimaryVariables = typename GridVariables::PrimaryVariables;
+    using Scalar = typename GridVariables::Scalar;
+
+    using FVGridGeometry = typename GridVariables::GridGeometry;
+    using FVElementGeometry = typename FVGridGeometry::LocalView;
+    using SubControlVolume = typename FVGridGeometry::SubControlVolume;
+    using SubControlVolumeFace = typename FVGridGeometry::SubControlVolumeFace;
+    using GridView = typename FVGridGeometry::GridView;
+    using Element = typename GridView::template Codim<0>::Entity;
+    using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
+
+    using BoundaryTypes = typename GET_PROP_TYPE(TypeTag, BoundaryTypes);
+    using CouplingManager = typename GET_PROP_TYPE(TypeTag, CouplingManager);
+
+public:
+    //! The constructor
+    OnePFacetProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
+                     std::shared_ptr<typename ParentType::SpatialParams> spatialParams,
+                     const std::string& paramGroup = "")
+    : ParentType(fvGridGeometry, spatialParams, paramGroup)
+    , aperture_(getParam<Scalar>("Extrusion.Aperture"))
+    {}
+
+    //! Specifies the kind of boundary condition at a boundary position
+    BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
+    {
+        BoundaryTypes values;
+        values.setAllNeumann();
+        return values;
+    }
+
+    /*!
+     * \brief Specifies which kind of interior boundary condition should be
+     *        used for which equation on a given sub-control volume face
+     *        that couples to a facet element.
+     *
+     * \param element The finite element the scvf is embedded in
+     * \param scvf The sub-control volume face
+     */
+    BoundaryTypes interiorBoundaryTypes(const Element& element, const SubControlVolumeFace& scvf) const
+    {
+        BoundaryTypes values;
+        values.setAllNeumann();
+        return values;
+    }
+
+    /*!
+     * \brief Evaluate the source term for all phases within a given
+     *        sub-control-volume.
+     */
+    NumEqVector source(const Element& element,
+                       const FVElementGeometry& fvGeometry,
+                       const ElementVolumeVariables& elemVolVars,
+                       const SubControlVolume& scv) const
+    {
+        // forward to solution independent, fully-implicit specific interface
+        auto source = couplingManagerPtr_->evalSourcesFromBulk(element, fvGeometry, elemVolVars, scv);
+        source /= scv.volume()*elemVolVars[scv].extrusionFactor();
+        return source;
+    }
+
+    //! Set the aperture as extrusion factor.
+    Scalar extrusionFactorAtPos(const GlobalPosition& globalPos) const
+    { return aperture_; }
+
+    //! Evaluate the initial conditions
+    PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
+    { return PrimaryVariables(1.0); }
+
+    //! Returns the temperature \f$\mathrm{[K]}\f$ for an isothermal problem.
+    Scalar temperature() const
+    { return 283.15; /*10°*/ }
+
+    //! Return const reference to the coupling manager.
+    const CouplingManager& couplingManager() const
+    { return *couplingManagerPtr_; }
+
+    //! sets the pointer to the coupling manager.
+    void setCouplingManager(std::shared_ptr<CouplingManager> cm)
+    { couplingManagerPtr_ = cm; }
+
+private:
+    Scalar aperture_;
+    std::shared_ptr<CouplingManager> couplingManagerPtr_;
+};
+
+} // end namespace Dumux
+
+#endif

test/multidomain/facet/1p_1p/threedomain/spatialparams.hh

+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program 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.                                     *
+ *                                                                           *
+ *   This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ * \ingroup OnePTests
+ * \brief The spatial params the single-phase facet coupling test
+ */
+#ifndef DUMUX_TEST_TPFAFACETCOUPLING_THREEDOMAIN_ONEP_SPATIALPARAMS_HH
+#define DUMUX_TEST_TPFAFACETCOUPLING_THREEDOMAIN_ONEP_SPATIALPARAMS_HH
+
+#include <dumux/material/spatialparams/fv1p.hh>
+
+namespace Dumux {
+
+/*!
+ * \ingroup OnePTests
+ * \brief The spatial params the single-phase facet coupling test
+ */
+template<class TypeTag>
+class OnePSpatialParams : public FVSpatialParamsOneP< typename GET_PROP_TYPE(TypeTag, FVGridGeometry),
+                                                      typename GET_PROP_TYPE(TypeTag, Scalar),
+                                                      OnePSpatialParams<TypeTag> >
+{
+    using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+    using FVGridGeometry = typename GET_PROP_TYPE(TypeTag, FVGridGeometry);
+    using GridView = typename FVGridGeometry::GridView;
+    using Element = typename GridView::template Codim<0>::Entity;
+    using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
+
+    using ThisType = OnePSpatialParams<TypeTag>;
+    using ParentType = FVSpatialParamsOneP<FVGridGeometry, Scalar, ThisType>;
+
+public:
+    //! export the type used for permeability
+    using PermeabilityType = Scalar;
+
+    //! the constructor
+    OnePSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry, const std::string& paramGroup = "")
+    : ParentType(fvGridGeometry)
+    {
+        permeability_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Permeability");
+    }
+
+    //! Function for defining the (intrinsic) permeability \f$[m^2]\f$.
+    PermeabilityType permeabilityAtPos(const GlobalPosition& globalPos) const
+    { return permeability_; }
+
+    //! Return the porosity
+    Scalar porosityAtPos(const GlobalPosition& globalPos) const
+    { return 1.0; }
+
+private:
+    Scalar permeability_;
+};
+
+} // end namespace Dumux
+
+#endif

test/multidomain/facet/1p_1p/threedomain/test_facetcoupling_1p1p_threedomain.input

+[Problem]
+EnableGravity = false
+
+[Grid]
+File = ./grids/grid.msh
+
+[Bulk]
+Problem.Name = test_1p_1p_bulk
+SpatialParams.Permeability = 1
+FacetCoupling.Xi = 1.0
+
+[Extrusion]
+Aperture = 1e-2
+
+[Facet]
+SpatialParams.Permeability = 1e3
+Problem.Name = test_1p_1p_facet
+Problem.Aperture = 1e-2
+FacetCoupling.Xi = 1.0
+
+[Edge]
+SpatialParams.Permeability = 1e5
+Problem.Name = test_1p_1p_edge
+Problem.Aperture = 1e-3

test/multidomain/facet/1p_1p/threedomain/test_facetcoupling_tpfa_1p1p_threedomain.cc

+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*****************************************************************************
+ *   See the file COPYING for full copying permissions.                      *
+ *                                                                           *
+ *   This program 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.                                     *
+ *                                                                           *
+ *   This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.   *
+ *****************************************************************************/
+/*!
+ * \file
+ *
+ * \brief test for the one-phase facet coupling model
+ */
+#include <config.h>
+
+#include <iostream>
+
+#include <dune/common/parallel/mpihelper.hh>
+
+#include <dumux/common/properties.hh>
+#include <dumux/common/parameters.hh>
+#include <dumux/common/dumuxmessage.hh>
+#include <dumux/common/defaultusagemessage.hh>
+
+#include "bulkproblem.hh"
+#include "facetproblem.hh"
+#include "edgeproblem.hh"
+
+#include <dumux/assembly/diffmethod.hh>
+
+#include <dumux/linear/seqsolverbackend.hh>
+#include <dumux/multidomain/newtonsolver.hh>
+#include <dumux/multidomain/fvassembler.hh>
+#include <dumux/multidomain/traits.hh>
+
+#include <dumux/multidomain/facet/gridcreator.hh>
+#include <dumux/multidomain/facet/couplingmapper.hh>
+#include <dumux/multidomain/facet/couplingmanager.hh>
+
+#include <dumux/io/vtkoutputmodule.hh>
+
+// obtain/define some types to be used below in the property definitions and in main
+class TestTraits
+{
+    using BulkFVG = typename GET_PROP_TYPE(TTAG(OnePBulkTpfa), FVGridGeometry);
+    using FacetFVG = typename GET_PROP_TYPE(TTAG(OnePFacetTpfa), FVGridGeometry);
+    using EdgeFVG = typename GET_PROP_TYPE(TTAG(OnePEdgeTpfa), FVGridGeometry);
+public:
+    using MDTraits = Dumux::MultiDomainTraits<TTAG(OnePBulkTpfa), TTAG(OnePFacetTpfa), TTAG(OnePEdgeTpfa)>;
+    using CouplingMapper = Dumux::FacetCouplingThreeDomainMapper<BulkFVG, FacetFVG, EdgeFVG>;
+    using CouplingManager = Dumux::FacetCouplingThreeDomainManager<MDTraits, CouplingMapper>;
+};
+
+// set the coupling manager property in the sub-problems
+namespace Dumux {
+namespace Properties {
+
+NEW_PROP_TAG(CouplingManager);
+SET_TYPE_PROP(OnePBulkTpfa, CouplingManager, typename TestTraits::CouplingManager);
+SET_TYPE_PROP(OnePFacetTpfa, CouplingManager, typename TestTraits::CouplingManager);
+SET_TYPE_PROP(OnePEdgeTpfa, CouplingManager, typename TestTraits::CouplingManager);
+
+} // end namespace Properties
+} // end namespace Dumux
+
+int main(int argc, char** argv) try
+{
+    using namespace Dumux;
+
+    ////////////////////////////////////////////////////////////
+    ////////////////////////////////////////////////////////////
+
+    // initialize MPI, finalize is done automatically on exit
+    const auto& mpiHelper = Dune::MPIHelper::instance(argc, argv);
+
+    // print dumux start message
+    if (mpiHelper.rank() == 0)
+        DumuxMessage::print(/*firstCall=*/true);
+
+    // initialize parameter tree
+    Parameters::init(argc, argv);
+
+    //////////////////////////////////////////////////////////////////////
+    // try to create a grid (from the given grid file or the input file)
+    /////////////////////////////////////////////////////////////////////
+    using BulkProblemTypeTag = TTAG(OnePBulkTpfa);
+    using FacetProblemTypeTag = TTAG(OnePFacetTpfa);
+    using EdgeProblemTypeTag = TTAG(OnePEdgeTpfa);
+    using BulkGrid = typename GET_PROP_TYPE(BulkProblemTypeTag, Grid);
+    using FacetGrid = typename GET_PROP_TYPE(FacetProblemTypeTag, Grid);
+    using EdgeGrid = typename GET_PROP_TYPE(EdgeProblemTypeTag, Grid);
+
+    using GridCreator = FacetCouplingGridCreator<BulkGrid, FacetGrid, EdgeGrid>;
+    GridCreator gridCreator;
+    gridCreator.makeGrids(getParam<std::string>("Grid.File"));
+    gridCreator.loadBalance();
+
+    ////////////////////////////////////////////////////////////
+    // run stationary non-linear problem on this grid
+    ////////////////////////////////////////////////////////////
+
+    // we compute on the leaf grid views
+    const auto& bulkGridView = gridCreator.template grid<0>().leafGridView();
+    const auto& facetGridView = gridCreator.template grid<1>().leafGridView();
+    const auto& edgeGridView = gridCreator.template grid<2>().leafGridView();
+
+    // create the finite volume grid geometries
+    using BulkFVGridGeometry = typename GET_PROP_TYPE(BulkProblemTypeTag, FVGridGeometry);
+    using FacetFVGridGeometry = typename GET_PROP_TYPE(FacetProblemTypeTag, FVGridGeometry);
+    using EdgeFVGridGeometry = typename GET_PROP_TYPE(EdgeProblemTypeTag, FVGridGeometry);
+    auto bulkFvGridGeometry = std::make_shared<BulkFVGridGeometry>(bulkGridView);
+    auto facetFvGridGeometry = std::make_shared<FacetFVGridGeometry>(facetGridView);
+    auto edgeFvGridGeometry = std::make_shared<EdgeFVGridGeometry>(edgeGridView);
+    bulkFvGridGeometry->update();
+    facetFvGridGeometry->update();
+    edgeFvGridGeometry->update();
+
+    // the problems (boundary conditions)
+    using BulkProblem = typename GET_PROP_TYPE(BulkProblemTypeTag, Problem);
+    using FacetProblem = typename GET_PROP_TYPE(FacetProblemTypeTag, Problem);
+    using EdgeProblem = typename GET_PROP_TYPE(EdgeProblemTypeTag, Problem);
+    auto bulkSpatialParams = std::make_shared<typename BulkProblem::SpatialParams>(bulkFvGridGeometry, "Bulk");
+    auto bulkProblem = std::make_shared<BulkProblem>(bulkFvGridGeometry, bulkSpatialParams, "Bulk");
+    auto facetSpatialParams = std::make_shared<typename FacetProblem::SpatialParams>(facetFvGridGeometry, "Facet");
+    auto facetProblem = std::make_shared<FacetProblem>(facetFvGridGeometry, facetSpatialParams, "Facet");
+    auto edgeSpatialParams = std::make_shared<typename EdgeProblem::SpatialParams>(edgeFvGridGeometry, "Edge");
+    auto edgeProblem = std::make_shared<EdgeProblem>(edgeFvGridGeometry, edgeSpatialParams, "Edge");
+
+    // the solution vector
+    using Traits = typename TestTraits::MDTraits;
+    using SolutionVector = typename Traits::SolutionVector;
+    SolutionVector x;
+
+    static const auto bulkId = Traits::template DomainIdx<0>();
+    static const auto facetId = Traits::template DomainIdx<1>();
+    static const auto edgeId = Traits::template DomainIdx<2>();
+    x[bulkId].resize(bulkFvGridGeometry->numDofs());
+    x[facetId].resize(facetFvGridGeometry->numDofs());
+    x[edgeId].resize(edgeFvGridGeometry->numDofs());
+    bulkProblem->applyInitialSolution(x[bulkId]);
+    facetProblem->applyInitialSolution(x[facetId]);
+    edgeProblem->applyInitialSolution(x[edgeId]);
+
+    // the coupling mapper
+    using CouplingMapper = typename TestTraits::CouplingMapper;
+    auto couplingMapper = std::make_shared<CouplingMapper>();
+    couplingMapper->update(*bulkFvGridGeometry, *facetFvGridGeometry, *edgeFvGridGeometry, gridCreator);
+
+    // the coupling manager
+    using CouplingManager = typename TestTraits::CouplingManager;
+    auto couplingManager = std::make_shared<CouplingManager>();
+    couplingManager->init(bulkProblem, facetProblem, edgeProblem, couplingMapper, x);
+
+    // set coupling manager pointer in sub-problems
+    bulkProblem->setCouplingManager(couplingManager);
+    facetProblem->setCouplingManager(couplingManager);
+    edgeProblem->setCouplingManager(couplingManager);
+
+    // the grid variables
+    using BulkGridVariables = typename GET_PROP_TYPE(BulkProblemTypeTag, GridVariables);
+    using FacetGridVariables = typename GET_PROP_TYPE(FacetProblemTypeTag, GridVariables);
+    using EdgeGridVariables = typename GET_PROP_TYPE(EdgeProblemTypeTag, GridVariables);
+    auto bulkGridVariables = std::make_shared<BulkGridVariables>(bulkProblem, bulkFvGridGeometry);
+    auto facetGridVariables = std::make_shared<FacetGridVariables>(facetProblem, facetFvGridGeometry);
+    auto edgeGridVariables = std::make_shared<EdgeGridVariables>(edgeProblem, edgeFvGridGeometry);
+    bulkGridVariables->init(x[bulkId]);
+    facetGridVariables->init(x[facetId]);
+    edgeGridVariables->init(x[edgeId]);
+
+    // intialize the vtk output module
+    VtkOutputModule<BulkProblemTypeTag> bulkVtkWriter(*bulkProblem, *bulkFvGridGeometry, *bulkGridVariables, x[bulkId], bulkProblem->name());
+    VtkOutputModule<FacetProblemTypeTag> facetVtkWriter(*facetProblem, *facetFvGridGeometry, *facetGridVariables, x[facetId], facetProblem->name());
+    VtkOutputModule<EdgeProblemTypeTag> edgeVtkWriter(*edgeProblem, *edgeFvGridGeometry, *edgeGridVariables, x[edgeId], edgeProblem->name());
+
+    // Add model specific output fields
+    using BulkVtkOutputFields = typename GET_PROP_TYPE(BulkProblemTypeTag, VtkOutputFields);
+    using FacetVtkOutputFields = typename GET_PROP_TYPE(FacetProblemTypeTag, VtkOutputFields);
+    using EdgeVtkOutputFields = typename GET_PROP_TYPE(EdgeProblemTypeTag, VtkOutputFields);
+    BulkVtkOutputFields::init(bulkVtkWriter);
+    FacetVtkOutputFields::init(facetVtkWriter);
+    EdgeVtkOutputFields::init(edgeVtkWriter);
+    bulkVtkWriter.write(0.0);
+    facetVtkWriter.write(0.0);
+    edgeVtkWriter.write(0.0);
+
+    // the assembler
+    using Assembler = MultiDomainFVAssembler<Traits, CouplingManager, DiffMethod::numeric, /*implicit?*/true>;
+    auto assembler = std::make_shared<Assembler>( std::make_tuple(bulkProblem, facetProblem, edgeProblem),
+                                                  std::make_tuple(bulkFvGridGeometry, facetFvGridGeometry, edgeFvGridGeometry),
+                                                  std::make_tuple(bulkGridVariables, facetGridVariables, edgeGridVariables),
+                                                  couplingManager);
+
+    // the linear solver
+    using LinearSolver = ILU0BiCGSTABBackend;
+    auto linearSolver = std::make_shared<LinearSolver>();
+
+    // the non-linear solver
+    using NewtonSolver = Dumux::MultiDomainNewtonSolver<Assembler, LinearSolver, CouplingManager>;
+    auto newtonSolver = std::make_shared<NewtonSolver>(assembler, linearSolver, couplingManager);
+
+    // linearize & solve
+    newtonSolver->solve(x);
+
+    // update grid variables for output
+    bulkGridVariables->update(x[bulkId]);
+    facetGridVariables->update(x[facetId]);
+    edgeGridVariables->update(x[edgeId]);
+
+    // write vtk output
+    bulkVtkWriter.write(1.0);
+    facetVtkWriter.write(1.0);
+    edgeVtkWriter.write(1.0);
+
+    ////////////////////////////////////////////////////////////
+    // finalize, print dumux message to say goodbye
+    ////////////////////////////////////////////////////////////
+    if (mpiHelper.rank() == 0)
+        DumuxMessage::print(/*firstCall=*/false);
+
+    return 0;
+
+}
+catch (Dumux::ParameterException &e)
+{
+    std::cerr << std::endl << e << " ---> Abort!" << std::endl;
+    return 1;
+}
+catch (Dune::DGFException & e)
+{
+    std::cerr << "DGF exception thrown (" << e <<
+                 "). Most likely, the DGF file name is wrong "
+                 "or the DGF file is corrupted, "
+                 "e.g. missing hash at end of file or wrong number (dimensions) of entries."
+                 << " ---> Abort!" << std::endl;
+    return 2;
+}
+catch (Dune::Exception &e)
+{
+    std::cerr << "Dune reported error: " << e << " ---> Abort!" << std::endl;
+    return 3;
+}
+catch (...)
+{
+    std::cerr << "Unknown exception thrown! ---> Abort!" << std::endl;
+    return 4;
+}

test/multidomain/facet/CMakeLists.txt

+add_subdirectory(1p_1p)
 dune_symlink_to_source_files(FILES "grid.msh" "grid2.msh" "test_gridcreator.input" "2d_grid.msh" "3d_grid.msh")
 
 dune_add_test(NAME test_facetgridcreator_alu

test/references/facetcoupling_1p1p_threedomain_edge.vtp

+<?xml version="1.0"?>
+<VTKFile type="PolyData" version="0.1" byte_order="LittleEndian">
+  <PolyData>
+    <Piece NumberOfLines="4" NumberOfPoints="5">
+      <CellData Scalars="pressure">
+        <DataArray type="Float32" Name="pressure" NumberOfComponents="1" format="ascii">
+          1.50645 1.5054 1.50401 1.50294
+        </DataArray>
+        <DataArray type="Float32" Name="process rank" NumberOfComponents="1" format="ascii">
+          0 0 0 0
+        </DataArray>
+      </CellData>
+      <Points>
+        <DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="ascii">
+          0.5 0.25 0.5 0.5 0.375 0.5 0.5 0.5 0.5 0.5 0.625 0.5
+          0.5 0.75 0.5
+        </DataArray>
+      </Points>
+      <Lines>
+        <DataArray type="Int32" Name="connectivity" NumberOfComponents="1" format="ascii">
+          0 1 1 2 2 3 3 4
+        </DataArray>
+        <DataArray type="Int32" Name="offsets" NumberOfComponents="1" format="ascii">
+          2 4 6 8
+        </DataArray>
+      </Lines>
+    </Piece>
+  </PolyData>
+</VTKFile>

test/references/facetcoupling_1p1p_threedomain_facet.vtu

+<?xml version="1.0"?>
+<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian">
+  <UnstructuredGrid>
+    <Piece NumberOfCells="88" NumberOfPoints="57">
+      <CellData Scalars="pressure">
+        <DataArray type="Float32" Name="pressure" NumberOfComponents="1" format="ascii">
+          1.5078 1.50688 1.50091 1.50238 1.50942 1.50914 1.4994 1.49969 1.50845 1.50035 1.50284 1.50155
+          1.50638 1.50463 1.50764 1.50558 1.50654 1.50375 1.50468 1.50261 1.50154 1.50767 1.50697 1.50807
+          1.50246 1.50121 1.50929 1.5096 1.4997 1.49998 1.50563 1.50664 1.5038 1.50475 1.50271 1.50063
+          1.50304 1.50177 1.50648 1.50476 1.50778 1.50868 1.50173 1.50783 1.50985 1.50955 1.49979 1.49947
+          1.50702 1.50819 1.50245 1.50116 1.50566 1.50671 1.5038 1.50478 1.50269 1.50886 1.5067 1.508
+          1.50304 1.50484 1.50171 1.5005 1.50167 1.508 1.50843 1.5071 1.50141 1.50252 1.50997 1.50965
+          1.49987 1.49958 1.50905 1.50678 1.5081 1.50311 1.50491 1.50179 1.50067 1.5057 1.50679 1.50384
+          1.50484 1.50277 1.50178 1.50813
+        </DataArray>
+        <DataArray type="Float32" Name="process rank" NumberOfComponents="1" format="ascii">
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0 0 0 0 0 0 0 0 0
+          0 0 0 0
+        </DataArray>
+      </CellData>
+      <Points>
+        <DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="ascii">
+          0.5 0.25 0.5 0.411458 0.338542 0.5 0.375 0.25 0.5 0.5 0.375 0.5
+          0.5 0.75 0.5 0.375 0.75 0.5 0.418846 0.671266 0.5 0.5 0.625 0.5
+          0.25 0.25 0.5 0.329409 0.326542 0.5 0.25 0.375 0.5 0.25 0.75 0.5
+          0.341841 0.661125 0.5 0.25 0.625 0.5 0.25 0.5 0.5 0.377392 0.570206 0.5
+          0.373141 0.425983 0.5 0.5 0.5 0.5 0.5 0.338542 0.588542 0.5 0.25 0.625
+          0.5 0.671266 0.581154 0.5 0.75 0.625 0.5 0.25 0.75 0.5 0.326542 0.670591
+          0.5 0.375 0.75 0.5 0.75 0.75 0.5 0.661125 0.658159 0.5 0.625 0.75
+          0.5 0.425983 0.626859 0.5 0.570206 0.622608 0.5 0.5 0.75 0.5 0.25 0.25
+          0.5 0.25 0.375 0.5 0.326542 0.329409 0.5 0.375 0.25 0.5 0.75 0.25
+          0.5 0.625 0.25 0.5 0.661125 0.341841 0.5 0.75 0.375 0.5 0.338542 0.411458
+          0.5 0.671266 0.418846 0.5 0.425983 0.373141 0.5 0.570206 0.377392 0.5 0.5 0.25
+          0.588542 0.338542 0.5 0.625 0.25 0.5 0.625 0.75 0.5 0.581154 0.671266 0.5
+          0.75 0.25 0.5 0.670591 0.326542 0.5 0.75 0.375 0.5 0.75 0.75 0.5
+          0.75 0.625 0.5 0.658159 0.661125 0.5 0.626859 0.425983 0.5 0.75 0.5 0.5
+          0.622608 0.570206 0.5
+        </DataArray>
+      </Points>
+      <Cells>
+        <DataArray type="Int32" Name="connectivity" NumberOfComponents="1" format="ascii">
+          0 1 2 0 3 1 4 5 6 4 6 7
+          8 2 9 8 9 10 11 12 5 11 13 12
+          2 1 9 5 12 6 13 14 15 13 15 12
+          14 10 16 14 16 15 10 9 16 3 17 16
+          3 16 1 17 7 15 17 15 16 7 6 15
+          15 6 12 16 9 1 0 3 18 0 18 19
+          4 20 7 4 21 20 22 23 24 22 19 23
+          25 26 21 25 27 26 3 17 28 3 28 18
+          17 7 29 17 29 28 7 20 29 21 26 20
+          27 30 29 27 29 26 30 24 28 30 28 29
+          24 23 28 19 18 23 29 20 26 28 23 18
+          31 32 33 31 33 34 35 36 37 35 37 38
+          0 3 39 0 39 32 4 40 7 4 38 40
+          3 17 41 3 41 39 17 7 42 17 42 41
+          7 40 42 32 39 33 34 41 43 34 33 41
+          43 42 36 43 41 42 36 42 37 38 37 40
+          42 40 37 41 33 39 0 44 45 0 3 44
+          4 46 47 4 47 7 48 45 49 48 49 50
+          51 52 53 51 53 46 45 44 49 50 54 55
+          50 49 54 55 56 52 55 54 56 52 56 53
+          46 53 47 3 17 54 3 54 44 17 7 56
+          17 56 54 7 47 56 56 47 53 54 49 44
+        </DataArray>
+        <DataArray type="Int32" Name="offsets" NumberOfComponents="1" format="ascii">
+          3 6 9 12 15 18 21 24 27 30 33 36
+          39 42 45 48 51 54 57 60 63 66 69 72
+          75 78 81 84 87 90 93 96 99 102 105 108
+          111 114 117 120 123 126 129 132 135 138 141 144
+          147 150 153 156 159 162 165 168 171 174 177 180
+          183 186 189 192 195 198 201 204 207 210 213 216
+          219 222 225 228 231 234 237 240 243 246 249 252
+          255 258 261 264
+        </DataArray>
+        <DataArray type="UInt8" Name="types" NumberOfComponents="1" format="ascii">
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5 5 5 5 5 5 5 5 5
+          5 5 5 5
+        </DataArray>
+      </Cells>
+    </Piece>
+  </UnstructuredGrid>
+</VTKFile>