96 struct MPI_Comm_Deleter;
116 using ModelParameters =
typename Model::ModelParameters;
117 using SolverParameters =
typename Solver::SolverParameters;
123 : simulator_(simulator)
126 simulator_.vanguard().eclState().getIOConfig(),
127 Parameters::Get<Parameters::SaveStep>(),
128 Parameters::Get<Parameters::LoadStep>(),
129 Parameters::Get<Parameters::SaveFile>(),
130 Parameters::Get<Parameters::LoadFile>())
134 this->terminalOutput_ =
false;
135 if (this->grid().comm().rank() == 0) {
139 [compNames =
typename Model::ComponentName{}](
const int compIdx)
140 {
return std::string_view { compNames.name(compIdx) }; }
143 if (!simulator_.vanguard().eclState().getIOConfig().initOnly()) {
144 this->convergence_output_.
145 startThread(this->simulator_.vanguard().eclState(),
147 R
"(OutputExtraConvergenceInfo (--output-extra-convergence-info))",
159 static void registerParameters()
161 ModelParameters::registerParameters();
162 SolverParameters::registerParameters();
163 TimeStepper::registerParameters();
164 detail::registerSimulatorParameters();
176#ifdef RESERVOIR_COUPLING_ENABLED
177 SimulatorReport
run(SimulatorTimer& timer,
int argc,
char** argv)
179 init(timer, argc, argv);
188 simulator_.model().invalidateAndUpdateIntensiveQuantities(0);
190 while (!timer.
done()) {
191 simulator_.problem().writeReports(timer);
192 bool continue_looping = runStep(timer);
193 if (!continue_looping)
break;
195 simulator_.problem().writeReports(timer);
197#ifdef RESERVOIR_COUPLING_ENABLED
200 if (this->reservoirCouplingMaster_) {
201 this->reservoirCouplingMaster_->sendTerminateAndDisconnect();
203 else if (this->reservoirCouplingSlave_ && !this->reservoirCouplingSlave_->terminated()) {
207 this->reservoirCouplingSlave_->receiveTerminateAndDisconnect();
214#ifdef RESERVOIR_COUPLING_ENABLED
219 bool checkRunningAsReservoirCouplingMaster()
221 for (std::size_t report_step = 0; report_step < this->schedule().size(); ++report_step) {
222 auto rescoup = this->schedule()[report_step].rescoup();
223 auto slave_count = rescoup.slaveCount();
224 auto master_group_count = rescoup.masterGroupCount();
228 if (slave_count > 0) {
231 else if (master_group_count > 0) {
233 throw ReservoirCouplingError(
234 "Inconsistent reservoir coupling master schedule: "
235 "Master group count is greater than 0 but slave count is 0"
243#ifdef RESERVOIR_COUPLING_ENABLED
245 void init(
const SimulatorTimer& timer,
int argc,
char** argv)
247 auto slave_mode = Parameters::Get<Parameters::Slave>();
249 this->reservoirCouplingSlave_ =
250 std::make_unique<ReservoirCouplingSlave<Scalar>>(
252 this->schedule(), timer
254 this->reservoirCouplingSlave_->sendAndReceiveInitialData();
255 this->simulator_.setReservoirCouplingSlave(this->reservoirCouplingSlave_.get());
256 wellModel_().setReservoirCouplingSlave(this->reservoirCouplingSlave_.get());
259 auto master_mode = checkRunningAsReservoirCouplingMaster();
261 this->reservoirCouplingMaster_ =
262 std::make_unique<ReservoirCouplingMaster<Scalar>>(
267 this->simulator_.setReservoirCouplingMaster(this->reservoirCouplingMaster_.get());
268 wellModel_().setReservoirCouplingMaster(this->reservoirCouplingMaster_.get());
272 void init(
const SimulatorTimer& timer)
275 simulator_.setEpisodeIndex(-1);
278 solverTimer_ = std::make_unique<time::StopWatch>();
279 totalTimer_ = std::make_unique<time::StopWatch>();
280 totalTimer_->start();
283 bool enableAdaptive = Parameters::Get<Parameters::EnableAdaptiveTimeStepping>();
284 bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
285 if (enableAdaptive) {
286 const UnitSystem& unitSystem = this->simulator_.vanguard().eclState().getUnits();
287 const auto& sched_state = schedule()[timer.currentStepNum()];
288 auto max_next_tstep = sched_state.max_next_tstep(enableTUNING);
290 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(max_next_tstep,
291 sched_state.tuning(),
292 unitSystem, report_, terminalOutput_);
295 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(unitSystem, report_, max_next_tstep, terminalOutput_);
300 adaptiveTimeStepping_->setSuggestedNextStep(simulator_.timeStepSize());
305 void updateTUNING(
const Tuning& tuning)
307 modelParam_.tolerance_cnv_ = tuning.TRGCNV;
308 modelParam_.tolerance_cnv_relaxed_ = tuning.XXXCNV;
309 modelParam_.tolerance_mb_ = tuning.TRGMBE;
310 modelParam_.tolerance_mb_relaxed_ = tuning.XXXMBE;
311 modelParam_.newton_max_iter_ = tuning.NEWTMX;
312 modelParam_.newton_min_iter_ = tuning.NEWTMN;
313 if (terminalOutput_) {
314 const auto msg = fmt::format(fmt::runtime(
"Tuning values: "
315 "MB: {:.2e}, CNV: {:.2e}, NEWTMN: {}, NEWTMX: {}"),
316 tuning.TRGMBE, tuning.TRGCNV, tuning.NEWTMN, tuning.NEWTMX);
318 if (tuning.TRGTTE_has_value) {
319 OpmLog::warning(
"Tuning item 2-1 (TRGTTE) is not supported.");
321 if (tuning.TRGLCV_has_value) {
322 OpmLog::warning(
"Tuning item 2-4 (TRGLCV) is not supported.");
324 if (tuning.XXXTTE_has_value) {
325 OpmLog::warning(
"Tuning item 2-5 (XXXTTE) is not supported.");
327 if (tuning.XXXLCV_has_value) {
328 OpmLog::warning(
"Tuning item 2-8 (XXXLCV) is not supported.");
330 if (tuning.XXXWFL_has_value) {
331 OpmLog::warning(
"Tuning item 2-9 (XXXWFL) is not supported.");
333 if (tuning.TRGFIP_has_value) {
334 OpmLog::warning(
"Tuning item 2-10 (TRGFIP) is not supported.");
336 if (tuning.TRGSFT_has_value) {
337 OpmLog::warning(
"Tuning item 2-11 (TRGSFT) is not supported.");
339 if (tuning.THIONX_has_value) {
340 OpmLog::warning(
"Tuning item 2-12 (THIONX) is not supported.");
342 if (tuning.TRWGHT_has_value) {
343 OpmLog::warning(
"Tuning item 2-13 (TRWGHT) is not supported.");
345 if (tuning.LITMAX_has_value) {
346 OpmLog::warning(
"Tuning item 3-3 (LITMAX) is not supported.");
348 if (tuning.LITMIN_has_value) {
349 OpmLog::warning(
"Tuning item 3-4 (LITMIN) is not supported.");
351 if (tuning.MXWSIT_has_value) {
352 OpmLog::warning(
"Tuning item 3-5 (MXWSIT) is not supported.");
354 if (tuning.MXWPIT_has_value) {
355 OpmLog::warning(
"Tuning item 3-6 (MXWPIT) is not supported.");
357 if (tuning.DDPLIM_has_value) {
358 OpmLog::warning(
"Tuning item 3-7 (DDPLIM) is not supported.");
360 if (tuning.DDSLIM_has_value) {
361 OpmLog::warning(
"Tuning item 3-8 (DDSLIM) is not supported.");
363 if (tuning.TRGDPR_has_value) {
364 OpmLog::warning(
"Tuning item 3-9 (TRGDPR) is not supported.");
366 if (tuning.XXXDPR_has_value) {
367 OpmLog::warning(
"Tuning item 3-10 (XXXDPR) is not supported.");
369 if (tuning.MNWRFP_has_value) {
370 OpmLog::warning(
"Tuning item 3-11 (MNWRFP) is not supported.");
375 void updateTUNINGDP(
const TuningDp& tuning_dp)
378 modelParam_.tolerance_max_dp_ = tuning_dp.TRGDDP;
379 modelParam_.tolerance_max_ds_ = tuning_dp.TRGDDS;
380 modelParam_.tolerance_max_drs_ = tuning_dp.TRGDDRS;
381 modelParam_.tolerance_max_drv_ = tuning_dp.TRGDDRV;
384 if (terminalOutput_) {
386 if (tuning_dp.TRGLCV_has_value) {
387 OpmLog::warning(
"TUNINGDP item 1 (TRGLCV) is not supported.");
389 if (tuning_dp.XXXLCV_has_value) {
390 OpmLog::warning(
"TUNINGDP item 2 (XXXLCV) is not supported.");
395 bool runStep(SimulatorTimer& timer)
397 if (schedule().exitStatus().has_value()) {
398 if (terminalOutput_) {
399 OpmLog::info(
"Stopping simulation since EXIT was triggered by an action keyword.");
401 report_.success.exit_status = schedule().exitStatus().value();
405 if (serializer_.shouldLoad()) {
406 serializer_.loadTimerInfo(timer);
410 if (terminalOutput_) {
411 std::ostringstream ss;
413 OpmLog::debug(ss.str());
414 details::outputReportStep(timer);
418 if (timer.initialStep()) {
419 Dune::Timer perfTimer;
422 simulator_.setEpisodeIndex(-1);
423 simulator_.setEpisodeLength(0.0);
424 simulator_.setTimeStepSize(0.0);
425 wellModel_().beginReportStep(timer.currentStepNum());
426 simulator_.problem().writeOutput(
true);
428 report_.success.output_write_time += perfTimer.stop();
432 solverTimer_->start();
435 solver_ = createSolver(wellModel_());
438 simulator_.startNextEpisode(
439 simulator_.startTime()
440 + schedule().seconds(timer.currentStepNum()),
441 timer.currentStepLength());
442 simulator_.setEpisodeIndex(timer.currentStepNum());
444 if (serializer_.shouldLoad()) {
445 wellModel_().prepareDeserialize(serializer_.loadStep() - 1);
446 serializer_.loadState();
447 simulator_.model().invalidateAndUpdateIntensiveQuantities(0);
450 this->solver_->model().beginReportStep();
452 const bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
459 if (adaptiveTimeStepping_) {
460 auto tuningUpdater = [enableTUNING,
this,
461 reportStep = timer.currentStepNum()](
const double curr_time,
462 double substep_length,
463 const int sub_step_number)
465 auto& schedule = this->simulator_.vanguard().schedule();
466 auto& events = this->schedule()[reportStep].events();
469 if (events.hasEvent(ScheduleEvents::TUNING_CHANGE)) {
471 schedule.clear_event(ScheduleEvents::TUNING_CHANGE, reportStep);
472 const auto& sched_state = schedule[reportStep];
473 const auto& max_next_tstep = sched_state.max_next_tstep(enableTUNING);
474 const auto& tuning = sched_state.tuning();
477 adaptiveTimeStepping_->updateTUNING(max_next_tstep, tuning);
480 solver_->model().updateTUNING(tuning);
481 this->updateTUNING(tuning);
482 substep_length = this->adaptiveTimeStepping_->suggestedNextStep();
484 substep_length = max_next_tstep;
485 this->adaptiveTimeStepping_->updateNEXTSTEP(max_next_tstep);
487 result = max_next_tstep > 0;
490 if (events.hasEvent(ScheduleEvents::TUNINGDP_CHANGE)) {
492 schedule.clear_event(ScheduleEvents::TUNINGDP_CHANGE, reportStep);
497 const auto& sched_state = schedule[reportStep];
498 const auto& tuning_dp = sched_state.tuning_dp();
499 solver_->model().updateTUNINGDP(tuning_dp);
500 this->updateTUNINGDP(tuning_dp);
503 const auto& wcycle = schedule[reportStep].wcycle.get();
504 if (wcycle.empty()) {
508 const auto& wmatcher = schedule.wellMatcher(reportStep);
509 double wcycle_time_step =
510 wcycle.nextTimeStep(curr_time,
513 this->wellModel_().wellOpenTimes(),
514 this->wellModel_().wellCloseTimes(),
516 &wg_events = this->wellModel_().reportStepStartEvents()]
517 (
const std::string& name)
519 if (sub_step_number != 0) {
522 return wg_events.hasEvent(name, ScheduleEvents::REQUEST_OPEN_WELL);
525 wcycle_time_step = this->grid().comm().min(wcycle_time_step);
526 if (substep_length != wcycle_time_step) {
527 this->adaptiveTimeStepping_->updateNEXTSTEP(wcycle_time_step);
534 tuningUpdater(timer.simulationTimeElapsed(),
535 this->adaptiveTimeStepping_->suggestedNextStep(), 0);
537#ifdef RESERVOIR_COUPLING_ENABLED
538 if (this->reservoirCouplingMaster_) {
539 this->reservoirCouplingMaster_->maybeSpawnSlaveProcesses(timer.currentStepNum());
540 this->reservoirCouplingMaster_->maybeActivate(timer.currentStepNum());
542 else if (this->reservoirCouplingSlave_) {
543 this->reservoirCouplingSlave_->maybeActivate(timer.currentStepNum());
546 const auto& events = schedule()[timer.currentStepNum()].events();
547 bool event = events.hasEvent(ScheduleEvents::NEW_WELL) ||
548 events.hasEvent(ScheduleEvents::INJECTION_TYPE_CHANGED) ||
549 events.hasEvent(ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER) ||
550 events.hasEvent(ScheduleEvents::PRODUCTION_UPDATE) ||
551 events.hasEvent(ScheduleEvents::INJECTION_UPDATE) ||
552 events.hasEvent(ScheduleEvents::WELL_STATUS_CHANGE);
553 auto stepReport = adaptiveTimeStepping_->step(timer, *solver_, event, tuningUpdater);
554 report_ += stepReport;
557 auto stepReport = solver_->step(timer,
nullptr);
558 report_ += stepReport;
560 simulator_.problem().setSubStepReport(stepReport);
561 simulator_.problem().setSimulationReport(report_);
562 simulator_.problem().endTimeStep();
563 if (terminalOutput_) {
564 std::ostringstream ss;
565 stepReport.reportStep(ss);
566 OpmLog::info(ss.str());
571 Dune::Timer perfTimer;
573 const double nextstep = adaptiveTimeStepping_ ? adaptiveTimeStepping_->suggestedNextStep() : -1.0;
574 simulator_.problem().setNextTimeStepSize(nextstep);
575 simulator_.problem().writeOutput(
true);
576 report_.success.output_write_time += perfTimer.stop();
578 solver_->model().endReportStep();
581 solverTimer_->stop();
584 report_.success.solver_time += solverTimer_->secsSinceStart();
586 if (this->grid().comm().rank() == 0) {
589 const auto& reps = this->solver_->model().stepReports();
590 convergence_output_.write(reps);
596 if (terminalOutput_) {
598 "Time step took " + std::to_string(solverTimer_->secsSinceStart()) +
" seconds; "
599 "total solver time " + std::to_string(report_.success.solver_time) +
" seconds.";
603 serializer_.save(timer);
608 SimulatorReport finalize()
612 Dune::Timer finalOutputTimer;
613 finalOutputTimer.start();
615 simulator_.problem().finalizeOutput();
616 report_.success.output_write_time += finalOutputTimer.stop();
621 report_.success.total_time = totalTimer_->secsSinceStart();
622 report_.success.converged =
true;
627 const Grid& grid()
const
628 {
return simulator_.vanguard().grid(); }
630 template<
class Serializer>
631 void serializeOp(Serializer& serializer)
633 serializer(simulator_);
635 serializer(adaptiveTimeStepping_);
638 const Model& model()
const
639 {
return solver_->model(); }
644 [[maybe_unused]]
const std::string& groupName)
override
647 serializer.read(*
this, groupName,
"simulator_data");
653 [[maybe_unused]]
const std::string& groupName)
const override
656 serializer.write(*
this, groupName,
"simulator_data");
663 std::ostringstream str;
668 simulator_.vanguard().caseName(),
675 return simulator_.vanguard().globalCell();
678 std::unique_ptr<Solver> createSolver(WellModel& wellModel)
680 auto model = std::make_unique<Model>(simulator_,
685 if (this->modelParam_.write_partitions_) {
686 const auto& iocfg = this->eclState().cfg().io();
688 const auto odir = iocfg.getOutputDir()
689 / std::filesystem::path {
"partition" }
690 / iocfg.getBaseName();
692 if (this->grid().comm().rank() == 0) {
693 create_directories(odir);
696 this->grid().comm().barrier();
698 model->writePartitions(odir);
700 this->modelParam_.write_partitions_ =
false;
703 return std::make_unique<Solver>(solverParam_, std::move(model));
706 const EclipseState& eclState()
const
707 {
return simulator_.vanguard().eclState(); }
710 const Schedule& schedule()
const
711 {
return simulator_.vanguard().schedule(); }
713 bool isRestart()
const
715 const auto& initconfig = eclState().getInitConfig();
716 return initconfig.restartRequested();
719 WellModel& wellModel_()
720 {
return simulator_.problem().wellModel(); }
722 const WellModel& wellModel_()
const
723 {
return simulator_.problem().wellModel(); }
726 Simulator& simulator_;
728 ModelParameters modelParam_;
729 SolverParameters solverParam_;
731 std::unique_ptr<Solver> solver_;
735 bool terminalOutput_;
737 SimulatorReport report_;
738 std::unique_ptr<time::StopWatch> solverTimer_;
739 std::unique_ptr<time::StopWatch> totalTimer_;
740 std::unique_ptr<TimeStepper> adaptiveTimeStepping_;
742 SimulatorConvergenceOutput convergence_output_{};
744#ifdef RESERVOIR_COUPLING_ENABLED
745 bool slaveMode_{
false};
746 std::unique_ptr<ReservoirCouplingMaster<Scalar>> reservoirCouplingMaster_{
nullptr};
747 std::unique_ptr<ReservoirCouplingSlave<Scalar>> reservoirCouplingSlave_{
nullptr};
750 SimulatorSerializer serializer_;