pnmvtkoutputmodule.hh 6 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
// -*- 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 3 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 PoreNetworkModels
 * \copydoc Dumux::PoreNetwork::VtkOutputModule
 */

#ifndef DUMUX_PNM_VTK_OUTPUT_MODULE_HH
#define DUMUX_PNM_VTK_OUTPUT_MODULE_HH

#include <dumux/io/vtkoutputmodule.hh>
#include "velocityoutput.hh"

namespace Dumux::PoreNetwork {

/*!
 * \ingroup PoreNetworkModels
 * \brief Adds vtk output fields specific to pore-network models.
 */
template<class GridVariables, class FluxVariables, class SolutionVector>
class VtkOutputModule : public Dumux::VtkOutputModule<GridVariables, SolutionVector>
{
    using ParentType = Dumux::VtkOutputModule<GridVariables, SolutionVector>;
    using Scalar = typename GridVariables::Scalar;
    using FluxVarsCache = typename GridVariables::GridFluxVariablesCache::FluxVariablesCache;

    struct ThroatFluxDataInfo { std::function<Scalar(const FluxVariables&, const FluxVarsCache&)> get; std::string name; };

public:

    //! The constructor
    VtkOutputModule(const GridVariables& gridVariables,
                    const SolutionVector& sol,
                    const std::string& name,
                    const std::string& paramGroup = "",
                    Dune::VTK::DataMode dm = Dune::VTK::conforming,
                    bool verbose = true)
    : ParentType(gridVariables, sol, name, paramGroup, dm, verbose)
    {
        if constexpr (GridVariables::VolumeVariables::numFluidPhases() >= 1)
        {
            // enable velocity output per default
            using VelocityOutput = VelocityOutput<GridVariables, FluxVariables>;
            this->addVelocityOutput(std::make_shared<VelocityOutput>(gridVariables));
        }
    }

    //! Output a scalar flux variable related to pore throats. This is basically a wrapper for the ParentType's addField method.
    //! \param f A function taking a Problem, FluxVariables and FluxVarsCache object and returning the desired scalar
    //! \param name The name of the vtk field
    void addFluxVariable(std::function<Scalar(const FluxVariables&, const FluxVarsCache&)>&& f, const std::string& name)
    {
        throatFluxDataInfo_.push_back(ThroatFluxDataInfo{f, name});
        const auto numElems = problem().gridGeometry().gridView().size(0);
        throatFluxData_.push_back(std::vector<Scalar>(numElems));
73
        ParentType::addField(throatFluxData_.back(), throatFluxDataInfo_.back().name, Vtk::FieldType::element);
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
    }

    //! Gather and process all required data and write them to a vtk file.
    void write(double time, Dune::VTK::OutputType type = Dune::VTK::ascii)
    {
        const auto gridView = problem().gridGeometry().gridView();
        const auto numElems = gridView.size(0);

        // resize all fields
        for (auto& data : throatFluxData_)
            data.resize(numElems);

        // iterate over all elements
        for (const auto& element : elements(gridView, Dune::Partitions::interior))
        {
            // make sure FVElementGeometry & vol vars are bound to the element
            auto fvElementGeometry = localView(problem().gridGeometry());
            fvElementGeometry.bindElement(element);

            const auto eIdx = problem().gridGeometry().elementMapper().index(element);

            auto elemVolVars = localView(this->gridVariables().curGridVolVars());
            auto elemFluxVarsCache = localView(this->gridVariables().gridFluxVarsCache());

            elemVolVars.bind(element, fvElementGeometry, this->sol());
            elemFluxVarsCache.bind(element, fvElementGeometry, elemVolVars);

            // treat the throat flux related data
            std::size_t dataIdx = 0;
            for (auto&& scvf : scvfs(fvElementGeometry))
            {
                if (!scvf.boundary())
                {
                    FluxVariables fluxVars;
                    fluxVars.init(problem(), element, fvElementGeometry, elemVolVars, scvf, elemFluxVarsCache);

                    dataIdx = 0;
                    for(auto& data : throatFluxData_)
                        data[eIdx] = throatFluxDataInfo_[dataIdx++].get(fluxVars, elemFluxVarsCache[scvf]);
                }
            }
        }

        // call the ParentType's write method to write out all data
        ParentType::write(time, type);

        // empty the data containers in order to save some memory
        auto clearAndShrink = [] (auto& data)
        {
            data.clear();
            data.shrink_to_fit();
        };

        for (auto& data : throatFluxData_)
            clearAndShrink(data);
    }

    //! Return a reference to the problem
    const auto& problem() const { return ParentType::problem(); }

private:
    std::vector<ThroatFluxDataInfo> throatFluxDataInfo_;
    std::list<std::vector<Scalar>> throatFluxData_;
};

} //namespace Dumux::PoreNetwork

#endif