Merge remote-tracking branch 'upstream/master' into master-refactor-for-cpgrid-support

Removed conflicts in
	opm/core/wells/WellsManager.cpp
that were due to the change
```diff
-                    pd.well_index = WellsManagerDetail::computeWellIndex(radius, cubical, cell_perm, completion->getDiameter());
+                    pd.well_index = WellsManagerDetail::computeWellIndex(radius, cubical, cell_perm, completion->getSkinFactor());
```
in  WellsManager::createWellsFromSpecs which moved from WellsManager.cpp to WellsManager_impl.hpp file in a previous commit.
This commit is contained in:
Markus Blatt 2014-04-04 21:21:22 +02:00
commit d83374e6b7
9 changed files with 309 additions and 129 deletions

View File

@ -175,6 +175,7 @@ list (APPEND TEST_SOURCE_FILES
tests/test_wellcontrols.cpp tests/test_wellcontrols.cpp
tests/test_wellsgroup.cpp tests/test_wellsgroup.cpp
tests/test_wellcollection.cpp tests/test_wellcollection.cpp
tests/test_timer.cpp
) )
# originally generated with the command: # originally generated with the command:
@ -190,6 +191,7 @@ list (APPEND TEST_DATA_FILES
tests/wells_manager_data_expanded.data tests/wells_manager_data_expanded.data
tests/wells_manager_data_wellSTOP.data tests/wells_manager_data_wellSTOP.data
tests/wells_group.data tests/wells_group.data
tests/TESTTIMER.DATA
) )
# originally generated with the command: # originally generated with the command:

View File

