OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-01-24 01:26:26 -06:00
Code Issues Packages Projects Releases Wiki Activity

collecting the NaN and too large well residuals

Browse Source 
make sure all the processes will throw if there is any of the processes
found abnormal residual values.
This commit is contained in:
Kai Bao 2017-08-22 14:49:30 +02:00
parent 8abe48a693
commit 355be6c26c
5 changed files with 102 additions and 25 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

8
opm/autodiff/StandardWell.hpp
View File

@ -58,6 +58,8 @@ namespace Opm
};
using typename Base::Scalar;
using typename Base::ConvergenceReport;
using Base::numEq;
using Base::has_solvent;
@ -123,9 +125,9 @@ namespace Opm
wellhelpers::WellSwitchingLogger& logger) const;
/// check whether the well equations get converged for this well
virtual bool getWellConvergence(Simulator& ebosSimulator,
const std::vector<double>& B_avg,
const ModelParameters& param) const;
virtual ConvergenceReport getWellConvergence(Simulator& ebosSimulator,
const std::vector<double>& B_avg,
const ModelParameters& param) const;
/// computing the accumulation term for later use in well mass equations
virtual void computeAccumWell();

33
opm/autodiff/StandardWell_impl.hpp
View File

@ -1457,7 +1457,7 @@ namespace Opm
template<typename TypeTag>
bool
typename StandardWell<TypeTag>::ConvergenceReport
StandardWell<TypeTag>::
getWellConvergence(Simulator& ebosSimulator,
const std::vector<double>& B_avg,
@ -1485,30 +1485,43 @@ namespace Opm
}
Vector well_flux_residual(numComp);
bool converged_Well = true;
// Finish computation
for ( int compIdx = 0; compIdx < numComp; ++compIdx )
{
well_flux_residual[compIdx] = B_avg[compIdx] * res[compIdx];
converged_Well = converged_Well && (well_flux_residual[compIdx] < tol_wells);
}
// if one of the residuals is NaN, throw exception, so that the solver can be restarted
ConvergenceReport report;
// checking if any NaN or too large residuals found
// TODO: not understand why phase here while component in other places.
for (int phaseIdx = 0; phaseIdx < np; ++phaseIdx) {
const auto& phaseName = FluidSystem::phaseName(flowPhaseToEbosPhaseIdx(phaseIdx));
if (std::isnan(well_flux_residual[phaseIdx])) {
OPM_THROW(Opm::NumericalProblem, "NaN residual for phase " << phaseName << " for well " << name());
}
if (well_flux_residual[phaseIdx] > maxResidualAllowed) {
OPM_THROW(Opm::NumericalProblem, "Too large residual for phase " << phaseName << " for well " << name());
report.nan_residual_found = true;
const typename ConvergenceReport::ProblemWell problem_well = {name(), phaseName};
report.nan_residual_wells.push_back(problem_well);
} else {
if (well_flux_residual[phaseIdx] > maxResidualAllowed) {
report.too_large_residual_found = true;
const typename ConvergenceReport::ProblemWell problem_well = {name(), phaseName};
report.too_large_residual_wells.push_back(problem_well);
}
}
}
return converged_Well;
if ( !(report.nan_residual_found || report.too_large_residual_found) ) { // no abnormal residual value found
// check convergence
for ( int compIdx = 0; compIdx < numComp; ++compIdx )
{
report.converged = report.converged && (well_flux_residual[compIdx] < tol_wells);
}
} else { // abnormal values found and no need to check the convergence
report.converged = false;
}
return report;
}

2
opm/autodiff/StandardWellsDense.hpp
View File

