Merge pull request #734 from qilicun/format-step-output

Format step output
This commit is contained in:
Atgeirr Flø Rasmussen
2016-06-30 09:19:54 +02:00
committed by GitHub
14 changed files with 129 additions and 67 deletions

View File

@@ -82,6 +82,7 @@ namespace Opm {
bool failed;
bool converged;
int linear_iterations;
int well_iterations;
};
@@ -196,9 +197,11 @@ namespace Opm {
/// \param[in] reservoir_state reservoir state variables
/// \param[in, out] well_state well state variables
/// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep
void assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly);
/// \return well iterations.
IterationReport
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly);
/// \brief Compute the residual norms of the mass balance for each phase,
/// the well flux, and the well equation.
@@ -373,7 +376,7 @@ namespace Opm {
assembleMassBalanceEq(const SolutionState& state);
bool
IterationReport
solveWellEq(const std::vector<ADB>& mob_perfcells,
const std::vector<ADB>& b_perfcells,
SolutionState& state,

View File

@@ -288,7 +288,7 @@ namespace detail {
current_relaxation_ = 1.0;
dx_old_ = V::Zero(sizeNonLinear());
}
asImpl().assemble(reservoir_state, well_state, iteration == 0);
IterationReport iter_report = asImpl().assemble(reservoir_state, well_state, iteration == 0);
residual_norms_history_.push_back(asImpl().computeResidualNorms());
const bool converged = asImpl().getConvergence(dt, iteration);
const bool must_solve = (iteration < nonlinear_solver.minIter()) || (!converged);
@@ -322,7 +322,7 @@ namespace detail {
}
const bool failed = false; // Not needed in this model.
const int linear_iters = must_solve ? asImpl().linearIterationsLastSolve() : 0;
return IterationReport{ failed, converged, linear_iters };
return IterationReport{ failed, converged, linear_iters , iter_report.well_iterations};
}
@@ -727,7 +727,7 @@ namespace detail {
template <class Grid, class WellModel, class Implementation>
void
IterationReport
BlackoilModelBase<Grid, WellModel, Implementation>::
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
@@ -775,9 +775,9 @@ namespace detail {
asImpl().assembleMassBalanceEq(state);
// -------- Well equations ----------
IterationReport iter_report = {false, false, 0, std::numeric_limits<int>::min()};
if ( ! wellsActive() ) {
return;
return iter_report;
}
std::vector<ADB> mob_perfcells;
@@ -785,7 +785,7 @@ namespace detail {
asImpl().wellModel().extractWellPerfProperties(state, rq_, mob_perfcells, b_perfcells);
if (param_.solve_welleq_initially_ && initial_assembly) {
// solve the well equations as a pre-processing step
asImpl().solveWellEq(mob_perfcells, b_perfcells, state, well_state);
iter_report = asImpl().solveWellEq(mob_perfcells, b_perfcells, state, well_state);
}
V aliveWells;
std::vector<ADB> cq_s;
@@ -800,7 +800,7 @@ namespace detail {
asImpl().makeConstantState(state0);
asImpl().wellModel().computeWellPotentials(mob_perfcells, b_perfcells, state0, well_state);
}
return iter_report;
}
@@ -971,7 +971,7 @@ namespace detail {
template <class Grid, class WellModel, class Implementation>
bool
IterationReport
BlackoilModelBase<Grid, WellModel, Implementation>::
solveWellEq(const std::vector<ADB>& mob_perfcells,
const std::vector<ADB>& b_perfcells,
@@ -1041,9 +1041,7 @@ namespace detail {
} while (it < 15);
if (converged) {
if ( terminal_output_ ) {
OpmLog::info("well converged iter: " + std::to_string(it));
}
OpmLog::note("well converged iter: " + std::to_string(it));
const int nw = wells().number_of_wells;
{
// We will set the bhp primary variable to the new ones,
@@ -1069,8 +1067,9 @@ namespace detail {
if (!converged) {
well_state = well_state0;
}
return converged;
const bool failed = false; // Not needed in this method.
const int linear_iters = 0; // Not needed in this method
return IterationReport{failed, converged, linear_iters, it};
}
@@ -1863,7 +1862,7 @@ namespace detail {
msg += " W-FLUX(" + materialName(idx).substr(0, 1) + ")";
}
// std::cout << " WELL-CONT ";
OpmLog::info(msg);
OpmLog::note(msg);
}
std::ostringstream ss;
const std::streamsize oprec = ss.precision(3);
@@ -1881,7 +1880,7 @@ namespace detail {
// std::cout << std::setw(11) << residualWell;
ss.precision(oprec);
ss.flags(oflags);
OpmLog::info(ss.str());
OpmLog::note(ss.str());
}
for (int idx = 0; idx < nm; ++idx) {
@@ -1969,7 +1968,7 @@ namespace detail {
for (int idx = 0; idx < np; ++idx) {
msg += " W-FLUX(" + materialName(idx).substr(0, 1) + ")";
}
OpmLog::info(msg);
OpmLog::note(msg);
}
std::ostringstream ss;
const std::streamsize oprec = ss.precision(3);
@@ -1980,7 +1979,7 @@ namespace detail {
}
ss.precision(oprec);
ss.flags(oflags);
OpmLog::info(ss.str());
OpmLog::note(ss.str());
}
return converged;
}

View File

@@ -105,9 +105,10 @@ namespace Opm {
/// \param[in] reservoir_state reservoir state variables
/// \param[in, out] well_state well state variables
/// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep
void assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly);
IterationReport
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly);
using Base::numPhases;
using Base::numMaterials;
@@ -164,7 +165,7 @@ namespace Opm {
const MultisegmentWells::MultisegmentWellOps& msWellOps() const { return well_model_.wellOps(); }
bool
IterationReport
solveWellEq(const std::vector<ADB>& mob_perfcells,
const std::vector<ADB>& b_perfcells,
SolutionState& state,

View File

@@ -126,7 +126,7 @@ namespace Opm {
template <class Grid>
void
IterationReport
BlackoilMultiSegmentModel<Grid>::
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
@@ -182,9 +182,9 @@ namespace Opm {
asImpl().assembleMassBalanceEq(state);
// -------- Well equations ----------
IterationReport iter_report = {false, false, 0, std::numeric_limits<int>::min()};
if ( ! wellsActive() ) {
return;
return iter_report;
}
wellModel().computeSegmentFluidProperties(state);
@@ -197,7 +197,7 @@ namespace Opm {
wellModel().extractWellPerfProperties(state, rq_, mob_perfcells, b_perfcells);
if (param_.solve_welleq_initially_ && initial_assembly) {
// solve the well equations as a pre-processing step
asImpl().solveWellEq(mob_perfcells, b_perfcells, state, well_state);
iter_report = asImpl().solveWellEq(mob_perfcells, b_perfcells, state, well_state);
}
// the perforation flux here are different
@@ -209,6 +209,7 @@ namespace Opm {
wellModel().addWellFluxEq(cq_s, state, residual_);
asImpl().addWellContributionToMassBalanceEq(cq_s, state, well_state);
wellModel().addWellControlEq(state, well_state, aliveWells, residual_);
return iter_report;
}
@@ -216,14 +217,15 @@ namespace Opm {
template <class Grid>
bool BlackoilMultiSegmentModel<Grid>::solveWellEq(const std::vector<ADB>& mob_perfcells,
const std::vector<ADB>& b_perfcells,
SolutionState& state,
WellState& well_state)
IterationReport
BlackoilMultiSegmentModel<Grid>::solveWellEq(const std::vector<ADB>& mob_perfcells,
const std::vector<ADB>& b_perfcells,
SolutionState& state,
WellState& well_state)
{
const bool converged = Base::solveWellEq(mob_perfcells, b_perfcells, state, well_state);
IterationReport iter_report = Base::solveWellEq(mob_perfcells, b_perfcells, state, well_state);
if (converged) {
if (iter_report.converged) {
// We must now update the state.segp and state.segqs members,
// that the base version does not know about.
const int np = numPhases();
@@ -252,7 +254,7 @@ namespace Opm {
asImpl().computeWellConnectionPressures(state, well_state);
}
return converged;
return iter_report;
}

View File

@@ -147,12 +147,12 @@ namespace Opm {
void assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly)
IterationReport
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly)
{
Base::assemble(reservoir_state, well_state, initial_assembly);
IterationReport iter_report = Base::assemble(reservoir_state, well_state, initial_assembly);
if (initial_assembly) {
}
@@ -184,6 +184,7 @@ namespace Opm {
assert(int(well_state.perfRates().size()) == wflux.size());
std::copy_n(wflux.data(), wflux.size(), well_state.perfRates().begin());
}
return iter_report;
}

View File

@@ -150,7 +150,7 @@ namespace Opm {
}
// Report and return.
return IterationReport { false, true, pressure_liniter + transport_liniter };
return IterationReport { false, true, pressure_liniter + transport_liniter, 0 };
} else {
// Iterate to fully implicit solution.
// This call is just for a single iteration (one pressure and one transport solve),
@@ -179,7 +179,7 @@ namespace Opm {
// Report and return.
const bool converged = iteration >= 3; // TODO: replace this with a proper convergence check
return IterationReport { false, converged, pressure_liniter + transport_liniter };
return IterationReport { false, converged, pressure_liniter + transport_liniter, 0 };
}
}

View File

@@ -82,10 +82,10 @@ namespace Opm {
asImpl().makeConstantState(state0_);
}
void assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly)
IterationReport
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly)
{
using namespace Opm::AutoDiffGrid;
@@ -124,9 +124,9 @@ namespace Opm {
asImpl().assembleMassBalanceEq(state);
// -------- Well equations ----------
IterationReport iter_report = {false, false, 0, std::numeric_limits<int>::min()};
if ( ! wellsActive() ) {
return;
return iter_report;
}
std::vector<ADB> mob_perfcells;
@@ -134,7 +134,7 @@ namespace Opm {
asImpl().wellModel().extractWellPerfProperties(state, rq_, mob_perfcells, b_perfcells);
if (param_.solve_welleq_initially_ && initial_assembly) {
// solve the well equations as a pre-processing step
asImpl().solveWellEq(mob_perfcells, b_perfcells, state, well_state);
iter_report = asImpl().solveWellEq(mob_perfcells, b_perfcells, state, well_state);
}
V aliveWells;
std::vector<ADB> cq_s;
@@ -153,6 +153,7 @@ namespace Opm {
asImpl().makeConstantState(state0);
asImpl().wellModel().computeWellPotentials(mob_perfcells, b_perfcells, state0, well_state);
}
return iter_report;
}

View File

@@ -108,12 +108,18 @@ namespace Opm {
/// Number of linear solver iterations used in all calls to step().
unsigned int linearIterations() const;
/// Number of well iterations used in all calls to step().
unsigned int wellIterations() const;
/// Number of nonlinear solver iterations used in the last call to step().
unsigned int nonlinearIterationsLastStep() const;
/// Number of linear solver iterations used in the last call to step().
unsigned int linearIterationsLastStep() const;
/// Number of well iterations used in the last call to step().
unsigned int wellIterationsLastStep() const;
/// Reference to physical model.
const PhysicalModel& model() const;
@@ -151,8 +157,10 @@ namespace Opm {
std::unique_ptr<PhysicalModel> model_;
unsigned int nonlinearIterations_;
unsigned int linearIterations_;
unsigned int wellIterations_;
unsigned int nonlinearIterationsLast_;
unsigned int linearIterationsLast_;
unsigned int wellIterationsLast_;
};
} // namespace Opm

View File

@@ -34,8 +34,10 @@ namespace Opm
model_(std::move(model_arg)),
nonlinearIterations_(0),
linearIterations_(0),
wellIterations_(0),
nonlinearIterationsLast_(0),
linearIterationsLast_(0)
linearIterationsLast_(0),
wellIterationsLast_(0)
{
if (!model_) {
OPM_THROW(std::logic_error, "Must provide a non-null model argument for NonlinearSolver.");
@@ -54,6 +56,12 @@ namespace Opm
return linearIterations_;
}
template <class PhysicalModel>
unsigned int NonlinearSolver<PhysicalModel>::wellIterations() const
{
return wellIterations_;
}
template <class PhysicalModel>
const PhysicalModel& NonlinearSolver<PhysicalModel>::model() const
{
@@ -78,6 +86,12 @@ namespace Opm
return linearIterationsLast_;
}
template <class PhysicalModel>
unsigned int NonlinearSolver<PhysicalModel>::wellIterationsLastStep() const
{
return wellIterationsLast_;
}
template <class PhysicalModel>
int
@@ -110,6 +124,7 @@ namespace Opm
// Set up for main solver loop.
int linIters = 0;
bool converged = false;
int wellIters = 0;
// ---------- Main nonlinear solver loop ----------
do {
@@ -122,6 +137,7 @@ namespace Opm
}
converged = report.converged;
linIters += report.linear_iterations;
wellIters += report.well_iterations;
++iteration;
} while ( (!converged && (iteration <= maxIter())) || (iteration <= minIter()));
@@ -134,8 +150,10 @@ namespace Opm
linearIterations_ += linIters;
nonlinearIterations_ += iteration - 1; // Since the last one will always be trivial.
wellIterations_ += wellIters;
linearIterationsLast_ = linIters;
nonlinearIterationsLast_ = iteration;
wellIterationsLast_ = wellIters;
// Do model-specific post-step actions.
model_->afterStep(dt, reservoir_state, well_state);

View File

@@ -20,7 +20,7 @@
*/
#include <algorithm>
#include <locale>
#include <opm/parser/eclipse/EclipseState/Schedule/Events.hpp>
#include <opm/core/utility/initHydroCarbonState.hpp>
#include <opm/core/well_controls.h>
@@ -137,7 +137,7 @@ namespace Opm
{
std::ostringstream ss;
timer.report(ss);
OpmLog::info(ss.str());
OpmLog::note(ss.str());
}
// Create wells and well state.
@@ -176,6 +176,15 @@ namespace Opm
const WellModel well_model(wells);
auto solver = asImpl().createSolver(well_model);
std::ostringstream step_msg;
boost::posix_time::time_facet* facet = new boost::posix_time::time_facet("%d-%b-%Y");
step_msg.imbue(std::locale(std::locale::classic(), facet));
step_msg << "\nTime step " << std::setw(4) <<timer.currentStepNum()
<< " at day " << (double)unit::convert::to(timer.simulationTimeElapsed(), unit::day)
<< "/" << (double)unit::convert::to(timer.totalTime(), unit::day)
<< ", date = " << timer.currentDateTime()
<< "\n";
OpmLog::info(step_msg.str());
// If sub stepping is enabled allow the solver to sub cycle
// in case the report steps are too large for the solver to converge
@@ -188,6 +197,15 @@ namespace Opm
else {
// solve for complete report step
solver->step(timer.currentStepLength(), state, well_state);
std::ostringstream iter_msg;
iter_msg << "Stepsize " << (double)unit::convert::to(timer.currentStepLength(), unit::day);
if (solver->wellIterations() != std::numeric_limits<int>::min()) {
iter_msg << " days well iterations = " << solver->wellIterations() << ", ";
}
iter_msg << "non-linear iterations = " << solver->nonlinearIterations()
<< ", total linear iterations = " << solver->linearIterations()
<< "\n";
OpmLog::info(iter_msg.str());
}
// update the derived geology (transmissibilities, pore volumes, etc) if the
@@ -221,7 +239,7 @@ namespace Opm
{
std::string msg;
msg = "Fully implicit solver took: " + std::to_string(st) + " seconds. Total solver time taken: " + std::to_string(stime) + " seconds.";
OpmLog::info(msg);
OpmLog::note(msg);
}
if ( output_writer_.output() ) {

View File

@@ -120,9 +120,10 @@ namespace Opm {
/// \param[in] reservoir_state reservoir state variables
/// \param[in, out] well_state well state variables
/// \param[in] initial_assembly pass true if this is the first call to assemble() in this timestep
void assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly);
IterationReport
assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly);
protected:

View File

@@ -489,7 +489,7 @@ namespace Opm {
template <class Grid>
void
IterationReport
BlackoilPolymerModel<Grid>::assemble(const ReservoirState& reservoir_state,
WellState& well_state,
const bool initial_assembly)
@@ -526,10 +526,10 @@ namespace Opm {
// -------- Mass balance equations --------
assembleMassBalanceEq(state);
IterationReport iter_report = {false, false, 0, std::numeric_limits<int>::min()};
// -------- Well equations ----------
if ( ! wellsActive() ) {
return;
return iter_report;
}
std::vector<ADB> mob_perfcells;
@@ -567,6 +567,7 @@ namespace Opm {
wellModel().addWellFluxEq(cq_s, state, residual_);
addWellContributionToMassBalanceEq(cq_s, state, well_state);
wellModel().addWellControlEq(state, well_state, aliveWells, residual_);
return iter_report;
}

View File

@@ -219,6 +219,8 @@ namespace {
const double r0 = residualNorm();
const double r_polymer = residual_.material_balance_eq[2].value().matrix().lpNorm<Eigen::Infinity>();
int it = 0;
// this solver does not solve well equations first.
wellIterations_ = 0;
std::cout << "\nIteration Residual Polymer Res\n"
<< std::setw(9) << it << std::setprecision(9)
<< std::setw(18) << r0 << std::setprecision(9)
@@ -242,7 +244,7 @@ namespace {
newtonIterations_ += 1;
std::cout << std::setw(9) << it << std::setprecision(9)
<< std::setw(18) << r << std::setprecision(9)
<< std::setw(18) << rr_polymer << std::endl;
<< std::setw(18) << rr_polymer << std::endl;
}
if (resTooLarge) {
@@ -267,6 +269,11 @@ namespace {
return linearIterations_;
}
int FullyImplicitCompressiblePolymerSolver::wellIterations() const
{
return wellIterations_;
}
FullyImplicitCompressiblePolymerSolver::ReservoirResidualQuant::ReservoirResidualQuant()
: accum(2, ADB::null())
, mflux( ADB::null())

View File

@@ -64,12 +64,12 @@ namespace Opm {
/// \param[in] wells well structure
/// \param[in] linsolver linear solver
FullyImplicitCompressiblePolymerSolver(const UnstructuredGrid& grid ,
const BlackoilPropsAdInterface& fluid,
const DerivedGeology& geo ,
const RockCompressibility* rock_comp_props,
const PolymerPropsAd& polymer_props_ad,
const Wells& wells,
const NewtonIterationBlackoilInterface& linsolver);
const BlackoilPropsAdInterface& fluid,
const DerivedGeology& geo ,
const RockCompressibility* rock_comp_props,
const PolymerPropsAd& polymer_props_ad,
const Wells& wells,
const NewtonIterationBlackoilInterface& linsolver);
/// Take a single forward step, modifiying
/// state.pressure()
@@ -88,6 +88,7 @@ namespace Opm {
int nonlinearIterations() const;
int linearIterations() const;
int wellIterations() const;
/// Not used by this class except to satisfy interface requirements.
typedef parameter::ParameterGroup SolverParameters;
@@ -159,6 +160,7 @@ namespace Opm {
unsigned int newtonIterations_;
unsigned int linearIterations_;
unsigned int wellIterations_;
// Private methods.
SolutionState