My Project
Loading...
Searching...
No Matches
Schedule.hpp
1/*
2 Copyright 2013 Statoil ASA.
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 3 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#ifndef SCHEDULE_HPP
20#define SCHEDULE_HPP
21
22#include <cstddef>
23#include <ctime>
24#include <functional>
25#include <iosfwd>
26#include <map>
27#include <memory>
28#include <optional>
29#include <set>
30#include <string>
31#include <unordered_map>
32#include <utility>
33#include <vector>
34
35#include <opm/input/eclipse/EclipseState/Runspec.hpp>
36#include <opm/input/eclipse/Schedule/Action/WGNames.hpp>
37#include <opm/input/eclipse/Schedule/CompletedCells.hpp>
38#include <opm/input/eclipse/Schedule/Group/Group.hpp>
39#include <opm/input/eclipse/Schedule/MessageLimits.hpp>
40#include <opm/input/eclipse/Schedule/ScheduleDeck.hpp>
41#include <opm/input/eclipse/Schedule/ScheduleState.hpp>
42#include <opm/input/eclipse/Schedule/Well/PAvg.hpp>
43#include <opm/input/eclipse/Schedule/WriteRestartFileEvents.hpp>
44#include <opm/input/eclipse/Units/UnitSystem.hpp>
45
46namespace Opm
47{
48 namespace Action {
49 class ActionX;
50 class PyAction;
51 class State;
52 }
53 class ActiveGridCells;
54 class Deck;
55 class DeckKeyword;
56 class DeckRecord;
57 class EclipseState;
58 class ErrorGuard;
59 class FieldPropsManager;
60 class GasLiftOpt;
61 class GTNode;
62 class GuideRateConfig;
63 class GuideRateModel;
64 enum class InputErrorAction;
65 class ParseContext;
66 class Python;
67 class RPTConfig;
68 class SCHEDULESection;
69 class SegmentMatcher;
70 struct SimulatorUpdate;
71 class SummaryState;
72 class TracerConfig;
73 class UDQConfig;
74 class Well;
75 enum class WellGasInflowEquation;
76 class WellMatcher;
77 enum class WellProducerCMode;
78 enum class WellStatus;
79 class WellTestConfig;
80
81 namespace RestartIO { struct RstState; }
82
83
85 std::shared_ptr<const Python> m_python_handle;
86 std::string m_input_path;
87 ScheduleRestartInfo rst_info;
88 MessageLimits m_deck_message_limits;
89 UnitSystem m_unit_system;
90 Runspec m_runspec;
91 RSTConfig rst_config;
92 std::optional<int> output_interval;
93 double sumthin{-1.0};
94 bool rptonly{false};
95 bool gaslift_opt_active{false};
96 std::optional<OilVaporizationProperties> oilVap;
97
98 ScheduleStatic() = default;
99
100 explicit ScheduleStatic(std::shared_ptr<const Python> python_handle) :
101 m_python_handle(python_handle)
102 {}
103
104 ScheduleStatic(std::shared_ptr<const Python> python_handle,
105 const ScheduleRestartInfo& restart_info,
106 const Deck& deck,
107 const Runspec& runspec,
108 const std::optional<int>& output_interval_,
109 const ParseContext& parseContext,
110 ErrorGuard& errors);
111
112 template<class Serializer>
113 void serializeOp(Serializer& serializer)
114 {
115 serializer(m_deck_message_limits);
116 serializer(this->rst_info);
117 serializer(m_runspec);
118 serializer(m_unit_system);
119 serializer(this->m_input_path);
120 serializer(rst_info);
121 serializer(rst_config);
122 serializer(this->output_interval);
123 serializer(this->gaslift_opt_active);
124 }
125
126
127 static ScheduleStatic serializationTestObject();
128
129 bool operator==(const ScheduleStatic& other) const;
130 };
131
132
133 class Schedule {
134 public:
135
136 struct PairComp
137 {
138 bool operator()(const std::pair<std::string,KeywordLocation>& pair,
139 const std::string& str) const
140 {
141 return std::get<0>(pair) < str;
142 }
143 bool operator()(const std::pair<std::string,KeywordLocation>& pair1,
144 const std::pair<std::string,KeywordLocation>& pair2) const
145 {
146 return std::get<0>(pair1) < std::get<0>(pair2);
147 }
148 bool operator()(const std::string& str,
149 const std::pair<std::string,KeywordLocation>& pair) const
150 {
151 return str < std::get<0>(pair);
152 }
153 };
154
155 using WelSegsSet = std::set<std::pair<std::string,KeywordLocation>,PairComp>;
156
157 Schedule() = default;
158 explicit Schedule(std::shared_ptr<const Python> python_handle);
159 Schedule(const Deck& deck,
160 const EclipseGrid& grid,
161 const FieldPropsManager& fp,
162 const Runspec &runspec,
163 const ParseContext& parseContext,
164 ErrorGuard& errors,
165 std::shared_ptr<const Python> python,
166 const std::optional<int>& output_interval = {},
167 const RestartIO::RstState* rst = nullptr,
168 const TracerConfig* tracer_config = nullptr);
169
170 template<typename T>
171 Schedule(const Deck& deck,
172 const EclipseGrid& grid,
173 const FieldPropsManager& fp,
174 const Runspec &runspec,
175 const ParseContext& parseContext,
176 T&& errors,
177 std::shared_ptr<const Python> python,
178 const std::optional<int>& output_interval = {},
179 const RestartIO::RstState* rst = nullptr,
180 const TracerConfig* tracer_config = nullptr);
181
182 Schedule(const Deck& deck,
183 const EclipseGrid& grid,
184 const FieldPropsManager& fp,
185 const Runspec &runspec,
186 std::shared_ptr<const Python> python,
187 const std::optional<int>& output_interval = {},
188 const RestartIO::RstState* rst = nullptr,
189 const TracerConfig* tracer_config = nullptr);
190
191 Schedule(const Deck& deck,
192 const EclipseState& es,
193 const ParseContext& parseContext,
194 ErrorGuard& errors,
195 std::shared_ptr<const Python> python,
196 const std::optional<int>& output_interval = {},
197 const RestartIO::RstState* rst = nullptr);
198
199 template <typename T>
200 Schedule(const Deck& deck,
201 const EclipseState& es,
202 const ParseContext& parseContext,
203 T&& errors,
204 std::shared_ptr<const Python> python,
205 const std::optional<int>& output_interval = {},
206 const RestartIO::RstState* rst = nullptr);
207
208 Schedule(const Deck& deck,
209 const EclipseState& es,
210 std::shared_ptr<const Python> python,
211 const std::optional<int>& output_interval = {},
212 const RestartIO::RstState* rst = nullptr);
213
214 // The constructor *without* the Python arg should really only be used from Python itself
215 Schedule(const Deck& deck,
216 const EclipseState& es,
217 const std::optional<int>& output_interval = {},
218 const RestartIO::RstState* rst = nullptr);
219
220 static Schedule serializationTestObject();
221
222 /*
223 * If the input deck does not specify a start time, Eclipse's 1. Jan
224 * 1983 is defaulted
225 */
226 std::time_t getStartTime() const;
227 std::time_t posixStartTime() const;
228 std::time_t posixEndTime() const;
229 std::time_t simTime(std::size_t timeStep) const;
230 double seconds(std::size_t timeStep) const;
231 double stepLength(std::size_t timeStep) const;
232 std::optional<int> exitStatus() const;
233 const UnitSystem& getUnits() const { return this->m_static.m_unit_system; }
234 const Runspec& runspec() const { return this->m_static.m_runspec; }
235
236 std::size_t numWells() const;
237 std::size_t numWells(std::size_t timestep) const;
238 bool hasWell(const std::string& wellName) const;
239 bool hasWell(const std::string& wellName, std::size_t timeStep) const;
240
241 WellMatcher wellMatcher(std::size_t report_step) const;
242 std::function<std::unique_ptr<SegmentMatcher>()> segmentMatcherFactory(std::size_t report_step) const;
243 std::vector<std::string> wellNames(const std::string& pattern, std::size_t timeStep, const std::vector<std::string>& matching_wells = {}) const;
244 std::vector<std::string> wellNames(const std::string& pattern) const;
245 std::vector<std::string> wellNames(std::size_t timeStep) const;
246 std::vector<std::string> wellNames() const;
247
248 bool hasGroup(const std::string& groupName, std::size_t timeStep) const;
249 std::vector<std::string> groupNames(const std::string& pattern, std::size_t timeStep) const;
250 std::vector<std::string> groupNames(std::size_t timeStep) const;
251 std::vector<std::string> groupNames(const std::string& pattern) const;
252 std::vector<std::string> groupNames() const;
253 /*
254 The restart_groups function returns a vector of groups pointers which
255 is organized as follows:
256
257 1. The number of elements is WELLDIMS::MAXGROUPS + 1
258 2. The elements are sorted according to group.insert_index().
259 3. If there are less than WELLDIMS::MAXGROUPS nullptr is used.
260 4. The very last element corresponds to the FIELD group.
261 */
262 std::vector<const Group*> restart_groups(std::size_t timeStep) const;
263
264 std::vector<std::string> changed_wells(std::size_t reportStep) const;
265 const Well& getWell(std::size_t well_index, std::size_t timeStep) const;
266 const Well& getWell(const std::string& wellName, std::size_t timeStep) const;
267 const Well& getWellatEnd(const std::string& well_name) const;
268 // get the list of the constant flux aquifer specified in the whole schedule
269 std::unordered_set<int> getAquiferFluxSchedule() const;
270 std::vector<Well> getWells(std::size_t timeStep) const;
271 std::vector<Well> getWellsatEnd() const;
272 void shut_well(const std::string& well_name, std::size_t report_step);
273 void stop_well(const std::string& well_name, std::size_t report_step);
274 void open_well(const std::string& well_name, std::size_t report_step);
275 void applyWellProdIndexScaling(const std::string& well_name, const std::size_t reportStep, const double scalingFactor);
276
277 std::vector<const Group*> getChildGroups2(const std::string& group_name, std::size_t timeStep) const;
278 std::vector<Well> getChildWells2(const std::string& group_name, std::size_t timeStep) const;
279 WellProducerCMode getGlobalWhistctlMmode(std::size_t timestep) const;
280
281 const UDQConfig& getUDQConfig(std::size_t timeStep) const;
282 void evalAction(const SummaryState& summary_state, std::size_t timeStep);
283
284 GTNode groupTree(std::size_t report_step) const;
285 GTNode groupTree(const std::string& root_node, std::size_t report_step) const;
286 const Group& getGroup(const std::string& groupName, std::size_t timeStep) const;
287
288 std::optional<std::size_t> first_RFT() const;
289 /*
290 Will remove all completions which are connected to cell which is not
291 active. Will scan through all wells and all timesteps.
292 */
293 void filterConnections(const ActiveGridCells& grid);
294 std::size_t size() const;
295
296 bool write_rst_file(std::size_t report_step) const;
297 const std::map< std::string, int >& rst_keywords( size_t timestep ) const;
298
299 /*
300 The applyAction() is invoked from the simulator *after* an ACTIONX has
301 evaluated to true. The return value is a small structure with
302 'information' which the simulator should take into account when
303 updating internal datastructures after the ACTIONX keywords have been
304 applied.
305 */
306 SimulatorUpdate applyAction(std::size_t reportStep, const Action::ActionX& action, const std::vector<std::string>& matching_wells, const std::unordered_map<std::string, double>& wellpi);
307 /*
308 The runPyAction() will run the Python script in a PYACTION keyword. In
309 the case of Schedule updates the recommended way of doing that from
310 PYACTION is to invoke a "normal" ACTIONX keyword internally from the
311 Python code. he return value from runPyAction() comes from such a
312 internal ACTIONX.
313 */
314 SimulatorUpdate runPyAction(std::size_t reportStep, const Action::PyAction& pyaction, Action::State& action_state, EclipseState& ecl_state, SummaryState& summary_state);
315
316
317 const GasLiftOpt& glo(std::size_t report_step) const;
318
319 bool operator==(const Schedule& data) const;
320 std::shared_ptr<const Python> python() const;
321
322
323 const ScheduleState& back() const;
324 const ScheduleState& operator[](std::size_t index) const;
325 std::vector<ScheduleState>::const_iterator begin() const;
326 std::vector<ScheduleState>::const_iterator end() const;
327 void create_next(const time_point& start_time, const std::optional<time_point>& end_time);
328 void create_next(const ScheduleBlock& block);
329 void create_first(const time_point& start_time, const std::optional<time_point>& end_time);
330
331
332 /*
333 The cmp() function compares two schedule instances in a context aware
334 manner. Floating point numbers are compared with a tolerance. The
335 purpose of this comparison function is to implement regression tests
336 for the schedule instances created by loading a restart file.
337 */
338 static bool cmp(const Schedule& sched1, const Schedule& sched2, std::size_t report_step);
339 void applyKeywords(std::vector<DeckKeyword*>& keywords, std::size_t timeStep);
340
341 template<class Serializer>
342 void serializeOp(Serializer& serializer)
343 {
344 serializer(this->m_static);
345 serializer(this->m_sched_deck);
346 serializer(this->action_wgnames);
347 serializer(this->exit_status);
348 serializer(this->snapshots);
349 serializer(this->restart_output);
350 serializer(this->completed_cells);
351
352 this->template pack_unpack<PAvg>(serializer);
353 this->template pack_unpack<WellTestConfig>(serializer);
354 this->template pack_unpack<GConSale>(serializer);
355 this->template pack_unpack<GConSump>(serializer);
356 this->template pack_unpack<GroupEconProductionLimits>(serializer);
357 this->template pack_unpack<WListManager>(serializer);
358 this->template pack_unpack<Network::ExtNetwork>(serializer);
359 this->template pack_unpack<Network::Balance>(serializer);
360 this->template pack_unpack<RPTConfig>(serializer);
361 this->template pack_unpack<Action::Actions>(serializer);
362 this->template pack_unpack<UDQActive>(serializer);
363 this->template pack_unpack<UDQConfig>(serializer);
364 this->template pack_unpack<NameOrder>(serializer);
365 this->template pack_unpack<GroupOrder>(serializer);
366 this->template pack_unpack<GuideRateConfig>(serializer);
367 this->template pack_unpack<GasLiftOpt>(serializer);
368 this->template pack_unpack<RFTConfig>(serializer);
369 this->template pack_unpack<RSTConfig>(serializer);
370
371 this->template pack_unpack_map<int, VFPProdTable>(serializer);
372 this->template pack_unpack_map<int, VFPInjTable>(serializer);
373 this->template pack_unpack_map<std::string, Group>(serializer);
374 this->template pack_unpack_map<std::string, Well>(serializer);
375 }
376
377 template <typename T, class Serializer>
378 void pack_unpack(Serializer& serializer) {
379 std::vector<T> value_list;
380 std::vector<std::size_t> index_list;
381
382 if (serializer.isSerializing())
383 this->template pack_state<T>(value_list, index_list);
384
385 serializer(value_list);
386 serializer(index_list);
387
388 if (!serializer.isSerializing())
389 this->template unpack_state<T>(value_list, index_list);
390 }
391
392 template <typename T>
393 std::vector<std::pair<std::size_t, T>> unique() const {
394 std::vector<std::pair<std::size_t, T>> values;
395 for (std::size_t index = 0; index < this->snapshots.size(); index++) {
396 const auto& member = this->snapshots[index].get<T>();
397 const auto& value = member.get();
398 if (values.empty() || !(value == values.back().second))
399 values.push_back( std::make_pair(index, value));
400 }
401 return values;
402 }
403
404
405 template <typename T>
406 void pack_state(std::vector<T>& value_list, std::vector<std::size_t>& index_list) const {
407 auto unique_values = this->template unique<T>();
408 for (auto& [index, value] : unique_values) {
409 value_list.push_back( std::move(value) );
410 index_list.push_back( index );
411 }
412 }
413
414
415 template <typename T>
416 void unpack_state(const std::vector<T>& value_list, const std::vector<std::size_t>& index_list) {
417 std::size_t unique_index = 0;
418 while (unique_index < value_list.size()) {
419 const auto& value = value_list[unique_index];
420 const auto& first_index = index_list[unique_index];
421 auto last_index = this->snapshots.size();
422 if (unique_index < (value_list.size() - 1))
423 last_index = index_list[unique_index + 1];
424
425 auto& target_state = this->snapshots[first_index];
426 target_state.get<T>().update( std::move(value) );
427 for (std::size_t index=first_index + 1; index < last_index; index++)
428 this->snapshots[index].get<T>().update( target_state.get<T>() );
429
430 unique_index++;
431 }
432 }
433
434
435 template <typename K, typename T, class Serializer>
436 void pack_unpack_map(Serializer& serializer) {
437 std::vector<T> value_list;
438 std::vector<std::size_t> index_list;
439
440 if (serializer.isSerializing())
441 pack_map<K,T>(value_list, index_list);
442
443 serializer(value_list);
444 serializer(index_list);
445
446 if (!serializer.isSerializing())
447 unpack_map<K,T>(value_list, index_list);
448 }
449
450
451 template <typename K, typename T>
452 void pack_map(std::vector<T>& value_list,
453 std::vector<std::size_t>& index_list) {
454
455 const auto& last_map = this->snapshots.back().get_map<K,T>();
456 std::vector<K> key_list{ last_map.keys() };
457 std::unordered_map<K,T> current_value;
458
459 for (std::size_t index = 0; index < this->snapshots.size(); index++) {
460 auto& state = this->snapshots[index];
461 const auto& current_map = state.template get_map<K,T>();
462 for (const auto& key : key_list) {
463 auto& value = current_map.get_ptr(key);
464 if (value) {
465 auto it = current_value.find(key);
466 if (it == current_value.end() || !(*value == it->second)) {
467 value_list.push_back( *value );
468 index_list.push_back( index );
469
470 current_value[key] = *value;
471 }
472 }
473 }
474 }
475 }
476
477
478 template <typename K, typename T>
479 void unpack_map(const std::vector<T>& value_list,
480 const std::vector<std::size_t>& index_list) {
481
482 std::unordered_map<K, std::vector<std::pair<std::size_t, T>>> storage;
483 for (std::size_t storage_index = 0; storage_index < value_list.size(); storage_index++) {
484 const auto& value = value_list[storage_index];
485 const auto& time_index = index_list[storage_index];
486
487 storage[ value.name() ].emplace_back( time_index, value );
488 }
489
490 for (const auto& [key, values] : storage) {
491 for (std::size_t unique_index = 0; unique_index < values.size(); unique_index++) {
492 const auto& [time_index, value] = values[unique_index];
493 auto last_index = this->snapshots.size();
494 if (unique_index < (values.size() - 1))
495 last_index = values[unique_index + 1].first;
496
497 auto& map_value = this->snapshots[time_index].template get_map<K,T>();
498 map_value.update(std::move(value));
499
500 for (std::size_t index=time_index + 1; index < last_index; index++) {
501 auto& forward_map = this->snapshots[index].template get_map<K,T>();
502 forward_map.update( key, map_value );
503 }
504 }
505 }
506 }
507
508 friend std::ostream& operator<<(std::ostream& os, const Schedule& sched);
509 void dump_deck(std::ostream& os) const;
510
511 private:
512 struct HandlerContext {
513
514 const ScheduleBlock& block;
515 const DeckKeyword& keyword;
516 const std::size_t currentStep;
517 const std::vector<std::string>& matching_wells;
518 const bool actionx_mode;
519 const ParseContext& parseContext;
520 ErrorGuard& errors;
521 SimulatorUpdate* sim_update{nullptr};
522 const std::unordered_map<std::string, double>* target_wellpi{nullptr};
523 std::unordered_map<std::string, double>* wpimult_global_factor{nullptr};
524 WelSegsSet* welsegs_wells{nullptr};
525 std::set<std::string>* compsegs_wells{nullptr};
526 const ScheduleGrid& grid;
527
530 HandlerContext(const ScheduleBlock& block_,
531 const DeckKeyword& keyword_,
532 const ScheduleGrid& grid_,
533 const std::size_t currentStep_,
534 const std::vector<std::string>& matching_wells_,
535 bool actionx_mode_,
536 const ParseContext& parseContext_,
537 ErrorGuard& errors_,
538 SimulatorUpdate* sim_update_,
539 const std::unordered_map<std::string, double>* target_wellpi_,
540 std::unordered_map<std::string, double>* wpimult_global_factor_,
541 WelSegsSet* welsegs_wells_,
542 std::set<std::string>* compsegs_wells_)
543 : block(block_)
544 , keyword(keyword_)
545 , currentStep(currentStep_)
546 , matching_wells(matching_wells_)
547 , actionx_mode(actionx_mode_)
548 , parseContext(parseContext_)
549 , errors(errors_)
550 , sim_update(sim_update_)
551 , target_wellpi(target_wellpi_)
552 , wpimult_global_factor(wpimult_global_factor_)
553 , welsegs_wells(welsegs_wells_)
554 , compsegs_wells(compsegs_wells_)
555 , grid(grid_)
556 {}
557
558 void affected_well(const std::string& well_name);
559 void record_well_structure_change();
560
562 void welsegs_handled(const std::string& well_name)
563 {
564 if (welsegs_wells)
565 welsegs_wells->insert({well_name, keyword.location()});
566 }
567
569 void compsegs_handled(const std::string& well_name)
570 {
571 if (compsegs_wells)
572 compsegs_wells->insert(well_name);
573 }
574
575 };
576
577 // Please update the member functions
578 // - operator==(const Schedule&) const
579 // - serializationTestObject()
580 // - serializeOp(Serializer&)
581 // when you update/change this list of data members.
582 ScheduleStatic m_static;
583 ScheduleDeck m_sched_deck;
584 Action::WGNames action_wgnames;
585 std::optional<int> exit_status;
586 std::vector<ScheduleState> snapshots;
587 WriteRestartFileEvents restart_output;
588 CompletedCells completed_cells;
589
590 void load_rst(const RestartIO::RstState& rst,
591 const TracerConfig& tracer_config,
592 const ScheduleGrid& grid,
593 const FieldPropsManager& fp);
594 void addWell(Well well);
595 void addWell(const std::string& wellName,
596 const std::string& group,
597 int headI,
598 int headJ,
599 Phase preferredPhase,
600 const std::optional<double>& refDepth,
601 double drainageRadius,
602 bool allowCrossFlow,
603 bool automaticShutIn,
604 int pvt_table,
605 WellGasInflowEquation gas_inflow,
606 std::size_t timeStep,
607 Connection::Order wellConnectionOrder);
608 bool updateWPAVE(const std::string& wname, std::size_t report_step, const PAvg& pavg);
609
610 void updateGuideRateModel(const GuideRateModel& new_model, std::size_t report_step);
611 GTNode groupTree(const std::string& root_node, std::size_t report_step, std::size_t level, const std::optional<std::string>& parent_name) const;
612 bool checkGroups(const ParseContext& parseContext, ErrorGuard& errors);
613 bool updateWellStatus( const std::string& well, std::size_t reportStep, WellStatus status, std::optional<KeywordLocation> = {});
614 void addWellToGroup( const std::string& group_name, const std::string& well_name , std::size_t timeStep);
615 void iterateScheduleSection(std::size_t load_start,
616 std::size_t load_end,
617 const ParseContext& parseContext,
618 ErrorGuard& errors,
619 const ScheduleGrid& grid,
620 const std::unordered_map<std::string, double> * target_wellpi,
621 const std::string& prefix,
622 const bool log_to_debug = false);
623 void addACTIONX(const Action::ActionX& action);
624 void addGroupToGroup( const std::string& parent_group, const std::string& child_group);
625 void addGroup(const std::string& groupName , std::size_t timeStep);
626 void addGroup(Group group);
627 void addGroup(const RestartIO::RstGroup& rst_group, std::size_t timeStep);
628 void addWell(const std::string& wellName, const DeckRecord& record, std::size_t timeStep, Connection::Order connection_order);
629 void checkIfAllConnectionsIsShut(std::size_t currentStep);
630 void end_report(std::size_t report_step);
633 void handleKeyword(std::size_t currentStep,
634 const ScheduleBlock& block,
635 const DeckKeyword& keyword,
636 const ParseContext& parseContext,
637 ErrorGuard& errors,
638 const ScheduleGrid& grid,
639 const std::vector<std::string>& matching_wells,
640 bool actionx_mode,
641 SimulatorUpdate* sim_update,
642 const std::unordered_map<std::string, double>* target_wellpi,
643 std::unordered_map<std::string, double>* wpimult_global_factor = nullptr,
644 WelSegsSet* welsegs_wells = nullptr,
645 std::set<std::string>* compsegs_wells = nullptr);
646
647 void prefetch_cell_properties(const ScheduleGrid& grid, const DeckKeyword& keyword);
648 void store_wgnames(const DeckKeyword& keyword);
649 std::vector<std::string> wellNames(const std::string& pattern, const HandlerContext& context);
650 std::vector<std::string> wellNames(const std::string& pattern, std::size_t timeStep, const std::vector<std::string>& matching_wells, InputErrorAction error_action, ErrorGuard& errors, const KeywordLocation& location) const;
651 void invalidNamePattern( const std::string& namePattern, const HandlerContext& context) const;
652 static std::string formatDate(std::time_t t);
653 std::string simulationDays(std::size_t currentStep) const;
654 void applyGlobalWPIMULT( const std::unordered_map<std::string, double>& wpimult_global_factor);
655
656 bool must_write_rst_file(std::size_t report_step) const;
657
658 void applyEXIT(const DeckKeyword&, std::size_t currentStep);
659 SimulatorUpdate applyAction(std::size_t reportStep, const std::string& action_name, const std::vector<std::string>& matching_wells);
660
676 bool handleNormalKeyword(HandlerContext& handlerContext);
677
678 // Keyword Handlers
679 void handlePYACTION(const DeckKeyword&);
680 void handleWELPIRuntime(HandlerContext&);
681
682 // Normal keyword handlers -- in KeywordHandlers.cpp
683
684 void handleAQUCT (HandlerContext&);
685 void handleAQUFETP (HandlerContext&);
686 void handleAQUFLUX (HandlerContext&);
687 void handleBCProp (HandlerContext&);
688 void handleBRANPROP (HandlerContext&);
689 void handleCOMPDAT (HandlerContext&);
690 void handleCOMPLUMP (HandlerContext&);
691 void handleCOMPORD (HandlerContext&);
692 void handleCOMPSEGS (HandlerContext&);
693 void handleCOMPTRAJ (HandlerContext&);
694 void handleCSKIN (HandlerContext&);
695 void handleDRSDT (HandlerContext&);
696 void handleDRSDTCON (HandlerContext&);
697 void handleDRSDTR (HandlerContext&);
698 void handleDRVDT (HandlerContext&);
699 void handleDRVDTR (HandlerContext&);
700 void handleEXIT (HandlerContext&);
701 void handleGCONINJE (HandlerContext&);
702 void handleGCONPROD (HandlerContext&);
703 void handleGCONSALE (HandlerContext&);
704 void handleGCONSUMP (HandlerContext&);
705 void handleGECON (HandlerContext&);
706 void handleGEFAC (HandlerContext&);
707 void handleGEOKeyword(HandlerContext&);
708 void handleGLIFTOPT (HandlerContext&);
709 void handleGPMAINT (HandlerContext&);
710 void handleGRUPNET (HandlerContext&);
711 void handleGRUPTREE (HandlerContext&);
712 void handleGUIDERAT (HandlerContext&);
713 void handleLIFTOPT (HandlerContext&);
714 void handleLINCOM (HandlerContext&);
715 void handleMESSAGES (HandlerContext&);
716 void handleMXUNSUPP (HandlerContext&);
717 void handleNETBALAN (HandlerContext&);
718 void handleNEXTSTEP (HandlerContext&);
719 void handleNODEPROP (HandlerContext&);
720 void handleNUPCOL (HandlerContext&);
721 void handleRPTONLY (HandlerContext&);
722 void handleRPTONLYO (HandlerContext&);
723 void handleRPTRST (HandlerContext&);
724 void handleRPTSCHED (HandlerContext&);
725 void handleTUNING (HandlerContext&);
726 void handleSAVE (HandlerContext&);
727 void handleSUMTHIN (HandlerContext&);
728 void handleUDQ (HandlerContext&);
729 void handleVAPPARS (HandlerContext&);
730 void handleVFPINJ (HandlerContext&);
731 void handleVFPPROD (HandlerContext&);
732 void handleWCONHIST (HandlerContext&);
733 void handleWCONINJE (HandlerContext&);
734 void handleWCONINJH (HandlerContext&);
735 void handleWCONPROD (HandlerContext&);
736 void handleWECON (HandlerContext&);
737 void handleWEFAC (HandlerContext&);
738 void handleWELOPEN (HandlerContext&);
739 void handleWELPI (HandlerContext&);
740 void handleWELSEGS (HandlerContext&);
741 void handleWELSPECS (HandlerContext&);
742 void handleWELTARG (HandlerContext&);
743 void handleWELTRAJ (HandlerContext&);
744 void handleWFOAM (HandlerContext&);
745 void handleWGRUPCON (HandlerContext&);
746 void handleWHISTCTL (HandlerContext&);
747 void handleWINJTEMP (HandlerContext&);
748 void handleWLIFTOPT (HandlerContext&);
749 void handleWLIST (HandlerContext&);
750 void handleWMICP (HandlerContext&);
751 void handleWPAVE (HandlerContext&);
752 void handleWPAVEDEP (HandlerContext&);
753 void handleWVFPDP (HandlerContext&);
754 void handleWVFPEXP (HandlerContext&);
755 void handleWWPAVE (HandlerContext&);
756 void handleWPIMULT (HandlerContext&);
757 void handleWINJCLN (HandlerContext&);
758 void handleWINJDAM (HandlerContext&);
759 void handleWINJFCNC (HandlerContext&);
760 void handleWINJMULT (HandlerContext&);
761 void handleWPMITAB (HandlerContext&);
762 void handleWPOLYMER (HandlerContext&);
763 void handleWRFT (HandlerContext&);
764 void handleWRFTPLT (HandlerContext&);
765 void handleWSALT (HandlerContext&);
766 void handleWSEGITER (HandlerContext&);
767 void handleWSEGSICD (HandlerContext&);
768 void handleWSEGAICD (HandlerContext&);
769 void handleWSEGVALV (HandlerContext&);
770 void handleWSKPTAB (HandlerContext&);
771 void handleWSOLVENT (HandlerContext&);
772 void handleWTEMP (HandlerContext&);
773 void handleWTEST (HandlerContext&);
774 void handleWTMULT (HandlerContext&);
775 void handleWTRACER (HandlerContext&);
776 };
777}
778
779#endif
Definition Deck.hpp:49
About cell information and dimension: The actual grid information is held in a pointer to an ERT ecl_...
Definition EclipseGrid.hpp:55
Definition ErrorGuard.hpp:29
Definition FieldPropsManager.hpp:41
Definition MessageLimits.hpp:28
Definition ParseContext.hpp:84
Definition RSTConfig.hpp:197
Definition Runspec.hpp:449
Definition Schedule.hpp:133
Class for (de-)serializing.
Definition Serializer.hpp:84
Definition TracerConfig.hpp:33
Definition UnitSystem.hpp:33
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition state.hpp:54
Definition ScheduleDeck.hpp:92
Definition Schedule.hpp:84
Definition Schedule.hpp:137