93 struct MPI_Comm_Deleter;
113 using ModelParameters =
typename Model::ModelParameters;
114 using SolverParameters =
typename Solver::SolverParameters;
120 : simulator_(simulator)
123 simulator_.vanguard().eclState().getIOConfig(),
124 Parameters::Get<Parameters::SaveStep>(),
125 Parameters::Get<Parameters::LoadStep>(),
126 Parameters::Get<Parameters::SaveFile>(),
127 Parameters::Get<Parameters::LoadFile>())
131 this->terminalOutput_ =
false;
132 if (this->grid().comm().rank() == 0) {
136 [compNames =
typename Model::ComponentName{}](
const int compIdx)
137 {
return std::string_view { compNames.name(compIdx) }; }
140 if (!simulator_.vanguard().eclState().getIOConfig().initOnly()) {
141 this->convergence_output_.
142 startThread(this->simulator_.vanguard().eclState(),
144 R
"(OutputExtraConvergenceInfo (--output-extra-convergence-info))",
156 static void registerParameters()
158 ModelParameters::registerParameters();
159 SolverParameters::registerParameters();
160 TimeStepper::registerParameters();
161 detail::registerSimulatorParameters();
170#ifdef RESERVOIR_COUPLING_ENABLED
171 SimulatorReport
run(SimulatorTimer& timer,
int argc,
char** argv)
173 init(timer, argc, argv);
182 simulator_.model().invalidateAndUpdateIntensiveQuantities(0);
184 while (!timer.
done()) {
185 simulator_.problem().writeReports(timer);
186 bool continue_looping = runStep(timer);
187 if (!continue_looping)
break;
189 simulator_.problem().writeReports(timer);
191#ifdef RESERVOIR_COUPLING_ENABLED
194 if (this->reservoirCouplingMaster_) {
195 this->reservoirCouplingMaster_->sendTerminateAndDisconnect();
197 else if (this->reservoirCouplingSlave_ && !this->reservoirCouplingSlave_->terminated()) {
201 this->reservoirCouplingSlave_->receiveTerminateAndDisconnect();
208#ifdef RESERVOIR_COUPLING_ENABLED
213 bool checkRunningAsReservoirCouplingMaster()
215 for (std::size_t report_step = 0; report_step < this->schedule().size(); ++report_step) {
216 auto rescoup = this->schedule()[report_step].rescoup();
217 auto slave_count = rescoup.slaveCount();
218 auto master_group_count = rescoup.masterGroupCount();
222 if (slave_count > 0) {
225 else if (master_group_count > 0) {
227 throw ReservoirCouplingError(
228 "Inconsistent reservoir coupling master schedule: "
229 "Master group count is greater than 0 but slave count is 0"
237#ifdef RESERVOIR_COUPLING_ENABLED
239 void init(
const SimulatorTimer& timer,
int argc,
char** argv)
241 auto slave_mode = Parameters::Get<Parameters::Slave>();
243 this->reservoirCouplingSlave_ =
244 std::make_unique<ReservoirCouplingSlave<Scalar>>(
246 this->schedule(), timer
248 this->reservoirCouplingSlave_->sendAndReceiveInitialData();
249 this->simulator_.setReservoirCouplingSlave(this->reservoirCouplingSlave_.get());
250 wellModel_().setReservoirCouplingSlave(this->reservoirCouplingSlave_.get());
253 auto master_mode = checkRunningAsReservoirCouplingMaster();
255 this->reservoirCouplingMaster_ =
256 std::make_unique<ReservoirCouplingMaster<Scalar>>(
261 this->simulator_.setReservoirCouplingMaster(this->reservoirCouplingMaster_.get());
262 wellModel_().setReservoirCouplingMaster(this->reservoirCouplingMaster_.get());
266 void init(
const SimulatorTimer& timer)
269 simulator_.setEpisodeIndex(-1);
272 solverTimer_ = std::make_unique<time::StopWatch>();
273 totalTimer_ = std::make_unique<time::StopWatch>();
274 totalTimer_->start();
277 bool enableAdaptive = Parameters::Get<Parameters::EnableAdaptiveTimeStepping>();
278 bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
279 if (enableAdaptive) {
280 const UnitSystem& unitSystem = this->simulator_.vanguard().eclState().getUnits();
281 const auto& sched_state = schedule()[timer.currentStepNum()];
282 auto max_next_tstep = sched_state.max_next_tstep(enableTUNING);
284 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(max_next_tstep,
285 sched_state.tuning(),
286 unitSystem, report_, terminalOutput_);
289 adaptiveTimeStepping_ = std::make_unique<TimeStepper>(unitSystem, report_, max_next_tstep, terminalOutput_);
294 adaptiveTimeStepping_->setSuggestedNextStep(simulator_.timeStepSize());
299 void updateTUNING(
const Tuning& tuning)
301 modelParam_.tolerance_cnv_ = tuning.TRGCNV;
302 modelParam_.tolerance_cnv_relaxed_ = tuning.XXXCNV;
303 modelParam_.tolerance_mb_ = tuning.TRGMBE;
304 modelParam_.tolerance_mb_relaxed_ = tuning.XXXMBE;
305 modelParam_.newton_max_iter_ = tuning.NEWTMX;
306 modelParam_.newton_min_iter_ = tuning.NEWTMN;
307 if (terminalOutput_) {
308 const auto msg = fmt::format(fmt::runtime(
"Tuning values: "
309 "MB: {:.2e}, CNV: {:.2e}, NEWTMN: {}, NEWTMX: {}"),
310 tuning.TRGMBE, tuning.TRGCNV, tuning.NEWTMN, tuning.NEWTMX);
312 if (tuning.TRGTTE_has_value) {
313 OpmLog::warning(
"Tuning item 2-1 (TRGTTE) is not supported.");
315 if (tuning.TRGLCV_has_value) {
316 OpmLog::warning(
"Tuning item 2-4 (TRGLCV) is not supported.");
318 if (tuning.XXXTTE_has_value) {
319 OpmLog::warning(
"Tuning item 2-5 (XXXTTE) is not supported.");
321 if (tuning.XXXLCV_has_value) {
322 OpmLog::warning(
"Tuning item 2-8 (XXXLCV) is not supported.");
324 if (tuning.XXXWFL_has_value) {
325 OpmLog::warning(
"Tuning item 2-9 (XXXWFL) is not supported.");
327 if (tuning.TRGFIP_has_value) {
328 OpmLog::warning(
"Tuning item 2-10 (TRGFIP) is not supported.");
330 if (tuning.TRGSFT_has_value) {
331 OpmLog::warning(
"Tuning item 2-11 (TRGSFT) is not supported.");
333 if (tuning.THIONX_has_value) {
334 OpmLog::warning(
"Tuning item 2-12 (THIONX) is not supported.");
336 if (tuning.TRWGHT_has_value) {
337 OpmLog::warning(
"Tuning item 2-13 (TRWGHT) is not supported.");
339 if (tuning.LITMAX_has_value) {
340 OpmLog::warning(
"Tuning item 3-3 (LITMAX) is not supported.");
342 if (tuning.LITMIN_has_value) {
343 OpmLog::warning(
"Tuning item 3-4 (LITMIN) is not supported.");
345 if (tuning.MXWSIT_has_value) {
346 OpmLog::warning(
"Tuning item 3-5 (MXWSIT) is not supported.");
348 if (tuning.MXWPIT_has_value) {
349 OpmLog::warning(
"Tuning item 3-6 (MXWPIT) is not supported.");
351 if (tuning.DDPLIM_has_value) {
352 OpmLog::warning(
"Tuning item 3-7 (DDPLIM) is not supported.");
354 if (tuning.DDSLIM_has_value) {
355 OpmLog::warning(
"Tuning item 3-8 (DDSLIM) is not supported.");
357 if (tuning.TRGDPR_has_value) {
358 OpmLog::warning(
"Tuning item 3-9 (TRGDPR) is not supported.");
360 if (tuning.XXXDPR_has_value) {
361 OpmLog::warning(
"Tuning item 3-10 (XXXDPR) is not supported.");
363 if (tuning.MNWRFP_has_value) {
364 OpmLog::warning(
"Tuning item 3-11 (MNWRFP) is not supported.");
369 void updateTUNINGDP(
const TuningDp& tuning_dp)
372 modelParam_.tolerance_max_dp_ = tuning_dp.TRGDDP;
373 modelParam_.tolerance_max_ds_ = tuning_dp.TRGDDS;
374 modelParam_.tolerance_max_drs_ = tuning_dp.TRGDDRS;
375 modelParam_.tolerance_max_drv_ = tuning_dp.TRGDDRV;
378 if (terminalOutput_) {
380 if (tuning_dp.TRGLCV_has_value) {
381 OpmLog::warning(
"TUNINGDP item 1 (TRGLCV) is not supported.");
383 if (tuning_dp.XXXLCV_has_value) {
384 OpmLog::warning(
"TUNINGDP item 2 (XXXLCV) is not supported.");
389 bool runStep(SimulatorTimer& timer)
391 if (schedule().exitStatus().has_value()) {
392 if (terminalOutput_) {
393 OpmLog::info(
"Stopping simulation since EXIT was triggered by an action keyword.");
395 report_.success.exit_status = schedule().exitStatus().value();
399 if (serializer_.shouldLoad()) {
400 serializer_.loadTimerInfo(timer);
404 if (terminalOutput_) {
405 std::ostringstream ss;
407 OpmLog::debug(ss.str());
408 details::outputReportStep(timer);
412 if (timer.initialStep()) {
413 Dune::Timer perfTimer;
416 simulator_.setEpisodeIndex(-1);
417 simulator_.setEpisodeLength(0.0);
418 simulator_.setTimeStepSize(0.0);
419 wellModel_().beginReportStep(timer.currentStepNum());
420 simulator_.problem().writeOutput(
true);
422 report_.success.output_write_time += perfTimer.stop();
426 solverTimer_->start();
429 solver_ = createSolver(wellModel_());
432 simulator_.startNextEpisode(
433 simulator_.startTime()
434 + schedule().seconds(timer.currentStepNum()),
435 timer.currentStepLength());
436 simulator_.setEpisodeIndex(timer.currentStepNum());
438 if (serializer_.shouldLoad()) {
439 wellModel_().prepareDeserialize(serializer_.loadStep() - 1);
440 serializer_.loadState();
441 simulator_.model().invalidateAndUpdateIntensiveQuantities(0);
444 this->solver_->model().beginReportStep();
446 const bool enableTUNING = Parameters::Get<Parameters::EnableTuning>();
453 if (adaptiveTimeStepping_) {
454 auto tuningUpdater = [enableTUNING,
this,
455 reportStep = timer.currentStepNum()](
const double curr_time,
456 double dt,
const int timeStep)
458 auto& schedule = this->simulator_.vanguard().schedule();
459 auto& events = this->schedule()[reportStep].events();
462 if (events.hasEvent(ScheduleEvents::TUNING_CHANGE)) {
464 schedule.clear_event(ScheduleEvents::TUNING_CHANGE, reportStep);
465 const auto& sched_state = schedule[reportStep];
466 const auto& max_next_tstep = sched_state.max_next_tstep(enableTUNING);
467 const auto& tuning = sched_state.tuning();
470 adaptiveTimeStepping_->updateTUNING(max_next_tstep, tuning);
473 solver_->model().updateTUNING(tuning);
474 this->updateTUNING(tuning);
475 dt = this->adaptiveTimeStepping_->suggestedNextStep();
478 this->adaptiveTimeStepping_->updateNEXTSTEP(max_next_tstep);
480 result = max_next_tstep > 0;
483 if (events.hasEvent(ScheduleEvents::TUNINGDP_CHANGE)) {
485 schedule.clear_event(ScheduleEvents::TUNINGDP_CHANGE, reportStep);
490 const auto& sched_state = schedule[reportStep];
491 const auto& tuning_dp = sched_state.tuning_dp();
492 solver_->model().updateTUNINGDP(tuning_dp);
493 this->updateTUNINGDP(tuning_dp);
496 const auto& wcycle = schedule[reportStep].wcycle.get();
497 if (wcycle.empty()) {
501 const auto& wmatcher = schedule.wellMatcher(reportStep);
502 double wcycle_time_step =
503 wcycle.nextTimeStep(curr_time,
506 this->wellModel_().wellOpenTimes(),
507 this->wellModel_().wellCloseTimes(),
509 &wg_events = this->wellModel_().reportStepStartEvents()]
510 (
const std::string& name)
515 return wg_events.hasEvent(name, ScheduleEvents::REQUEST_OPEN_WELL);
518 wcycle_time_step = this->grid().comm().min(wcycle_time_step);
519 if (dt != wcycle_time_step) {
520 this->adaptiveTimeStepping_->updateNEXTSTEP(wcycle_time_step);
527 tuningUpdater(timer.simulationTimeElapsed(),
528 this->adaptiveTimeStepping_->suggestedNextStep(), 0);
530#ifdef RESERVOIR_COUPLING_ENABLED
531 if (this->reservoirCouplingMaster_) {
532 this->reservoirCouplingMaster_->maybeSpawnSlaveProcesses(timer.currentStepNum());
533 this->reservoirCouplingMaster_->maybeActivate(timer.currentStepNum());
535 else if (this->reservoirCouplingSlave_) {
536 this->reservoirCouplingSlave_->maybeActivate(timer.currentStepNum());
539 const auto& events = schedule()[timer.currentStepNum()].events();
540 bool event = events.hasEvent(ScheduleEvents::NEW_WELL) ||
541 events.hasEvent(ScheduleEvents::INJECTION_TYPE_CHANGED) ||
542 events.hasEvent(ScheduleEvents::WELL_SWITCHED_INJECTOR_PRODUCER) ||
543 events.hasEvent(ScheduleEvents::PRODUCTION_UPDATE) ||
544 events.hasEvent(ScheduleEvents::INJECTION_UPDATE) ||
545 events.hasEvent(ScheduleEvents::WELL_STATUS_CHANGE);
546 auto stepReport = adaptiveTimeStepping_->step(timer, *solver_, event, tuningUpdater);
547 report_ += stepReport;
550 auto stepReport = solver_->step(timer,
nullptr);
551 report_ += stepReport;
553 simulator_.problem().setSubStepReport(stepReport);
554 simulator_.problem().setSimulationReport(report_);
555 simulator_.problem().endTimeStep();
556 if (terminalOutput_) {
557 std::ostringstream ss;
558 stepReport.reportStep(ss);
559 OpmLog::info(ss.str());
564 Dune::Timer perfTimer;
566 const double nextstep = adaptiveTimeStepping_ ? adaptiveTimeStepping_->suggestedNextStep() : -1.0;
567 simulator_.problem().setNextTimeStepSize(nextstep);
568 simulator_.problem().writeOutput(
true);
569 report_.success.output_write_time += perfTimer.stop();
571 solver_->model().endReportStep();
574 solverTimer_->stop();
577 report_.success.solver_time += solverTimer_->secsSinceStart();
579 if (this->grid().comm().rank() == 0) {
582 const auto& reps = this->solver_->model().stepReports();
583 convergence_output_.write(reps);
589 if (terminalOutput_) {
591 "Time step took " + std::to_string(solverTimer_->secsSinceStart()) +
" seconds; "
592 "total solver time " + std::to_string(report_.success.solver_time) +
" seconds.";
596 serializer_.save(timer);
601 SimulatorReport finalize()
605 Dune::Timer finalOutputTimer;
606 finalOutputTimer.start();
608 simulator_.problem().finalizeOutput();
609 report_.success.output_write_time += finalOutputTimer.stop();
614 report_.success.total_time = totalTimer_->secsSinceStart();
615 report_.success.converged =
true;
620 const Grid& grid()
const
621 {
return simulator_.vanguard().grid(); }
623 template<
class Serializer>
624 void serializeOp(Serializer& serializer)
626 serializer(simulator_);
628 serializer(adaptiveTimeStepping_);
631 const Model& model()
const
632 {
return solver_->model(); }
637 [[maybe_unused]]
const std::string& groupName)
override
640 serializer.read(*
this, groupName,
"simulator_data");
646 [[maybe_unused]]
const std::string& groupName)
const override
649 serializer.write(*
this, groupName,
"simulator_data");
656 std::ostringstream str;
661 simulator_.vanguard().caseName(),
668 return simulator_.vanguard().globalCell();
671 std::unique_ptr<Solver> createSolver(WellModel& wellModel)
673 auto model = std::make_unique<Model>(simulator_,
678 if (this->modelParam_.write_partitions_) {
679 const auto& iocfg = this->eclState().cfg().io();
681 const auto odir = iocfg.getOutputDir()
682 / std::filesystem::path {
"partition" }
683 / iocfg.getBaseName();
685 if (this->grid().comm().rank() == 0) {
686 create_directories(odir);
689 this->grid().comm().barrier();
691 model->writePartitions(odir);
693 this->modelParam_.write_partitions_ =
false;
696 return std::make_unique<Solver>(solverParam_, std::move(model));
699 const EclipseState& eclState()
const
700 {
return simulator_.vanguard().eclState(); }
703 const Schedule& schedule()
const
704 {
return simulator_.vanguard().schedule(); }
706 bool isRestart()
const
708 const auto& initconfig = eclState().getInitConfig();
709 return initconfig.restartRequested();
712 WellModel& wellModel_()
713 {
return simulator_.problem().wellModel(); }
715 const WellModel& wellModel_()
const
716 {
return simulator_.problem().wellModel(); }
719 Simulator& simulator_;
721 ModelParameters modelParam_;
722 SolverParameters solverParam_;
724 std::unique_ptr<Solver> solver_;
728 bool terminalOutput_;
730 SimulatorReport report_;
731 std::unique_ptr<time::StopWatch> solverTimer_;
732 std::unique_ptr<time::StopWatch> totalTimer_;
733 std::unique_ptr<TimeStepper> adaptiveTimeStepping_;
735 SimulatorConvergenceOutput convergence_output_{};
737#ifdef RESERVOIR_COUPLING_ENABLED
738 bool slaveMode_{
false};
739 std::unique_ptr<ReservoirCouplingMaster<Scalar>> reservoirCouplingMaster_{
nullptr};
740 std::unique_ptr<ReservoirCouplingSlave<Scalar>> reservoirCouplingSlave_{
nullptr};
743 SimulatorSerializer serializer_;