My Project
Loading...
Searching...
No Matches
FlowGenericVanguard.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*/
27#ifndef OPM_FLOW_GENERIC_VANGUARD_HPP
28#define OPM_FLOW_GENERIC_VANGUARD_HPP
29
30#include <dune/common/parallel/communication.hh>
31
32#include <opm/common/OpmLog/OpmLog.hpp>
33
34#include <opm/grid/common/GridEnums.hpp>
35
36#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
37
38#include <opm/simulators/utils/ParallelCommunication.hpp>
39
40#include <cassert>
41#include <memory>
42#include <optional>
43#include <string>
44#include <unordered_map>
45#include <utility>
46#include <vector>
47
48namespace Opm::Parameters {
49
50struct AllowDistributedWells { static constexpr bool value = false; };
51struct AllowSplittingInactiveWells { static constexpr bool value = true; };
52
53struct EclOutputInterval { static constexpr int value = -1; };
54struct EdgeWeightsMethod { static constexpr int value = 1; };
55struct EnableDryRun { static constexpr auto value = "auto"; };
56struct EnableEclOutput { static constexpr auto value = true; };
57struct EnableOpmRstFile { static constexpr bool value = false; };
58struct ExternalPartition { static constexpr auto* value = ""; };
59
60template<class Scalar>
61struct ImbalanceTol { static constexpr Scalar value = 1.1; };
62
63struct IgnoreKeywords { static constexpr auto value = ""; };
64struct InputSkipMode { static constexpr auto value = "100"; };
65struct MetisParams { static constexpr auto value = "default"; };
66
67#if HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
68struct NumJacobiBlocks { static constexpr int value = 0; };
69#endif // HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
70
71struct OwnerCellsFirst { static constexpr bool value = true; };
72struct ParsingStrictness { static constexpr auto value = "normal"; };
73struct ActionParsingStrictness { static constexpr auto value = "normal"; };
74
77struct PartitionMethod { static constexpr int value = 3; };
78struct AddCorners { static constexpr bool value = false; };
79struct NumOverlap { static constexpr int value = 1; };
80
81struct SchedRestart{ static constexpr bool value = false; };
82struct SerialPartitioning{ static constexpr bool value = false; };
83
84template<class Scalar>
85struct ZoltanImbalanceTol { static constexpr Scalar value = 1.1; };
86
87struct ZoltanPhgEdgeSizeThreshold { static constexpr auto value = 0.35; };
88
89struct ZoltanParams { static constexpr auto value = "graph"; };
90
91} // namespace Opm::Parameters
92
93namespace Opm {
94
95namespace Action { class State; }
96class Deck;
97class EclipseState;
99class ParseContext;
100class Schedule;
101class Python;
102class SummaryConfig;
103class SummaryState;
104class UDQState;
105class WellTestState;
106
108public:
109 using ParallelWellStruct = std::vector<std::pair<std::string,bool>>;
110
112 double setupTime_;
113 std::unique_ptr<UDQState> udqState_;
114 std::unique_ptr<Action::State> actionState_;
115 std::unique_ptr<SummaryState> summaryState_;
116 std::unique_ptr<WellTestState> wtestState_;
117 std::shared_ptr<EclipseState> eclState_;
118 std::shared_ptr<Schedule> eclSchedule_;
119 std::shared_ptr<SummaryConfig> eclSummaryConfig_;
120 };
121
122 static SimulationModelParams modelParams_;
123
130
136
137 static SimulationModelParams serializationTestParams();
138
145 static std::string canonicalDeckPath(const std::string& caseName);
146
150 double setupTime()
151 { return setupTime_; }
152
157 static void readDeck(const std::string& filename);
158
162 void defineSimulationModel(SimulationModelParams&& params);
163
167 const EclipseState& eclState() const
168 { return *eclState_; }
169
170 EclipseState& eclState()
171 { return *eclState_; }
172
176 const Schedule& schedule() const
177 { return *eclSchedule_; }
178
180 { return *eclSchedule_; }
181
187 { return *eclSummaryConfig_; }
188
197 { return *summaryState_; }
198
199 const SummaryState& summaryState() const
200 { return *summaryState_; }
201
207 Action::State& actionState()
208 { return *actionState_; }
209
210 const Action::State& actionState() const
211 { return *actionState_; }
212
219 { return *udqState_; }
220
221 const UDQState& udqState() const
222 { return *udqState_; }
223
224 std::unique_ptr<WellTestState> transferWTestState() {
225 return std::move(this->wtestState_);
226 }
227
228
235 const std::string& caseName() const
236 { return caseName_; }
237
241 Dune::EdgeWeightMethod edgeWeightsMethod() const
242 { return edgeWeightsMethod_; }
243
247 int numJacobiBlocks() const
248 {
249#if HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
250 return numJacobiBlocks_;
251#else
252 return 0;
253#endif // HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
254 }
255
259 bool ownersFirst() const
260 { return ownersFirst_; }
261
262#if HAVE_MPI
263 bool addCorners() const
264 { return addCorners_; }
265
266 int numOverlap() const
267 { return numOverlap_; }
268
272 Dune::PartitionMethod partitionMethod() const
273 { return partitionMethod_; }
274
279 bool serialPartitioning() const
280 { return serialPartitioning_; }
281
286 double imbalanceTol() const
287 {
289 OpmLog::info("The parameter --zoltan-imbalance-tol is deprecated "
290 "and has been renamed to --imbalance-tol, please "
291 "adjust your calls and scripts!");
292 return zoltanImbalanceTol_;
293 } else {
294 return imbalanceTol_;
295 }
296 }
297
298 const std::string& externalPartitionFile() const
299 {
300 return this->externalPartitionFile_;
301 }
302#endif // HAVE_MPI
303
308 { return enableDistributedWells_; }
309
314 bool enableEclOutput() const
315 { return enableEclOutput_; }
316
324 const ParallelWellStruct& parallelWells() const
325 { return parallelWells_; }
326
328 static void setCommunication(std::unique_ptr<Opm::Parallel::Communication> comm)
329 { comm_ = std::move(comm); }
330
332 static Parallel::Communication& comm()
333 {
334 assert(comm_ != nullptr);
335 return *comm_;
336 }
337
338 // Private to avoid pulling schedule in header.
339 // Static state is not serialized, only use for restart.
340 template<class Serializer>
341 void serializeOp(Serializer& serializer);
342
343 // Only compares dynamic state.
344 bool operator==(const FlowGenericVanguard& rhs) const;
345
346protected:
347 void updateOutputDir_(std::string outputDir,
349
350 void updateNOSIM_(std::string_view enableDryRun);
351
352 bool drsdtconEnabled() const;
353
354 std::unordered_map<std::size_t, const NumericalAquiferCell*> allAquiferCells() const;
355
356 void init();
357
358 template<class Scalar>
359 static void registerParameters_();
360
361 double setupTime_;
362
363 // These variables may be owned by both Python and the simulator
364 static std::unique_ptr<Parallel::Communication> comm_;
365
366 std::string caseName_;
367 std::string fileName_;
368 Dune::EdgeWeightMethod edgeWeightsMethod_;
369
370#if HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
371 int numJacobiBlocks_{0};
372#endif // HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
373
374 bool ownersFirst_;
375#if HAVE_MPI
376 bool addCorners_;
377 int numOverlap_;
378
379 Dune::PartitionMethod partitionMethod_;
381 double imbalanceTol_;
382
384 double zoltanImbalanceTol_;
386 std::string zoltanParams_;
387
388 std::string metisParams_;
389
390 std::string externalPartitionFile_{};
391#endif // HAVE_MPI
392
393 bool enableDistributedWells_;
394 bool enableEclOutput_;
395 bool allow_splitting_inactive_wells_;
396
397 std::string ignoredKeywords_;
398 std::optional<int> outputInterval_;
399 bool useMultisegmentWell_;
400 bool enableExperiments_;
401
402 std::unique_ptr<SummaryState> summaryState_;
403 std::unique_ptr<UDQState> udqState_;
404 std::unique_ptr<Action::State> actionState_;
405
406 // Observe that this instance is handled differently from the other state
407 // variables, it will only be initialized for a restart run. While
408 // initializing a restarted run this instance is transferred to the WGState
409 // member in the well model.
410 std::unique_ptr<WellTestState> wtestState_;
411
412 // these attributes point either to the internal or to the external version of the
413 // parser objects.
414 std::shared_ptr<Python> python;
415 // These variables may be owned by both Python and the simulator
416 std::shared_ptr<EclipseState> eclState_;
417 std::shared_ptr<Schedule> eclSchedule_;
418 std::shared_ptr<SummaryConfig> eclSummaryConfig_;
419
425 ParallelWellStruct parallelWells_;
426};
427
428} // namespace Opm
429
430#endif // OPM_FLOW_GENERIC_VANGUARD_HPP
Definition FlowGenericVanguard.hpp:107
static Parallel::Communication & comm()
Obtain global communicator.
Definition FlowGenericVanguard.hpp:332
int numJacobiBlocks() const
Number of blocks in the Block-Jacobi preconditioner.
Definition FlowGenericVanguard.hpp:247
const Schedule & schedule() const
Return a reference to the object that managages the ECL schedule.
Definition FlowGenericVanguard.hpp:176
double setupTime()
Returns the wall time required to set up the simulator before it was born.
Definition FlowGenericVanguard.hpp:150
static void readDeck(const std::string &filename)
Read a deck.
Definition FlowGenericVanguard.cpp:174
const SummaryConfig & summaryConfig() const
Return a reference to the object that determines which quantities ought to be put into the ECL summar...
Definition FlowGenericVanguard.hpp:186
ParallelWellStruct parallelWells_
Information about wells in parallel.
Definition FlowGenericVanguard.hpp:425
static void setCommunication(std::unique_ptr< Opm::Parallel::Communication > comm)
Set global communication.
Definition FlowGenericVanguard.hpp:328
bool enableDistributedWells() const
Whether perforations of a well might be distributed.
Definition FlowGenericVanguard.hpp:307
static std::string canonicalDeckPath(const std::string &caseName)
Returns the canonical path to a deck file.
Definition FlowGenericVanguard.cpp:191
~FlowGenericVanguard()
Destructor.
Action::State & actionState()
Returns the action state.
Definition FlowGenericVanguard.hpp:207
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition FlowGenericVanguard.hpp:167
FlowGenericVanguard()
Constructor.
Definition FlowGenericVanguard.cpp:103
void defineSimulationModel(SimulationModelParams &&params)
Set the simulation configuration objects.
Definition FlowGenericVanguard.cpp:162
bool enableEclOutput() const
Whether or not to emit result files that are compatible with a commercial reservoir simulator.
Definition FlowGenericVanguard.hpp:314
const std::string & caseName() const
Returns the name of the case.
Definition FlowGenericVanguard.hpp:235
bool ownersFirst() const
Parameter that decide if cells owned by rank are ordered before ghost cells.
Definition FlowGenericVanguard.hpp:259
Dune::EdgeWeightMethod edgeWeightsMethod() const
Parameter deciding the edge-weight strategy of the load balancer.
Definition FlowGenericVanguard.hpp:241
UDQState & udqState()
Returns the udq state.
Definition FlowGenericVanguard.hpp:218
const ParallelWellStruct & parallelWells() const
Retrieve collection (a vector of pairs) of well names and whether or not the corresponding well objec...
Definition FlowGenericVanguard.hpp:324
SummaryState & summaryState()
Returns the summary state.
Definition FlowGenericVanguard.hpp:196
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
Definition FlowGenericVanguard.hpp:111
Definition FlowGenericVanguard.hpp:73
Definition FlowGenericVanguard.hpp:78
Definition FlowGenericVanguard.hpp:50
Definition FlowGenericVanguard.hpp:51
Definition FlowGenericVanguard.hpp:53
Definition FlowGenericVanguard.hpp:54
Definition FlowGenericVanguard.hpp:55
Definition FlowGenericVanguard.hpp:56
Definition FlowGenericVanguard.hpp:57
Definition FlowGenericVanguard.hpp:58
Definition FlowGenericVanguard.hpp:63
Definition FlowGenericVanguard.hpp:61
Definition FlowGenericVanguard.hpp:64
Definition FlowGenericVanguard.hpp:65
Definition FlowGenericVanguard.hpp:79
Definition FlowGenericVanguard.hpp:71
Definition FlowGenericVanguard.hpp:72
0: simple, 1: Zoltan, 2: METIS, 3: Zoltan with a all cells of a well represented by one vertex in the...
Definition FlowGenericVanguard.hpp:77
Definition FlowGenericVanguard.hpp:81
Definition FlowGenericVanguard.hpp:82
Definition FlowGenericVanguard.hpp:85
Definition FlowGenericVanguard.hpp:89
Definition FlowGenericVanguard.hpp:87