@ -1178,6 +1178,22 @@ struct EclipseWellBhp : public EclipseWellReport {
"Pascal") "Pascal")
{ } { }
EclipseWellBhp (const EclipseSummary& summary,
Opm::DeckConstPtr newParserDeck,
int whichWell,
PhaseUsage uses,
BlackoilPhases::PhaseIndex phase,
WellType type)
: EclipseWellReport (summary,
newParserDeck,
whichWell,
uses,
phase,
type,
'B',
"Pascal")
{ }
virtual double update (const SimulatorTimer& /*timer*/, virtual double update (const SimulatorTimer& /*timer*/,
const WellState& wellState) const WellState& wellState)
{ {
@ -1306,6 +1322,27 @@ EclipseSummary::addWells (Opm::DeckConstPtr newParserDeck,
} }
} }
} }
// Add BHP monitors
for (int whichWell = 0; whichWell != numWells; ++whichWell) {
// In the call below: uses, phase and the well type arguments
// are not used, except to set up an index that stores the
// well indirectly. For details see the implementation of the
// EclipseWellReport constructor, and the method
// EclipseWellReport::bhp().
BlackoilPhases::PhaseIndex phase = BlackoilPhases::Liquid;
if (!uses.phase_used[BlackoilPhases::Liquid]) {
phase = BlackoilPhases::Vapour;
}
add (std::unique_ptr <EclipseWellReport> (
new EclipseWellBhp (*this,
newParserDeck,
whichWell,
uses,
phase,
WELL_TYPES[0])));
}
} }
namespace Opm { namespace Opm {

View File

@ -39,6 +39,7 @@
#include <dune/istl/solvers.hh> #include <dune/istl/solvers.hh>
#include <dune/istl/paamg/amg.hh> #include <dune/istl/paamg/amg.hh>
#include <dune/istl/paamg/kamg.hh> #include <dune/istl/paamg/kamg.hh>
#include <dune/istl/paamg/pinfo.hh>
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3) #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3)
#include <dune/istl/paamg/fastamg.hh> #include <dune/istl/paamg/fastamg.hh>
@ -59,25 +60,30 @@ namespace Opm
typedef Dune::BlockVector<VectorBlockType> Vector; typedef Dune::BlockVector<VectorBlockType> Vector;
typedef Dune::MatrixAdapter<Mat,Vector,Vector> Operator; typedef Dune::MatrixAdapter<Mat,Vector,Vector> Operator;
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveCG_ILU0(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity); solveCG_ILU0(O& A, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity);
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveCG_AMG(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity, solveCG_AMG(O& A, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity,
double prolongateFactor, int smoothsteps); double prolongateFactor, int smoothsteps);
#if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3) #if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3)
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveKAMG(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity, solveKAMG(O& A, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity,
double prolongateFactor, int smoothsteps); double prolongateFactor, int smoothsteps);
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveFastAMG(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity, solveFastAMG(O& A, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity,
double prolongateFactor); double prolongateFactor);
#endif #endif
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveBiCGStab_ILU0(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity); solveBiCGStab_ILU0(O& A, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity);
} // anonymous namespace } // anonymous namespace
@ -144,17 +150,33 @@ namespace Opm
A[ri][ja[i]] = sa[i]; A[ri][ja[i]] = sa[i];
} }
} }
int maxit = linsolver_max_iterations_;
if (maxit == 0) {
maxit = 5000;
}
Dune::SeqScalarProduct<Vector> sp;
Dune::Amg::SequentialInformation comm;
Operator opA(A);
return solveSystem(opA, solution, rhs, sp, comm, maxit);
}
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport
LinearSolverIstl::solveSystem (O& opA, double* solution, const double* rhs,
S& sp, const C& comm, int maxit) const
{
// System RHS // System RHS
Vector b(size); Vector b(opA.getmat().N());
std::copy(rhs, rhs + size, b.begin()); std::copy(rhs, rhs+b.size(), b.begin());
// System solution // System solution
Vector x(size); Vector x(opA.getmat().M());
x = 0.0; x = 0.0;
if (linsolver_save_system_) if (linsolver_save_system_)
{ {
// Save system to files. // Save system to files.
writeMatrixToMatlab(A, linsolver_save_filename_ + "-mat"); writeMatrixToMatlab(opA.getmat(), linsolver_save_filename_ + "-mat");
std::string rhsfile(linsolver_save_filename_ + "-rhs"); std::string rhsfile(linsolver_save_filename_ + "-rhs");
std::ofstream rhsf(rhsfile.c_str()); std::ofstream rhsf(rhsfile.c_str());
rhsf.precision(15); rhsf.precision(15);
@ -163,23 +185,18 @@ namespace Opm
std::ostream_iterator<VectorBlockType>(rhsf, "\n")); std::ostream_iterator<VectorBlockType>(rhsf, "\n"));
} }
int maxit = linsolver_max_iterations_;
if (maxit == 0) {
maxit = 5000;
}
LinearSolverReport res; LinearSolverReport res;
switch (linsolver_type_) { switch (linsolver_type_) {
case CG_ILU0: case CG_ILU0:
res = solveCG_ILU0(A, x, b, linsolver_residual_tolerance_, maxit, linsolver_verbosity_); res = solveCG_ILU0(opA, x, b, sp, comm, linsolver_residual_tolerance_, maxit, linsolver_verbosity_);
break; break;
case CG_AMG: case CG_AMG:
res = solveCG_AMG(A, x, b, linsolver_residual_tolerance_, maxit, linsolver_verbosity_, res = solveCG_AMG(opA, x, b, sp, comm, linsolver_residual_tolerance_, maxit, linsolver_verbosity_,
linsolver_prolongate_factor_, linsolver_smooth_steps_); linsolver_prolongate_factor_, linsolver_smooth_steps_);
break; break;
case KAMG: case KAMG:
#if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3) #if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3)
res = solveKAMG(A, x, b, linsolver_residual_tolerance_, maxit, linsolver_verbosity_, res = solveKAMG(opA, x, b, sp, comm, linsolver_residual_tolerance_, maxit, linsolver_verbosity_,
linsolver_prolongate_factor_, linsolver_smooth_steps_); linsolver_prolongate_factor_, linsolver_smooth_steps_);
#else #else
throw std::runtime_error("KAMG not supported with this version of DUNE"); throw std::runtime_error("KAMG not supported with this version of DUNE");
@ -187,17 +204,17 @@ namespace Opm
break; break;
case FastAMG: case FastAMG:
#if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3) #if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3)
res = solveFastAMG(A, x, b, linsolver_residual_tolerance_, maxit, linsolver_verbosity_, res = solveFastAMG(opA, x, b, sp, comm, linsolver_residual_tolerance_, maxit, linsolver_verbosity_,
linsolver_prolongate_factor_); linsolver_prolongate_factor_);
#else #else
if(linsolver_verbosity_) if(linsolver_verbosity_)
std::cerr<<"Fast AMG is not available; falling back to CG preconditioned with the normal one"<<std::endl; std::cerr<<"Fast AMG is not available; falling back to CG preconditioned with the normal one"<<std::endl;
res = solveCG_AMG(A, x, b, linsolver_residual_tolerance_, maxit, linsolver_verbosity_, res = solveCG_AMG(opA, x, b, sp, comm, linsolver_residual_tolerance_, maxit, linsolver_verbosity_,
linsolver_prolongate_factor_, linsolver_smooth_steps_); linsolver_prolongate_factor_, linsolver_smooth_steps_);
#endif #endif
break; break;
case BiCGStab_ILU0: case BiCGStab_ILU0:
res = solveBiCGStab_ILU0(A, x, b, linsolver_residual_tolerance_, maxit, linsolver_verbosity_); res = solveBiCGStab_ILU0(opA, x, b, sp, comm, linsolver_residual_tolerance_, maxit, linsolver_verbosity_);
break; break;
default: default:
std::cerr << "Unknown linsolver_type: " << int(linsolver_type_) << '\n'; std::cerr << "Unknown linsolver_type: " << int(linsolver_type_) << '\n';
@ -221,17 +238,37 @@ namespace Opm
namespace namespace
{ {
template<class P, class O>
LinearSolverInterface::LinearSolverReport struct PreconditionerTraits
solveCG_ILU0(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity) {
typedef std::shared_ptr<P> PointerType;
};
template<class P, class O, class C>
typename PreconditionerTraits<P,O>::PointerType
makePreconditioner(O& opA, double relax, const C& comm, int iterations=1)
{
typename Dune::Amg::SmootherTraits<P>::Arguments args;
typename Dune::Amg::ConstructionTraits<P>::Arguments cargs;
cargs.setMatrix(opA.getmat());
args.iterations=iterations;
args.relaxationFactor=relax;
cargs.setArgs(args);
cargs.setComm(comm);
return std::shared_ptr<P>(Dune::Amg::ConstructionTraits<P>::construct(cargs));
}
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport
solveCG_ILU0(O& opA, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity)
{ {
Operator opA(A);
// Construct preconditioner. // Construct preconditioner.
Dune::SeqILU0<Mat,Vector,Vector> precond(A, 1.0); typedef Dune::SeqILU0<Mat,Vector,Vector> Preconditioner;
auto precond = makePreconditioner<Preconditioner>(opA, 1.0, comm);
// Construct linear solver. // Construct linear solver.
Dune::CGSolver<Vector> linsolve(opA, precond, tolerance, maxit, verbosity); Dune::CGSolver<Vector> linsolve(opA, sp, *precond, tolerance, maxit, verbosity);
// Solve system. // Solve system.
Dune::InverseOperatorResult result; Dune::InverseOperatorResult result;
@ -266,8 +303,9 @@ namespace Opm
criterion.setGamma(1); // V-cycle; this is the default criterion.setGamma(1); // V-cycle; this is the default
} }
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveCG_AMG(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity, solveCG_AMG(O& opA, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity,
double linsolver_prolongate_factor, int linsolver_smooth_steps) double linsolver_prolongate_factor, int linsolver_smooth_steps)
{ {
// Solve with AMG solver. // Solve with AMG solver.
@ -295,13 +333,12 @@ namespace Opm
// Construct preconditioner. // Construct preconditioner.
Criterion criterion; Criterion criterion;
Precond::SmootherArgs smootherArgs; Precond::SmootherArgs smootherArgs;
Operator opA(A);
setUpCriterion(criterion, linsolver_prolongate_factor, verbosity, setUpCriterion(criterion, linsolver_prolongate_factor, verbosity,
linsolver_smooth_steps); linsolver_smooth_steps);
Precond precond(opA, criterion, smootherArgs); Precond precond(opA, criterion, smootherArgs);
// Construct linear solver. // Construct linear solver.
Dune::CGSolver<Vector> linsolve(opA, precond, tolerance, maxit, verbosity); Dune::CGSolver<Vector> linsolve(opA, sp, precond, tolerance, maxit, verbosity);
// Solve system. // Solve system.
Dune::InverseOperatorResult result; Dune::InverseOperatorResult result;
@ -317,8 +354,9 @@ namespace Opm
#if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3) #if defined(HAS_DUNE_FAST_AMG) || DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3)
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveKAMG(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity, solveKAMG(O& opA, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity,
double linsolver_prolongate_factor, int linsolver_smooth_steps) double linsolver_prolongate_factor, int linsolver_smooth_steps)
{ {
// Solve with AMG solver. // Solve with AMG solver.
@ -344,7 +382,6 @@ namespace Opm
typedef Dune::Amg::KAMG<Operator,Vector,Smoother,Dune::Amg::SequentialInformation> Precond; typedef Dune::Amg::KAMG<Operator,Vector,Smoother,Dune::Amg::SequentialInformation> Precond;
// Construct preconditioner. // Construct preconditioner.
Operator opA(A);
Precond::SmootherArgs smootherArgs; Precond::SmootherArgs smootherArgs;
Criterion criterion; Criterion criterion;
setUpCriterion(criterion, linsolver_prolongate_factor, verbosity, setUpCriterion(criterion, linsolver_prolongate_factor, verbosity,
@ -352,7 +389,7 @@ namespace Opm
Precond precond(opA, criterion, smootherArgs); Precond precond(opA, criterion, smootherArgs);
// Construct linear solver. // Construct linear solver.
Dune::GeneralizedPCGSolver<Vector> linsolve(opA, precond, tolerance, maxit, verbosity); Dune::GeneralizedPCGSolver<Vector> linsolve(opA, sp, precond, tolerance, maxit, verbosity);
// Solve system. // Solve system.
Dune::InverseOperatorResult result; Dune::InverseOperatorResult result;
@ -366,8 +403,9 @@ namespace Opm
return res; return res;
} }
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveFastAMG(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity, solveFastAMG(O& opA, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity,
double linsolver_prolongate_factor) double linsolver_prolongate_factor)
{ {
// Solve with AMG solver. // Solve with AMG solver.
@ -388,7 +426,6 @@ namespace Opm
typedef Dune::Amg::FastAMG<Operator,Vector> Precond; typedef Dune::Amg::FastAMG<Operator,Vector> Precond;
// Construct preconditioner. // Construct preconditioner.
Operator opA(A);
Criterion criterion; Criterion criterion;
const int smooth_steps = 1; const int smooth_steps = 1;
setUpCriterion(criterion, linsolver_prolongate_factor, verbosity, smooth_steps); setUpCriterion(criterion, linsolver_prolongate_factor, verbosity, smooth_steps);
@ -400,7 +437,7 @@ namespace Opm
Precond precond(opA, criterion, parms); Precond precond(opA, criterion, parms);
// Construct linear solver. // Construct linear solver.
Dune::GeneralizedPCGSolver<Vector> linsolve(opA, precond, tolerance, maxit, verbosity); Dune::GeneralizedPCGSolver<Vector> linsolve(opA, sp, precond, tolerance, maxit, verbosity);
// Solve system. // Solve system.
Dune::InverseOperatorResult result; Dune::InverseOperatorResult result;
@ -415,17 +452,17 @@ namespace Opm
} }
#endif #endif
template<class O, class S, class C>
LinearSolverInterface::LinearSolverReport LinearSolverInterface::LinearSolverReport
solveBiCGStab_ILU0(const Mat& A, Vector& x, Vector& b, double tolerance, int maxit, int verbosity) solveBiCGStab_ILU0(O& opA, Vector& x, Vector& b, S& sp, const C& comm, double tolerance, int maxit, int verbosity)
{ {
Operator opA(A);
// Construct preconditioner. // Construct preconditioner.
Dune::SeqILU0<Mat,Vector,Vector> precond(A, 1.0); typedef Dune::SeqILU0<Mat,Vector,Vector> Preconditioner;
auto precond = makePreconditioner<Preconditioner>(opA, 1.0, comm);
// Construct linear solver. // Construct linear solver.
Dune::BiCGSTABSolver<Vector> linsolve(opA, precond, tolerance, maxit, verbosity); Dune::BiCGSTABSolver<Vector> linsolve(opA, sp, *precond, tolerance, maxit, verbosity);
// Solve system. // Solve system.
Dune::InverseOperatorResult result; Dune::InverseOperatorResult result;

View File

@ -86,6 +86,17 @@ namespace Opm
virtual double getTolerance() const; virtual double getTolerance() const;
private: private:
/// \brief Solve the linear system using ISTL
/// \param[in] opA The linear operator of the system to solve.
/// \param[out] solution C array for storing the solution vector.
/// \param[in] rhs C array containing the right hand side.
/// \param[in] sp The scalar product to use.
/// \param[in] comm The information about the parallel domain decomposition.
/// \param[in] maxit The maximum number of iterations allowed.
template<class O, class S, class C>
LinearSolverReport solveSystem(O& opA, double* solution, const double *rhs,
S& sp, const C& comm, int maxit) const;
double linsolver_residual_tolerance_; double linsolver_residual_tolerance_;
int linsolver_verbosity_; int linsolver_verbosity_;
enum LinsolverType { CG_ILU0 = 0, CG_AMG = 1, BiCGStab_ILU0 = 2, FastAMG=3, KAMG=4 }; enum LinsolverType { CG_ILU0 = 0, CG_AMG = 1, BiCGStab_ILU0 = 2, FastAMG=3, KAMG=4 };

View File

@ -58,63 +58,34 @@ namespace Opm
} }
/// Use the SimulatorTimer as a shim around opm-parser's Opm::TimeMap /// Use the SimulatorTimer as a shim around opm-parser's Opm::TimeMap
void SimulatorTimer::init(Opm::TimeMapConstPtr timeMap, void SimulatorTimer::init(Opm::TimeMapConstPtr timeMap)
size_t beginReportStepIdx,
size_t endReportStepIdx)
{ {
timeMap_ = timeMap;
current_step_ = 0; current_step_ = 0;
beginReportStepIdx_ = beginReportStepIdx; total_time_ = timeMap->getTotalTime();
endReportStepIdx_ = std::min(timeMap_->numTimesteps() + 1, endReportStepIdx); timesteps_.resize(timeMap->numTimesteps());
for ( size_t i = 0; i < timeMap->numTimesteps(); ++i ) {
timesteps_[i] = timeMap->getTimeStepLength(i);
}
boost::posix_time::ptime start_time = timeMap->getStartTime(0);
start_date_ = start_time.date();
} }
/// Total number of steps. /// Total number of steps.
int SimulatorTimer::numSteps() const int SimulatorTimer::numSteps() const
{ {
if (timeMap_)
return endReportStepIdx_ - beginReportStepIdx_;
else
return timesteps_.size(); return timesteps_.size();
} }
/// Index of the first considered simulation episode
size_t SimulatorTimer::beginReportStepIndex() const
{
if (!timeMap_) {
OPM_THROW(std::runtime_error, "indexFirstEpisode() is only implemented "
"for simulation timers which are based on Opm::TimeMap");
}
return beginReportStepIdx_;
}
/// Index of the last considered simulation episode
size_t SimulatorTimer::endReportStepIndex() const
{
if (!timeMap_) {
OPM_THROW(std::runtime_error, "indexLastEpisode() is only implemented "
"for simulation timers which are based on Opm::TimeMap");
}
return endReportStepIdx_;
}
/// Current step number. /// Current step number.
int SimulatorTimer::currentStepNum() const int SimulatorTimer::currentStepNum() const
{ {
return current_step_ + beginReportStepIdx_; return current_step_;
} }
/// Set current step number. /// Set current step number.
void SimulatorTimer::setCurrentStepNum(int step) void SimulatorTimer::setCurrentStepNum(int step)
{ {
if (current_step_ < 0 || current_step_ > int(numSteps())) {
// Note that we do allow current_step_ == timesteps_.size(),
// that is the done() state.
OPM_THROW(std::runtime_error, "Trying to set invalid step number: " << step);
}
current_step_ = step; current_step_ = step;
if (timeMap_)
current_time_ = std::accumulate(timesteps_.begin(), timesteps_.begin() + step, 0.0); current_time_ = std::accumulate(timesteps_.begin(), timesteps_.begin() + step, 0.0);
} }
@ -123,28 +94,18 @@ namespace Opm
double SimulatorTimer::currentStepLength() const double SimulatorTimer::currentStepLength() const
{ {
assert(!done()); assert(!done());
if (timeMap_)
return timeMap_->getTimeStepLength(beginReportStepIdx_ + current_step_);
else
return timesteps_[current_step_]; return timesteps_[current_step_];
} }
double SimulatorTimer::stepLengthTaken() const double SimulatorTimer::stepLengthTaken() const
{ {
assert(current_step_ > 0); assert(current_step_ > 0);
if (timeMap_)
return timeMap_->getTimeStepLength(beginReportStepIdx_ + current_step_ - 1);
else
return timesteps_[current_step_ - 1]; return timesteps_[current_step_ - 1];
} }
/// time elapsed since the start of the simulation [s]. /// time elapsed since the start of the simulation [s].
double SimulatorTimer::simulationTimeElapsed() const double SimulatorTimer::simulationTimeElapsed() const
{ {
if (timeMap_)
return
timeMap_->getTimePassedUntil(beginReportStepIdx_ + current_step_);
else
return current_time_; return current_time_;
} }
@ -157,9 +118,6 @@ namespace Opm
boost::posix_time::ptime SimulatorTimer::currentDateTime() const boost::posix_time::ptime SimulatorTimer::currentDateTime() const
{ {
if (timeMap_)
return timeMap_->getStartTime(beginReportStepIdx_ + current_step_);
else
return boost::posix_time::ptime(start_date_) + boost::posix_time::seconds( (int) current_time_ ); return boost::posix_time::ptime(start_date_) + boost::posix_time::seconds( (int) current_time_ );
} }
@ -168,10 +126,6 @@ namespace Opm
/// Total time. /// Total time.
double SimulatorTimer::totalTime() const double SimulatorTimer::totalTime() const
{ {
if (timeMap_)
return
timeMap_->getTotalTime();
else
return total_time_; return total_time_;
} }
@ -181,12 +135,6 @@ namespace Opm
/// access to later timesteps. /// access to later timesteps.
void SimulatorTimer::setTotalTime(double time) void SimulatorTimer::setTotalTime(double time)
{ {
if (timeMap_) {
// well, what can we do if we use opm-parser's TimeMap?
OPM_THROW(std::logic_error,
"Not implemented: SimulatorTimer::setTotalTime() if using a TimeMap.");
}
else
total_time_ = time; total_time_ = time;
} }
@ -204,7 +152,6 @@ namespace Opm
SimulatorTimer& SimulatorTimer::operator++() SimulatorTimer& SimulatorTimer::operator++()
{ {
assert(!done()); assert(!done());
if (!timeMap_)
current_time_ += timesteps_[current_step_]; current_time_ += timesteps_[current_step_];
++current_step_; ++current_step_;
return *this; return *this;
@ -213,9 +160,6 @@ namespace Opm
/// Return true if op++() has been called numSteps() times. /// Return true if op++() has been called numSteps() times.
bool SimulatorTimer::done() const bool SimulatorTimer::done() const
{ {
if (timeMap_)
return current_step_ > int(endReportStepIdx_ - beginReportStepIdx_ - 1);
else
return int(timesteps_.size()) == current_step_; return int(timesteps_.size()) == current_step_;
} }

View File

@ -50,19 +50,11 @@ namespace Opm
void init(const EclipseGridParser& deck); void init(const EclipseGridParser& deck);
/// Use the SimulatorTimer as a shim around opm-parser's Opm::TimeMap /// Use the SimulatorTimer as a shim around opm-parser's Opm::TimeMap
void init(TimeMapConstPtr timeMap, void init(TimeMapConstPtr timeMap);
size_t beginReportStepIdx = 0,
size_t endReportStepIdx = std::numeric_limits<size_t>::max());
/// Total number of steps. /// Total number of steps.
int numSteps() const; int numSteps() const;
/// Index of the first report step considered
size_t beginReportStepIndex() const;
/// Index of the next-after-last report step to be considered
size_t endReportStepIndex() const;
/// Current step number. This is the number of timesteps that /// Current step number. This is the number of timesteps that
/// has been completed from the start of the run. The time /// has been completed from the start of the run. The time
/// after initialization but before the simulation has started /// after initialization but before the simulation has started
@ -117,10 +109,7 @@ namespace Opm
bool done() const; bool done() const;
private: private:
Opm::TimeMapConstPtr timeMap_;
std::vector<double> timesteps_; std::vector<double> timesteps_;
size_t beginReportStepIdx_;
size_t endReportStepIdx_;
int current_step_; int current_step_;
double current_time_; double current_time_;
double total_time_; double total_time_;

View File

@ -144,7 +144,7 @@ void WellsManager::createWellsFromSpecs(std::vector<WellConstPtr>& wells, size_t
std::array<double, 3> cubical = WellsManagerDetail::getCubeDim(c2f, begin_face_centroids, std::array<double, 3> cubical = WellsManagerDetail::getCubeDim(c2f, begin_face_centroids,
dimensions, cell); dimensions, cell);
const double* cell_perm = &permeability[dimensions*dimensions*cell]; const double* cell_perm = &permeability[dimensions*dimensions*cell];
pd.well_index = WellsManagerDetail::computeWellIndex(radius, cubical, cell_perm, completion->getDiameter()); pd.well_index = WellsManagerDetail::computeWellIndex(radius, cubical, cell_perm, completion->getSkinFactor());
} }
wellperf_data[well_index].push_back(pd); wellperf_data[well_index].push_back(pd);
} }

49
tests/TESTTIMER.DATA Normal file
View File

@ -0,0 +1,49 @@
RUNSPEC
START
26 'MAR' 2014 /
TSTEP
1.0 2*5.0
/
TSTEP
7*10.0 14*25.0
/
TSTEP
19.0
/
TSTEP
18*13.0
/
TSTEP
17*10.0
/
TSTEP
13.0
/
TSTEP
18.0
/
TSTEP
11.0
/
TSTEP
17*5.0
/
TSTEP
19*6.0
/
TSTEP
21*5.0 /
END ==

111
tests/test_timer.cpp Normal file
View File

@ -0,0 +1,111 @@
/*
Copyright 2012 SINTEF ICT, Applied Mathematics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <config.h>
#if HAVE_DYNAMIC_BOOST_TEST
#define BOOST_TEST_DYN_LINK
#endif
#define NVERBOSE // Suppress own messages when throw()ing
#define BOOST_TEST_MODULE OPM-TimerTest
#include <boost/test/unit_test.hpp>
#include <opm/parser/eclipse/Parser/Parser.hpp>
#include <opm/core/simulator/SimulatorTimer.hpp>
#include <opm/core/utility/Units.hpp>
#include <string>
#include <iostream>
#include <vector>
#include <memory>
BOOST_AUTO_TEST_CASE(CreateTimer)
{
const std::string filename1 = "TESTTIMER.DATA";
Opm::ParserPtr parser(new Opm::Parser() );
Opm::DeckConstPtr parserDeck = parser->parseFile( filename1 );
Opm::TimeMapPtr timeMap(new Opm::TimeMap(parserDeck));
Opm::SimulatorTimer simtimer;
boost::gregorian::date defaultStartDate( 2012, 1, 1);
BOOST_CHECK_EQUAL( boost::posix_time::ptime(defaultStartDate), simtimer.currentDateTime() );
simtimer.init(timeMap);
boost::gregorian::date startDate( 2014, 3, 26);
BOOST_CHECK_EQUAL( boost::posix_time::ptime(startDate), simtimer.currentDateTime() );
BOOST_CHECK_EQUAL( 0, simtimer.currentStepNum() );
BOOST_CHECK_EQUAL( 0., simtimer.simulationTimeElapsed() );
BOOST_CHECK_EQUAL( 120, simtimer.numSteps() );
BOOST_CHECK_EQUAL( 1200., Opm::unit::convert::to(simtimer.totalTime(), Opm::unit::day) );
BOOST_CHECK_EQUAL( 0., Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::day) );
double testCurrentTime = 0.;
BOOST_CHECK_EQUAL( Opm::unit::convert::to(testCurrentTime, Opm::unit::day),
Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::day) );
for ( int i = 0; i < simtimer.numSteps(); ++i ) {
BOOST_CHECK_EQUAL( i, simtimer.currentStepNum() );
BOOST_CHECK_EQUAL( Opm::unit::convert::to(testCurrentTime, Opm::unit::minute),
Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::minute) );
testCurrentTime += simtimer.currentStepLength();
++simtimer;
}
for ( int i = 0; i <= simtimer.numSteps(); ++i ) {
simtimer.setCurrentStepNum(i);
BOOST_CHECK_EQUAL( i, simtimer.currentStepNum() );
}
BOOST_CHECK_EQUAL( true, simtimer.done() );
simtimer.setCurrentStepNum(0);
BOOST_CHECK_EQUAL( false, simtimer.done() );
BOOST_CHECK_EQUAL( 0., Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::day) );
simtimer.setCurrentStepNum(120);
BOOST_CHECK_EQUAL( Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::day),
Opm::unit::convert::to(simtimer.totalTime(), Opm::unit::day) );
int i = 0;
double testCurrentTime1 = 0.;
double testCurrentTime2 = 0.;
simtimer.setCurrentStepNum(0);
while (!simtimer.done()) {
testCurrentTime1 += simtimer.currentStepLength();
BOOST_CHECK_EQUAL( i, simtimer.currentStepNum() );
++i;
++simtimer;
testCurrentTime2 += simtimer.stepLengthTaken();
BOOST_CHECK_EQUAL( Opm::unit::convert::to(testCurrentTime1, Opm::unit::minute),
Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::minute) );
BOOST_CHECK_EQUAL( Opm::unit::convert::to(testCurrentTime2, Opm::unit::minute),
Opm::unit::convert::to(simtimer.simulationTimeElapsed(), Opm::unit::minute) );
}
BOOST_CHECK_EQUAL( true, simtimer.done() );
BOOST_CHECK_EQUAL( Opm::unit::convert::to(testCurrentTime1, Opm::unit::minute),
Opm::unit::convert::to(simtimer.totalTime(), Opm::unit::minute) );
BOOST_CHECK_EQUAL( Opm::unit::convert::to(testCurrentTime2, Opm::unit::minute),
Opm::unit::convert::to(simtimer.totalTime(), Opm::unit::minute) );
}