24#ifndef OPM_BLACKOILMODEL_NLDD_HEADER_INCLUDED
25#define OPM_BLACKOILMODEL_NLDD_HEADER_INCLUDED
27#include <dune/common/timer.hh>
29#include <opm/grid/common/SubGridPart.hpp>
52#include <fmt/format.h>
74template<
class TypeTag>
class BlackoilModel;
77template <
class TypeTag>
81 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
82 using Grid = GetPropType<TypeTag, Properties::Grid>;
83 using Indices = GetPropType<TypeTag, Properties::Indices>;
85 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
93 static constexpr int numEq = Indices::numEq;
100 : model_(model), rank_(model_.simulator().vanguard().grid().comm().rank())
103 const auto& [partition_vector, num_domains] = this->partitionCells();
106 std::vector<int> sizes(num_domains, 0);
107 for (
const auto& p : partition_vector) {
112 using EntitySeed =
typename Grid::template Codim<0>::EntitySeed;
113 std::vector<std::vector<EntitySeed>> seeds(num_domains);
114 std::vector<std::vector<int>> partitions(num_domains);
115 for (
int domain = 0; domain < num_domains; ++domain) {
116 seeds[domain].resize(sizes[domain]);
117 partitions[domain].resize(sizes[domain]);
122 const auto& grid = model_.simulator().vanguard().grid();
124 std::vector<int> count(num_domains, 0);
125 const auto& gridView = grid.leafGridView();
126 const auto beg = gridView.template begin<0, Dune::Interior_Partition>();
127 const auto end = gridView.template end<0, Dune::Interior_Partition>();
129 for (
auto it = beg; it != end; ++it, ++cell) {
130 const int p = partition_vector[cell];
131 seeds[p][count[p]] = it->seed();
132 partitions[p][count[p]] = cell;
135 assert(count == sizes);
138 for (
int index = 0; index < num_domains; ++index) {
139 std::vector<bool> interior(partition_vector.size(),
false);
140 for (
int ix : partitions[index]) {
144 Dune::SubGridPart<Grid> view{grid, std::move(seeds[index])};
146 this->domains_.emplace_back(index,
147 std::move(partitions[index]),
153 domain_matrices_.resize(num_domains);
156 for (
int index = 0; index < num_domains; ++index) {
160 const auto& eclState = model_.simulator().vanguard().eclState();
162 loc_param.template init<TypeTag>(eclState.getSimulationConfig().useCPR());
165 if (domains_[index].cells.size() < 200) {
172 const bool force_serial =
true;
173 domain_linsolvers_.emplace_back(model_.simulator(), loc_param, force_serial);
176 assert(
int(domains_.size()) == num_domains);
183 model_.wellModel().setupDomains(domains_);
187 template <
class NonlinearSolverType>
190 NonlinearSolverType& nonlinear_solver)
194 Dune::Timer perfTimer;
196 model_.initialLinearization(report, iteration, nonlinear_solver.minIter(), nonlinear_solver.maxIter(), timer);
204 auto& solution = model_.simulator().model().solution(0);
205 auto initial_solution = solution;
206 auto locally_solved = initial_solution;
209 const auto domain_order = this->getSubdomainOrder();
213 std::vector<SimulatorReportSingle> domain_reports(domains_.size());
214 for (
const int domain_index : domain_order) {
215 const auto& domain = domains_[domain_index];
218 switch (model_.param().local_solve_approach_) {
220 solveDomainJacobi(solution, locally_solved, local_report, logger,
221 iteration, timer, domain);
225 solveDomainGaussSeidel(solution, locally_solved, local_report, logger,
226 iteration, timer, domain);
239 logger.
debug(fmt::format(
"Convergence failure in domain {} on rank {}." , domain.index, rank_));
241 domain_reports[domain.index] = local_report;
246 global_logger.logMessages();
251 std::array<int, 4> counts{ 0, 0, 0,
static_cast<int>(domain_reports.size()) };
252 int& num_converged = counts[0];
253 int& num_converged_already = counts[1];
254 int& num_local_newtons = counts[2];
255 int& num_domains = counts[3];
258 for (
const auto& dr : domain_reports) {
261 if (dr.total_newton_iterations == 0) {
262 ++num_converged_already;
268 local_reports_accumulated_ += rep;
272 solution = locally_solved;
273 model_.simulator().model().invalidateAndUpdateIntensiveQuantities(0);
284 const auto& comm = model_.simulator().vanguard().grid().comm();
285 if (comm.size() > 1) {
286 const auto* ccomm = model_.simulator().model().newtonMethod().linearSolver().comm();
289 ccomm->copyOwnerToAll(solution, solution);
292 const std::size_t num = solution.size();
293 Dune::BlockVector<std::size_t> allmeanings(num);
294 for (std::size_t ii = 0; ii < num; ++ii) {
297 ccomm->copyOwnerToAll(allmeanings, allmeanings);
298 for (std::size_t ii = 0; ii < num; ++ii) {
303 model_.simulator().model().invalidateAndUpdateIntensiveQuantitiesOverlap(0);
306 comm.sum(counts.data(), counts.size());
310 const bool is_iorank = this->rank_ == 0;
312 OpmLog::debug(fmt::format(
"Local solves finished. Converged for {}/{} domains. {} domains did no work. {} total local Newton iterations.\n",
313 num_converged, num_domains, num_converged_already, num_local_newtons));
321 auto rep = model_.nonlinearIterationNewton(iteration + 100, timer, nonlinear_solver);
332 return local_reports_accumulated_;
337 const auto& elementMapper = this->model_.simulator().model().elementMapper();
338 const auto& cartMapper = this->model_.simulator().vanguard().cartesianIndexMapper();
340 const auto& grid = this->model_.simulator().vanguard().grid();
341 const auto& comm = grid.comm();
342 const auto nDigit = 1 +
static_cast<int>(std::floor(std::log10(comm.size())));
344 std::ofstream pfile { odir / fmt::format(
"{1:0>{0}}", nDigit, comm.rank()) };
346 const auto p = this->reconstitutePartitionVector();
348 for (
const auto& cell : elements(grid.leafGridView(), Dune::Partitions::interior)) {
349 pfile << comm.rank() <<
' '
350 << cartMapper.cartesianIndex(elementMapper.index(cell)) <<
' '
357 std::pair<SimulatorReportSingle, ConvergenceReport>
358 solveDomain(
const Domain& domain,
361 [[maybe_unused]]
const int global_iteration,
362 const bool initial_assembly_required)
364 auto& modelSimulator = model_.simulator();
367 Dune::Timer solveTimer;
369 Dune::Timer detailTimer;
371 modelSimulator.model().newtonMethod().setIterationIndex(0);
378 if (initial_assembly_required) {
380 modelSimulator.model().newtonMethod().setIterationIndex(iter);
383 modelSimulator.problem().wellModel().assembleDomain(modelSimulator.model().newtonMethod().numIterations(),
384 modelSimulator.timeStepSize(),
387 report += this->assembleReservoirDomain(domain);
392 std::vector<double> resnorms;
393 auto convreport = this->getDomainConvergence(domain, timer, 0, logger, resnorms);
394 if (convreport.converged()) {
397 return { report, convreport };
404 model_.wellModel().linearizeDomain(domain,
405 modelSimulator.model().linearizer().jacobian(),
406 modelSimulator.model().linearizer().residual());
407 const double tt1 = detailTimer.stop();
412 const int max_iter = model_.param().max_local_solve_iterations_;
413 const auto& grid = modelSimulator.vanguard().grid();
417 const int nc = grid.size(0);
421 this->solveJacobianSystemDomain(domain, x);
422 model_.wellModel().postSolveDomain(x, domain);
430 this->updateDomainSolution(domain, x);
437 modelSimulator.model().newtonMethod().setIterationIndex(iter);
441 modelSimulator.problem().wellModel().assembleDomain(modelSimulator.model().newtonMethod().numIterations(),
442 modelSimulator.timeStepSize(),
444 report += this->assembleReservoirDomain(domain);
450 convreport = this->getDomainConvergence(domain, timer, iter, logger, resnorms);
456 model_.wellModel().linearizeDomain(domain,
457 modelSimulator.model().linearizer().jacobian(),
458 modelSimulator.model().linearizer().residual());
459 const double tt2 = detailTimer.stop();
462 }
while (!convreport.converged() && iter <= max_iter);
464 modelSimulator.problem().endIteration();
466 report.
converged = convreport.converged();
471 return { report, convreport };
475 SimulatorReportSingle assembleReservoirDomain(
const Domain& domain)
478 model_.simulator().model().linearizer().linearizeDomain(domain);
479 return model_.wellModel().lastReport();
483 void solveJacobianSystemDomain(
const Domain& domain,
BVector& global_x)
485 const auto& modelSimulator = model_.simulator();
487 Dune::Timer perfTimer;
490 const Mat& main_matrix = modelSimulator.model().linearizer().jacobian().istlMatrix();
491 if (domain_matrices_[domain.index]) {
494 domain_matrices_[domain.index] = std::make_unique<Mat>(
Details::extractMatrix(main_matrix, domain.cells));
496 auto& jac = *domain_matrices_[domain.index];
505 auto& linsolver = domain_linsolvers_[domain.index];
507 linsolver.prepare(jac, res);
508 model_.linearSolveSetupTime() = perfTimer.stop();
509 linsolver.setResidual(res);
516 void updateDomainSolution(
const Domain& domain,
const BVector& dx)
518 auto& simulator = model_.simulator();
519 auto& newtonMethod = simulator.model().newtonMethod();
522 newtonMethod.update_(solution,
531 simulator.model().invalidateAndUpdateIntensiveQuantities(0, domain);
535 std::pair<double, double> localDomainConvergenceData(
const Domain& domain,
536 std::vector<Scalar>& R_sum,
537 std::vector<Scalar>& maxCoeff,
538 std::vector<Scalar>& B_avg,
539 std::vector<int>& maxCoeffCell)
541 const auto& modelSimulator = model_.simulator();
543 double pvSumLocal = 0.0;
544 double numAquiferPvSumLocal = 0.0;
545 const auto& model = modelSimulator.model();
546 const auto& problem = modelSimulator.problem();
548 const auto& modelResid = modelSimulator.model().linearizer().residual();
551 const auto& gridView = domain.view;
552 const auto& elemEndIt = gridView.template end<0>();
553 IsNumericalAquiferCell isNumericalAquiferCell(gridView.grid());
555 for (
auto elemIt = gridView.template begin</*codim=*/0>();
559 if (elemIt->partitionType() != Dune::InteriorEntity) {
562 const auto& elem = *elemIt;
563 elemCtx.updatePrimaryStencil(elem);
564 elemCtx.updatePrimaryIntensiveQuantities(0);
566 const unsigned cell_idx = elemCtx.globalSpaceIndex(0, 0);
567 const auto& intQuants = elemCtx.intensiveQuantities(0, 0);
568 const auto& fs = intQuants.fluidState();
570 const auto pvValue = problem.referencePorosity(cell_idx, 0) *
571 model.dofTotalVolume(cell_idx);
572 pvSumLocal += pvValue;
574 if (isNumericalAquiferCell(elem))
576 numAquiferPvSumLocal += pvValue;
579 model_.getMaxCoeff(cell_idx, intQuants, fs, modelResid, pvValue,
580 B_avg, R_sum, maxCoeff, maxCoeffCell);
584 const int bSize = B_avg.size();
585 for (
int i = 0; i<bSize; ++i )
587 B_avg[ i ] /=
Scalar(domain.cells.size());
590 return {pvSumLocal, numAquiferPvSumLocal};
593 ConvergenceReport getDomainReservoirConvergence(
const double reportTime,
597 DeferredLogger& logger,
598 std::vector<Scalar>& B_avg,
599 std::vector<Scalar>& residual_norms)
601 using Vector = std::vector<Scalar>;
603 const int numComp =
numEq;
604 Vector R_sum(numComp, 0.0 );
605 Vector maxCoeff(numComp, std::numeric_limits<Scalar>::lowest() );
606 std::vector<int> maxCoeffCell(numComp, -1);
607 const auto [ pvSum, numAquiferPvSum]
608 = this->localDomainConvergenceData(domain, R_sum, maxCoeff, B_avg, maxCoeffCell);
610 auto cnvErrorPvFraction = computeCnvErrorPvLocal(domain, B_avg, dt);
611 cnvErrorPvFraction /= (pvSum - numAquiferPvSum);
616 const bool use_relaxed_cnv = cnvErrorPvFraction < model_.param().relaxed_max_pv_fraction_ &&
617 iteration >= model_.param().min_strict_cnv_iter_;
620 const double tol_cnv = model_.param().local_tolerance_scaling_cnv_
621 * (use_relaxed_cnv ? model_.param().tolerance_cnv_relaxed_
622 : model_.param().tolerance_cnv_);
624 const bool use_relaxed_mb = iteration >= model_.param().min_strict_mb_iter_;
625 const double tol_mb = model_.param().local_tolerance_scaling_mb_
626 * (use_relaxed_mb ? model_.param().tolerance_mb_relaxed_ : model_.param().tolerance_mb_);
629 std::vector<Scalar> CNV(numComp);
630 std::vector<Scalar> mass_balance_residual(numComp);
631 for (
int compIdx = 0; compIdx < numComp; ++compIdx )
633 CNV[compIdx] = B_avg[compIdx] * dt * maxCoeff[compIdx];
634 mass_balance_residual[compIdx] = std::abs(B_avg[compIdx]*R_sum[compIdx]) * dt / pvSum;
635 residual_norms.push_back(CNV[compIdx]);
639 ConvergenceReport report{reportTime};
640 using CR = ConvergenceReport;
641 for (
int compIdx = 0; compIdx < numComp; ++compIdx) {
642 double res[2] = { mass_balance_residual[compIdx], CNV[compIdx] };
643 CR::ReservoirFailure::Type types[2] = { CR::ReservoirFailure::Type::MassBalance,
644 CR::ReservoirFailure::Type::Cnv };
645 double tol[2] = { tol_mb, tol_cnv };
646 for (
int ii : {0, 1}) {
647 if (std::isnan(res[ii])) {
648 report.setReservoirFailed({types[ii], CR::Severity::NotANumber, compIdx});
649 logger.debug(
"NaN residual for " + model_.compNames().name(compIdx) +
" equation.");
650 }
else if (res[ii] > model_.param().max_residual_allowed_) {
651 report.setReservoirFailed({types[ii], CR::Severity::TooLarge, compIdx});
652 logger.debug(
"Too large residual for " + model_.compNames().name(compIdx) +
" equation.");
653 }
else if (res[ii] < 0.0) {
654 report.setReservoirFailed({types[ii], CR::Severity::Normal, compIdx});
655 logger.debug(
"Negative residual for " + model_.compNames().name(compIdx) +
" equation.");
656 }
else if (res[ii] > tol[ii]) {
657 report.setReservoirFailed({types[ii], CR::Severity::Normal, compIdx});
659 report.setReservoirConvergenceMetric(types[ii], compIdx, res[ii]);
664 const bool converged_at_initial_state = (report.
converged() && iteration == 0);
665 if (!converged_at_initial_state) {
666 if (iteration == 0) {
668 std::string msg = fmt::format(
"Domain {} on rank {}, size {}, containing cell {}\n| Iter",
669 domain.index, this->rank_, domain.cells.size(), domain.cells[0]);
670 for (
int compIdx = 0; compIdx < numComp; ++compIdx) {
672 msg += model_.compNames().name(compIdx)[0];
675 for (
int compIdx = 0; compIdx < numComp; ++compIdx) {
677 msg += model_.compNames().name(compIdx)[0];
683 std::ostringstream ss;
685 const std::streamsize oprec = ss.precision(3);
686 const std::ios::fmtflags oflags = ss.setf(std::ios::scientific);
687 ss << std::setw(4) << iteration;
688 for (
int compIdx = 0; compIdx < numComp; ++compIdx) {
689 ss << std::setw(11) << mass_balance_residual[compIdx];
691 for (
int compIdx = 0; compIdx < numComp; ++compIdx) {
692 ss << std::setw(11) << CNV[compIdx];
696 logger.debug(ss.str());
702 ConvergenceReport getDomainConvergence(
const Domain& domain,
703 const SimulatorTimerInterface& timer,
705 DeferredLogger& logger,
706 std::vector<double>& residual_norms)
708 std::vector<Scalar> B_avg(
numEq, 0.0);
709 auto report = this->getDomainReservoirConvergence(timer.simulationTimeElapsed(),
710 timer.currentStepLength(),
716 report += model_.wellModel().getDomainWellConvergence(domain, B_avg, logger);
721 std::vector<int> getSubdomainOrder()
723 const auto& modelSimulator = model_.simulator();
724 const auto& solution = modelSimulator.model().solution(0);
726 std::vector<int> domain_order(domains_.size());
727 std::iota(domain_order.begin(), domain_order.end(), 0);
734 std::vector<double> measure_per_domain(domains_.size());
735 switch (model_.param().local_domain_ordering_) {
738 for (
const auto& domain : domains_) {
739 double press_sum = 0.0;
740 for (
const int c : domain.cells) {
741 press_sum += solution[c][Indices::pressureSwitchIdx];
743 const double avgpress = press_sum / domain.cells.size();
744 measure_per_domain[domain.index] = avgpress;
750 for (
const auto& domain : domains_) {
751 double maxpress = 0.0;
752 for (
const int c : domain.cells) {
753 maxpress = std::max(maxpress, solution[c][Indices::pressureSwitchIdx]);
755 measure_per_domain[domain.index] = maxpress;
761 const auto& residual = modelSimulator.model().linearizer().residual();
762 const int num_vars = residual[0].size();
763 for (
const auto& domain : domains_) {
765 for (
const int c : domain.cells) {
766 for (
int ii = 0; ii < num_vars; ++ii) {
767 maxres = std::max(maxres, std::fabs(residual[c][ii]));
770 measure_per_domain[domain.index] = maxres;
777 const auto& m = measure_per_domain;
778 std::stable_sort(domain_order.begin(), domain_order.end(),
779 [&m](
const int i1,
const int i2){ return m[i1] > m[i2]; });
782 throw std::logic_error(
"Domain solve approach must be Jacobi or Gauss-Seidel");
786 template<
class GlobalEqVector>
787 void solveDomainJacobi(GlobalEqVector& solution,
788 GlobalEqVector& locally_solved,
789 SimulatorReportSingle& local_report,
790 DeferredLogger& logger,
792 const SimulatorTimerInterface& timer,
795 auto initial_local_well_primary_vars = model_.wellModel().getPrimaryVarsDomain(domain);
797 auto res = solveDomain(domain, timer, logger, iteration,
false);
798 local_report = res.first;
799 if (local_report.converged) {
803 model_.simulator().model().invalidateAndUpdateIntensiveQuantities(0, domain);
805 model_.wellModel().setPrimaryVarsDomain(domain, initial_local_well_primary_vars);
807 model_.simulator().model().invalidateAndUpdateIntensiveQuantities(0, domain);
811 template<
class GlobalEqVector>
812 void solveDomainGaussSeidel(GlobalEqVector& solution,
813 GlobalEqVector& locally_solved,
814 SimulatorReportSingle& local_report,
815 DeferredLogger& logger,
817 const SimulatorTimerInterface& timer,
820 auto initial_local_well_primary_vars = model_.wellModel().getPrimaryVarsDomain(domain);
822 auto res = solveDomain(domain, timer, logger, iteration,
true);
823 local_report = res.first;
824 if (!local_report.converged) {
827 const auto& convrep = res.second;
829 if (!convrep.wellFailed()) {
832 double cnv_sum = 0.0;
833 for (
const auto& rc : convrep.reservoirConvergence()) {
835 mb_sum += rc.value();
837 cnv_sum += rc.value();
841 const double acceptable_local_mb_sum = 1e-3;
842 const double acceptable_local_cnv_sum = 1.0;
843 if (mb_sum < acceptable_local_mb_sum && cnv_sum < acceptable_local_cnv_sum) {
844 local_report.converged =
true;
845 logger.debug(fmt::format(
"Accepting solution in unconverged domain {} on rank {}.", domain.index, rank_));
846 logger.debug(fmt::format(
"Value of mb_sum: {} cnv_sum: {}", mb_sum, cnv_sum));
848 logger.debug(
"Unconverged local solution.");
851 logger.debug(
"Unconverged local solution with well convergence failures:");
852 for (
const auto& wf : convrep.wellFailures()) {
857 if (local_report.converged) {
861 model_.wellModel().setPrimaryVarsDomain(domain, initial_local_well_primary_vars);
863 model_.simulator().model().invalidateAndUpdateIntensiveQuantities(0, domain);
867 double computeCnvErrorPvLocal(
const Domain& domain,
868 const std::vector<Scalar>& B_avg,
double dt)
const
871 const auto& simulator = model_.simulator();
872 const auto& model = simulator.model();
873 const auto& problem = simulator.problem();
874 const auto& residual = simulator.model().linearizer().residual();
876 for (
const int cell_idx : domain.cells) {
877 const double pvValue = problem.referencePorosity(cell_idx, 0) *
878 model.dofTotalVolume(cell_idx);
879 const auto& cellResidual = residual[cell_idx];
880 bool cnvViolated =
false;
882 for (
unsigned eqIdx = 0; eqIdx < cellResidual.size(); ++eqIdx) {
884 Scalar CNV = cellResidual[eqIdx] * dt * B_avg[eqIdx] / pvValue;
885 cnvViolated = cnvViolated || (fabs(CNV) > model_.param().tolerance_cnv_);
895 decltype(
auto) partitionCells()
const
897 const auto& grid = this->model_.simulator().vanguard().grid();
899 using GridView = std::remove_cv_t<std::remove_reference_t<
decltype(grid.leafGridView())>>;
900 using Element = std::remove_cv_t<std::remove_reference_t<typename GridView::template Codim<0>::Entity>>;
902 const auto& param = this->model_.param();
904 auto zoltan_ctrl = ZoltanPartitioningControl<Element>{};
905 zoltan_ctrl.domain_imbalance = param.local_domain_partition_imbalance_;
907 [elementMapper = &this->model_.simulator().model().elementMapper()]
908 (
const Element& element)
910 return elementMapper->index(element);
912 zoltan_ctrl.local_to_global =
913 [cartMapper = &this->model_.simulator().vanguard().cartesianIndexMapper()]
916 return cartMapper->cartesianIndex(elemIdx);
920 const auto need_wells = param.local_domain_partition_method_ ==
"zoltan";
921 const auto wells = need_wells
922 ? this->model_.simulator().vanguard().schedule().getWellsatEnd()
923 : std::vector<Well>{};
926 constexpr int default_cells_per_domain = 1000;
928 const int num_domains = param.num_local_domains_ > 0
929 ? param.num_local_domains_
930 : num_cells / default_cells_per_domain;
934 grid.leafGridView(), wells, zoltan_ctrl);
937 std::vector<int> reconstitutePartitionVector()
const
939 const auto& grid = this->model_.simulator().vanguard().grid();
941 auto numD = std::vector<int>(grid.comm().size() + 1, 0);
942 numD[grid.comm().rank() + 1] =
static_cast<int>(this->domains_.size());
943 grid.comm().sum(numD.data(), numD.size());
944 std::partial_sum(numD.begin(), numD.end(), numD.begin());
946 auto p = std::vector<int>(grid.size(0));
947 auto maxCellIdx = std::numeric_limits<int>::min();
949 auto d = numD[grid.comm().rank()];
950 for (
const auto& domain : this->domains_) {
951 for (
const auto& cell : domain.cells) {
953 if (cell > maxCellIdx) {
961 p.erase(p.begin() + maxCellIdx + 1, p.end());
965 BlackoilModel<TypeTag>& model_;
966 std::vector<Domain> domains_;
967 std::vector<std::unique_ptr<Mat>> domain_matrices_;
968 std::vector<ISTLSolverType> domain_linsolvers_;
969 SimulatorReportSingle local_reports_accumulated_;
A NLDD implementation for three-phase black oil.
Definition: BlackoilModelNldd.hpp:78
GetPropType< TypeTag, Properties::SolutionVector > SolutionVector
Definition: BlackoilModelNldd.hpp:86
GetPropType< TypeTag, Properties::Indices > Indices
Definition: BlackoilModelNldd.hpp:83
void writePartitions(const std::filesystem::path &odir) const
Definition: BlackoilModelNldd.hpp:335
BlackoilModelNldd(BlackoilModel< TypeTag > &model)
The constructor sets up the subdomains.
Definition: BlackoilModelNldd.hpp:99
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilModelNldd.hpp:85
typename BlackoilModel< TypeTag >::Mat Mat
Definition: BlackoilModelNldd.hpp:91
SimulatorReportSingle nonlinearIterationNldd(const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &nonlinear_solver)
Do one non-linear NLDD iteration.
Definition: BlackoilModelNldd.hpp:188
SubDomain< Grid > Domain
Definition: BlackoilModelNldd.hpp:89
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: BlackoilModelNldd.hpp:80
const SimulatorReportSingle & localAccumulatedReports() const
return the statistics if the nonlinearIteration() method failed
Definition: BlackoilModelNldd.hpp:330
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: BlackoilModelNldd.hpp:81
void prepareStep()
Called before starting a time step.
Definition: BlackoilModelNldd.hpp:180
static constexpr int numEq
Definition: BlackoilModelNldd.hpp:93
GetPropType< TypeTag, Properties::Grid > Grid
Definition: BlackoilModelNldd.hpp:82
typename BlackoilModel< TypeTag >::BVector BVector
Definition: BlackoilModelNldd.hpp:88
Definition: BlackoilModel.hpp:164
typename SparseMatrixAdapter::IstlMatrix Mat
Definition: BlackoilModel.hpp:210
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilModel.hpp:211
Definition: DeferredLogger.hpp:57
void debug(const std::string &tag, const std::string &message)
Definition: ISTLSolver.hpp:143
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:34
Vector extractVector(const Vector &x, const std::vector< int > &indices)
Definition: extractMatrix.hpp:104
void copySubMatrix(const Matrix &A, const std::vector< int > &indices, Matrix &B)
Definition: extractMatrix.hpp:35
void setGlobal(const Vector &x, const std::vector< int > &indices, Vector &global_x)
Definition: extractMatrix.hpp:115
Matrix extractMatrix(const Matrix &m, const std::vector< int > &indices)
Definition: extractMatrix.hpp:47
void unPack(PV &privar, const std::size_t meanings)
Definition: priVarsPacking.hpp:41
std::size_t pack(const PV &privar)
Definition: priVarsPacking.hpp:31
std::size_t countGlobalCells(const Grid &grid)
Get the number of cells of a global grid.
Definition: countGlobalCells.hpp:82
Definition: BlackoilPhases.hpp:27
Opm::DeferredLogger gatherDeferredLogger(const Opm::DeferredLogger &local_deferredlogger, Parallel::Communication communicator)
Create a global log combining local logs.
std::pair< std::vector< int >, int > partitionCells(const std::string &method, const int num_local_domains, const GridView &grid_view, const std::vector< Well > &wells, const ZoltanPartitioningControl< Element > &zoltan_ctrl)
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:484
This class carries all parameters for the NewtonIterationBlackoilInterleaved class.
Definition: FlowLinearSolverParameters.hpp:237
double linear_solver_reduction_
Definition: FlowLinearSolverParameters.hpp:238
bool linear_solver_print_json_definition_
Definition: FlowLinearSolverParameters.hpp:252
std::string linsolver_
Definition: FlowLinearSolverParameters.hpp:251
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34
double linear_solve_time
Definition: SimulatorReport.hpp:43
double assemble_time
Definition: SimulatorReport.hpp:39
double assemble_time_well
Definition: SimulatorReport.hpp:41
bool converged
Definition: SimulatorReport.hpp:54
double total_time
Definition: SimulatorReport.hpp:37
double linear_solve_setup_time
Definition: SimulatorReport.hpp:42
unsigned int total_newton_iterations
Definition: SimulatorReport.hpp:49
double update_time
Definition: SimulatorReport.hpp:44
unsigned int total_linearizations
Definition: SimulatorReport.hpp:48
unsigned int total_linear_iterations
Definition: SimulatorReport.hpp:50
Definition: SubDomain.hpp:62