My Project
Loading...
Searching...
No Matches
EclGenericWriter.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef OPM_ECL_GENERIC_WRITER_HPP
29#define OPM_ECL_GENERIC_WRITER_HPP
30
32
33#include <opm/simulators/flow/CollectDataOnIORank.hpp>
35#include <opm/simulators/timestepping/SimulatorReport.hpp>
36
37#include <map>
38#include <memory>
39#include <string>
40#include <utility>
41#include <vector>
42
43namespace Opm {
44
45class EclipseIO;
46class EclipseState;
47class InterRegFlowMap;
48class Inplace;
49struct NNCdata;
50class Schedule;
51class SummaryConfig;
52class SummaryState;
53class UDQState;
54
55} // namespace Opm
56
57namespace Opm { namespace Action {
58class State;
59}} // namespace Opm::Action
60
61namespace Opm {
62
63template <class Grid, class EquilGrid, class GridView, class ElementMapper, class Scalar>
65{
70
71public:
72 EclGenericWriter(const Schedule& schedule,
73 const EclipseState& eclState,
74 const SummaryConfig& summaryConfig,
75 const Grid& grid,
76 const EquilGrid* equilGrid,
77 const GridView& gridView,
81 bool enableEsmry);
82
83 const EclipseIO& eclIO() const;
84
85 void writeInit();
86
87 void setTransmissibilities(const TransmissibilityType* globalTrans)
88 {
89 globalTrans_ = globalTrans;
90 }
91
92 void setSubStepReport(const SimulatorReportSingle& report)
93 {
94 sub_step_report_ = report;
95 }
96 void setSimulationReport(const SimulatorReport& report)
97 {
98 simulation_report_ = report;
99 }
100
101 const std::vector<NNCdata>& getOutputNnc() const
102 {
103 return outputNnc_;
104 }
105
106 const CollectDataOnIORankType& collectOnIORank() const
107 {
108 return collectOnIORank_;
109 }
110
111 void extractOutputTransAndNNC(const std::function<unsigned int(unsigned int)>& map);
112
113protected:
114 const TransmissibilityType& globalTrans() const;
115 unsigned int gridEquilIdxToGridIdx(unsigned int elemIndex) const;
116
117 void doWriteOutput(const int reportStepNum,
118 const std::optional<int> timeStepNum,
119 const bool isSubStep,
120 data::Solution&& localCellData,
121 data::Wells&& localWellData,
122 data::GroupAndNetworkValues&& localGroupAndNetworkData,
123 data::Aquifers&& localAquiferData,
125 const Action::State& actionState,
126 const UDQState& udqState,
127 const SummaryState& summaryState,
128 const std::vector<Scalar>& thresholdPressure,
129 Scalar curTime,
130 Scalar nextStepSize,
131 bool doublePrecision,
132 bool isFlowsn,
133 std::array<FlowsData<double>,3>&& flowsn,
134 bool isFloresn,
135 std::array<FlowsData<double>, 3>&& floresn);
136
137 void evalSummary(int reportStepNum,
138 Scalar curTime,
139 const data::Wells& localWellData,
140 const data::WellBlockAveragePressures& localWBPData,
141 const data::GroupAndNetworkValues& localGroupAndNetworkData,
142 const std::map<int,data::AquiferData>& localAquiferData,
143 const std::map<std::pair<std::string, int>, double>& blockData,
144 const std::map<std::string, double>& miscSummaryData,
145 const std::map<std::string, std::vector<double>>& regionData,
146 const Inplace& inplace,
147 const Inplace& initialInPlace,
149 SummaryState& summaryState,
150 UDQState& udqState);
151
152 CollectDataOnIORankType collectOnIORank_;
153 const Grid& grid_;
154 const GridView& gridView_;
155 const Schedule& schedule_;
156 const EclipseState& eclState_;
157 std::unique_ptr<EclipseIO> eclIO_;
158 std::unique_ptr<TaskletRunner> taskletRunner_;
159 Scalar restartTimeStepSize_;
160 const TransmissibilityType* globalTrans_ = nullptr;
161 const Dune::CartesianIndexMapper<Grid>& cartMapper_;
162 const Dune::CartesianIndexMapper<EquilGrid>* equilCartMapper_;
163 const EquilGrid* equilGrid_;
164 SimulatorReportSingle sub_step_report_;
165 SimulatorReport simulation_report_;
166 mutable std::vector<NNCdata> outputNnc_;
167 mutable std::unique_ptr<data::Solution> outputTrans_;
168
169private:
170 void computeTrans_(const std::unordered_map<int,int>& cartesianToActive, const std::function<unsigned int(unsigned int)>& map) const;
171 std::vector<NNCdata> exportNncStructure_(const std::unordered_map<int,int>& cartesianToActive, const std::function<unsigned int(unsigned int)>& map) const;
172};
173
174} // namespace Opm
175
176#endif // OPM_ECL_GENERIC_WRITER_HPP
Definition CollectDataOnIORank.hpp:51
Definition CollectDataOnIORank.hpp:58
Definition EclGenericWriter.hpp:65
Inter-region flow accumulation maps for all region definition arrays.
Definition InterRegFlows.hpp:179
Definition Transmissibility.hpp:54
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Simple container for FLOWS data.
Definition FlowsData.hpp:32
A struct for returning timing data from a simulator to its caller.
Definition SimulatorReport.hpp:34
Definition SimulatorReport.hpp:100
Provides a mechanism to dispatch work to separate threads.