@ -168,6 +168,8 @@ namespace Opm {
// later, might make share_ptr const later.
std::vector<WellInterfacePtr > well_container_;
using ConvergenceReport = typename WellInterface<TypeTag>::ConvergenceReport;
// create the well container
static std::vector<WellInterfacePtr > createWellContainer(const Wells* wells,
const std::vector<const Well*>& wells_ecl,

39
opm/autodiff/StandardWellsDense_impl.hpp
View File

@ -479,22 +479,51 @@ namespace Opm {
getWellConvergence(Simulator& ebosSimulator,
const std::vector<Scalar>& B_avg) const
{
bool converged_well = true;
ConvergenceReport report;
// currently, if there is any well not converged, we consider the well eqautions do not get converged
for (const auto& well : well_container_) {
if ( !well->getWellConvergence(ebosSimulator, B_avg, param_) ) {
converged_well = false;
report += well->getWellConvergence(ebosSimulator, B_avg, param_);
}
// checking NaN residuals
{
bool nan_residual_found = report.nan_residual_found;
const auto& grid = ebosSimulator.gridManager().grid();
int value = nan_residual_found ? 1 : 0;
nan_residual_found = grid.comm().max(value);
if (nan_residual_found) {
for (const auto& well : report.nan_residual_wells) {
OpmLog::debug("NaN residual found with phase " + well.phase_name + " for well " + well.well_name);
}
OPM_THROW(Opm::NumericalProblem, "NaN residual found!");
}
}
// checking too large residuals
{
bool too_large_residual_found = report.too_large_residual_found;
const auto& grid = ebosSimulator.gridManager().grid();
int value = too_large_residual_found ? 1 : 0;
too_large_residual_found = grid.comm().max(value);
if (too_large_residual_found) {
for (const auto& well : report.too_large_residual_wells) {
OpmLog::debug("Too large residual found with phase " + well.phase_name + " fow well " + well.well_name);
}
OPM_THROW(Opm::NumericalProblem, "Too large residual found!");
}
}
// checking convergence
bool converged_well = report.converged;
{
const auto& grid = ebosSimulator.gridManager().grid();
int value = converged_well ? 1 : 0;
converged_well = grid.comm().min(value);
}
// TODO: to think about the output here.
return converged_well;
}

45
opm/autodiff/WellInterface.hpp
View File

@ -109,9 +109,40 @@ namespace Opm
virtual void initPrimaryVariablesEvaluation() const = 0;
virtual bool getWellConvergence(Simulator& ebosSimulator,
const std::vector<double>& B_avg,
const ModelParameters& param) const = 0;
/// a struct to collect information about the convergence checking
struct ConvergenceReport {
struct ProblemWell {
std::string well_name;
std::string phase_name;
};
bool converged = true;
bool nan_residual_found = false;
std::vector<ProblemWell> nan_residual_wells;
// We consider Inf is large residual here
bool too_large_residual_found = false;
std::vector<ProblemWell> too_large_residual_wells;
ConvergenceReport& operator+=(const ConvergenceReport& rhs) {
converged = converged && rhs.converged;
nan_residual_found = nan_residual_found || rhs.nan_residual_found;
if (rhs.nan_residual_found) {
for (const ProblemWell& well : rhs.nan_residual_wells) {
nan_residual_wells.push_back(well);
}
}
too_large_residual_found = too_large_residual_found || rhs.too_large_residual_found;
if (rhs.too_large_residual_found) {
for (const ProblemWell& well : rhs.too_large_residual_wells) {
too_large_residual_wells.push_back(well);
}
}
return *this;
}
};
virtual ConvergenceReport getWellConvergence(Simulator& ebosSimulator,
const std::vector<double>& B_avg,
const ModelParameters& param) const = 0;
virtual void solveEqAndUpdateWellState(const ModelParameters& param,
WellState& well_state) = 0;
@ -128,15 +159,15 @@ namespace Opm
void computeRepRadiusPerfLength(const Grid& grid, const std::map<int, int>& cartesian_to_compressed);
// using the solution x to recover the solution xw for wells and applying
// xw to update Well State
/// using the solution x to recover the solution xw for wells and applying
/// xw to update Well State
virtual void recoverWellSolutionAndUpdateWellState(const BVector& x, const ModelParameters& param,
WellState& well_state) const = 0;
// Ax = Ax - C D^-1 B x
/// Ax = Ax - C D^-1 B x
virtual void apply(const BVector& x, BVector& Ax) const = 0;
// r = r - C D^-1 Rw
/// r = r - C D^-1 Rw
virtual void apply(BVector& r) const = 0;
virtual void computeWellPotentials(const Simulator& ebosSimulator,