OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #5373 from akva2/various_helpers_template_scalar

Browse Source 
Various helpers: template Scalar type
This commit is contained in:
Bård Skaflestad
2024-05-22 10:50:30 +02:00
committed by GitHub
parent 41a5e748d0 e0bcb314ea
commit 1d92aa952d
44 changed files with 570 additions and 469 deletions
Download Patch File Download Diff File Expand all files Collapse all files

6
opm/simulators/wells/BlackoilWellModelGeneric.cpp
View File

@@ -274,11 +274,11 @@ getLocalWells(const int timeStepIdx) const
}
template<class Scalar>
std::vector<std::reference_wrapper<ParallelWellInfo>>
std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>
BlackoilWellModelGeneric<Scalar>::
createLocalParallelWellInfo(const std::vector<Well>& wells)
{
std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info;
std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>> local_parallel_well_info;
local_parallel_well_info.reserve(wells.size());
for (const auto& well : wells)
{
@@ -318,7 +318,7 @@ initializeWellPerfData()
int connection_index = 0;
// INVALID_ECL_INDEX marks no above perf available
int connection_index_above = ParallelWellInfo::INVALID_ECL_INDEX;
int connection_index_above = ParallelWellInfo<Scalar>::INVALID_ECL_INDEX;
well_perf_data_[well_index].clear();
well_perf_data_[well_index].reserve(well.getConnections().size());

12
opm/simulators/wells/BlackoilWellModelGeneric.hpp
View File

@@ -58,7 +58,7 @@ namespace Opm {
template<class Scalar> class GasLiftWellState;
class Group;
class GuideRateConfig;
class ParallelWellInfo;
template<class Scalar> class ParallelWellInfo;
class RestartValue;
class Schedule;
struct SimulatorUpdate;
@@ -193,7 +193,7 @@ public:
bool forceShutWellByName(const std::string& wellname,
const double simulation_time);
const std::vector<PerforationData>& perfData(const int well_idx) const
const std::vector<PerforationData<Scalar>>& perfData(const int well_idx) const
{ return well_perf_data_[well_idx]; }
const Parallel::Communication& comm() const { return comm_; }
@@ -323,7 +323,7 @@ protected:
/// \brief Create the parallel well information
/// \param localWells The local wells from ECL schedule
std::vector<std::reference_wrapper<ParallelWellInfo>>
std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>
createLocalParallelWellInfo(const std::vector<Well>& wells);
void initializeWellProdIndCalculators();
@@ -450,7 +450,7 @@ protected:
std::optional<int> last_run_wellpi_{};
std::vector<Well> wells_ecl_;
std::vector<std::vector<PerforationData>> well_perf_data_;
std::vector<std::vector<PerforationData<Scalar>>> well_perf_data_;
/// Connection index mappings
class ConnectionIndexMap
@@ -546,8 +546,8 @@ protected:
std::vector<int> local_shut_wells_{};
std::vector<ParallelWellInfo> parallel_well_info_;
std::vector<std::reference_wrapper<ParallelWellInfo>> local_parallel_well_info_;
std::vector<ParallelWellInfo<Scalar>> parallel_well_info_;
std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>> local_parallel_well_info_;
std::vector<WellProdIndexCalculator> prod_index_calc_;
mutable ParallelWBPCalculation wbpCalculationService_;

4
opm/simulators/wells/BlackoilWellModelRestart.cpp
View File

@@ -71,7 +71,7 @@ template<class Scalar>
void BlackoilWellModelRestart<Scalar>::
loadRestartConnectionData(const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
const std::vector<PerforationData<Scalar>>& old_perf_data,
SingleWellState<Scalar>& ws) const
{
auto& perf_data = ws.perf_data;
@@ -130,7 +130,7 @@ loadRestartWellData(const std::string& well_name,
const bool handle_ms_well,
const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
const std::vector<PerforationData<Scalar>>& old_perf_data,
SingleWellState<Scalar>& ws) const
{
const auto np = phs.size();

6
opm/simulators/wells/BlackoilWellModelRestart.hpp
View File

@@ -38,7 +38,7 @@ class GroupAndNetworkValues;
template<class Scalar> class GroupState;
class GuideRate;
class GuideRateConfig;
struct PerforationData;
template<class Scalar> struct PerforationData;
struct PhaseUsage;
template<class Scalar> class SingleWellState;
template<class Scalar> class WellState;
@@ -76,7 +76,7 @@ private:
//! \brief Loads per-connection data from restart structures.
void loadRestartConnectionData(const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
const std::vector<PerforationData<Scalar>>& old_perf_data,
SingleWellState<Scalar>& ws) const;
//! \brief Loads per-segment data from restart structures.
@@ -90,7 +90,7 @@ private:
const bool handle_ms_well,
const std::vector<data::Rates::opt>& phs,
const data::Well& rst_well,
const std::vector<PerforationData>& old_perf_data,
const std::vector<PerforationData<Scalar>>& old_perf_data,
SingleWellState<Scalar>& ws) const;
//! \brief Loads per-group data from restart structures.

89
opm/simulators/wells/MSWellHelpers.cpp
View File

@@ -47,10 +47,10 @@
namespace {
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType haalandFormular(const ValueType& re,
const double diameter,
const double roughness)
const Scalar diameter,
const Scalar roughness)
{
const ValueType value = -3.6 * log10(6.9 / re + std::pow(roughness / (3.7 * diameter), 10. / 9.) );
@@ -62,10 +62,10 @@ ValueType haalandFormular(const ValueType& re,
// water in oil emulsion viscosity
// TODO: maybe it should be two different ValueTypes. When we calculate the viscosity for transitional zone
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType WIOEmulsionViscosity(const ValueType& oil_viscosity,
const ValueType& water_liquid_fraction,
const double max_visco_ratio)
const Scalar max_visco_ratio)
{
const ValueType temp_value = 1. / (1. - (0.8415 / 0.7480 * water_liquid_fraction) );
const ValueType viscosity_ratio = pow(temp_value, 2.5);
@@ -78,10 +78,10 @@ ValueType WIOEmulsionViscosity(const ValueType& oil_viscosity,
}
// oil in water emulsion viscosity
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType OIWEmulsionViscosity(const ValueType& water_viscosity,
const ValueType& water_liquid_fraction,
const double max_visco_ratio)
const Scalar max_visco_ratio)
{
const ValueType temp_value = 1. / (1. - (0.6019 / 0.6410) * (1. - water_liquid_fraction) );
const ValueType viscosity_ratio = pow(temp_value, 2.5);
@@ -95,10 +95,7 @@ ValueType OIWEmulsionViscosity(const ValueType& water_viscosity,
}
namespace Opm {
namespace mswellhelpers
{
namespace Opm::mswellhelpers {
/// Applies umfpack and checks for singularity
template <typename MatrixType, typename VectorType>
@@ -212,17 +209,17 @@ invDX(const MatrixType& D, VectorType x, DeferredLogger& deferred_logger)
return y;
}
template <typename ValueType>
ValueType frictionPressureLoss(const double l, const double diameter,
const double area, const double roughness,
template <typename ValueType, typename Scalar>
ValueType frictionPressureLoss(const Scalar l, const Scalar diameter,
const Scalar area, const Scalar roughness,
const ValueType& density,
const ValueType& w, const ValueType& mu)
{
// Reynolds number
const ValueType re = abs( diameter * w / (area * mu));
constexpr double re_value1 = 2000.;
constexpr double re_value2 = 4000.;
constexpr Scalar re_value1 = 2000.;
constexpr Scalar re_value2 = 4000.;
if (re < re_value1) {
// not using the formula directly because of the division with very small w
@@ -234,7 +231,7 @@ ValueType frictionPressureLoss(const double l, const double diameter,
if (re > re_value2) {
f = haalandFormular(re, diameter, roughness);
} else { // in between
constexpr double f1 = 16. / re_value1;
constexpr Scalar f1 = 16. / re_value1;
const ValueType f2 = haalandFormular(re_value2, diameter, roughness);
f = (f2 - f1) / (re_value2 - re_value1) * (re - re_value1) + f1;
}
@@ -242,19 +239,19 @@ ValueType frictionPressureLoss(const double l, const double diameter,
return 2. * f * l * w * w / (area * area * diameter * density);
}
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType valveContrictionPressureLoss(const ValueType& mass_rate,
const ValueType& density,
const double area_con, const double cv)
const Scalar area_con, const Scalar cv)
{
// the formulation is adjusted a little bit for convinience
// velocity = mass_rate / (density * area) is applied to the original formulation
const double area = (area_con > 1.e-10 ? area_con : 1.e-10);
const Scalar area = (area_con > 1.e-10 ? area_con : 1.e-10);
return mass_rate * mass_rate / (2. * density * cv * cv * area * area);
}
template <typename ValueType>
ValueType velocityHead(const double area, const ValueType& mass_rate,
template <typename ValueType, typename Scalar>
ValueType velocityHead(const Scalar area, const ValueType& mass_rate,
const ValueType& density)
{
// \Note: a factor of 2 is added to the formulation in order to match results from the
@@ -262,21 +259,21 @@ ValueType velocityHead(const double area, const ValueType& mass_rate,
return (mass_rate * mass_rate / (area * area * density));
}
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType emulsionViscosity(const ValueType& water_fraction,
const ValueType& water_viscosity,
const ValueType& oil_fraction,
const ValueType& oil_viscosity,
const SICD& sicd)
{
const double width_transition = sicd.widthTransitionRegion();
const Scalar width_transition = sicd.widthTransitionRegion();
// it is just for now, we should be able to treat it.
if (width_transition <= 0.) {
OPM_THROW(std::runtime_error, "Not handling non-positive transition width now");
}
const double critical_value = sicd.criticalValue();
const Scalar critical_value = sicd.criticalValue();
const ValueType transition_start_value = critical_value - width_transition / 2.0;
const ValueType transition_end_value = critical_value + width_transition / 2.0;
@@ -288,7 +285,7 @@ ValueType emulsionViscosity(const ValueType& water_fraction,
const ValueType water_liquid_fraction = water_fraction / liquid_fraction;
const double max_visco_ratio = sicd.maxViscosityRatio();
const Scalar max_visco_ratio = sicd.maxViscosityRatio();
if (water_liquid_fraction <= transition_start_value) {
return WIOEmulsionViscosity(oil_viscosity, water_liquid_fraction, max_visco_ratio);
} else if (water_liquid_fraction >= transition_end_value) {
@@ -317,30 +314,30 @@ INSTANCE_UMF(2)
INSTANCE_UMF(3)
INSTANCE_UMF(4)
#define INSTANCE_IMPL(...) \
#define INSTANCE_IMPL(T,...) \
template __VA_ARGS__ \
frictionPressureLoss<__VA_ARGS__>(const double, \
const double, \
const double, \
const double, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&); \
frictionPressureLoss(const T, \
const T, \
const T, \
const T, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&); \
template __VA_ARGS__ \
valveContrictionPressureLoss<__VA_ARGS__>(const __VA_ARGS__& mass_rate, \
const __VA_ARGS__& density, \
const double, const double); \
valveContrictionPressureLoss(const __VA_ARGS__& mass_rate, \
const __VA_ARGS__& density, \
const T, const T); \
template __VA_ARGS__ \
velocityHead<__VA_ARGS__>(const double, const __VA_ARGS__&, const __VA_ARGS__&); \
velocityHead(const T, const __VA_ARGS__&, const __VA_ARGS__&); \
template __VA_ARGS__ \
emulsionViscosity<__VA_ARGS__>(const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const SICD&);
emulsionViscosity<__VA_ARGS__,T>(const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const SICD&);
#define INSTANCE_EVAL(Dim) \
INSTANCE_IMPL(DenseAd::Evaluation<double,Dim>)
INSTANCE_IMPL(double, DenseAd::Evaluation<double,Dim>)
INSTANCE_EVAL(3)
INSTANCE_EVAL(4)
@@ -350,6 +347,4 @@ INSTANCE_EVAL(7)
INSTANCE_EVAL(8)
INSTANCE_EVAL(9)
} // namespace mswellhelpers
} // namespace Opm
} // namespace Opm::mswellhelpers

18
opm/simulators/wells/MSWellHelpers.hpp
View File

@@ -29,8 +29,6 @@ namespace Dune {
template<class Matrix> class UMFPack;
}
#include <memory>
namespace Opm {
class DeferredLogger;
@@ -69,26 +67,26 @@ namespace mswellhelpers
// density is density
// roughness is the absolute roughness
// mu is the average phase viscosity
template <typename ValueType>
ValueType frictionPressureLoss(const double l, const double diameter,
const double area, const double roughness,
template <typename ValueType, typename Scalar>
ValueType frictionPressureLoss(const Scalar l, const Scalar diameter,
const Scalar area, const Scalar roughness,
const ValueType& density,
const ValueType& w, const ValueType& mu);
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType valveContrictionPressureLoss(const ValueType& mass_rate,
const ValueType& density,
const double area_con, const double cv);
const Scalar area_con, const Scalar cv);
template <typename ValueType>
ValueType velocityHead(const double area, const ValueType& mass_rate,
template <typename ValueType, typename Scalar>
ValueType velocityHead(const Scalar area, const ValueType& mass_rate,
const ValueType& density);
// calculating the viscosity of oil-water emulsion at local conditons
template <typename ValueType>
template <typename ValueType, typename Scalar>
ValueType emulsionViscosity(const ValueType& water_fraction,
const ValueType& water_viscosity,
const ValueType& oil_fraction,

8
opm/simulators/wells/MultisegmentWell.hpp
View File

@@ -68,7 +68,7 @@ namespace Opm
using typename Base::PressureMatrix;
MultisegmentWell(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const ModelParameters& param,
const RateConverterType& rate_converter,
@@ -76,7 +76,7 @@ namespace Opm
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data);
const std::vector<PerforationData<Scalar>>& perf_data);
void init(const PhaseUsage* phase_usage_arg,
const std::vector<double>& depth_arg,
@@ -197,7 +197,7 @@ namespace Opm
const bool& allow_cf,
std::vector<Value>& cq_s,
Value& perf_press,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const;
template<class Value>
@@ -214,7 +214,7 @@ namespace Opm
const std::vector<Value>& cmix_s,
std::vector<Value>& cq_s,
Value& perf_press,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const;
// compute the fluid properties, such as densities, viscosities, and so on, in the segments

2
opm/simulators/wells/MultisegmentWellSegments.cpp
View File

@@ -624,7 +624,7 @@ pressureDropSpiralICD(const int seg,
// viscosity contribution from the liquid
const EvalWell liquid_viscosity_fraction = liquid_fraction < 1.e-30 ? oil_fraction * oil_viscosity + water_fraction * water_viscosity :
liquid_fraction * mswellhelpers::emulsionViscosity(water_fraction, water_viscosity, oil_fraction, oil_viscosity, sicd);
liquid_fraction * mswellhelpers::emulsionViscosity<EvalWell,Scalar>(water_fraction, water_viscosity, oil_fraction, oil_viscosity, sicd);
const EvalWell mixture_viscosity = liquid_viscosity_fraction + gas_fraction * gas_viscosity;

14
opm/simulators/wells/MultisegmentWell_impl.hpp
View File

@@ -57,7 +57,7 @@ namespace Opm
template <typename TypeTag>
MultisegmentWell<TypeTag>::
MultisegmentWell(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const ModelParameters& param,
const RateConverterType& rate_converter,
@@ -65,7 +65,7 @@ namespace Opm
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data)
const std::vector<PerforationData<Scalar>>& perf_data)
: Base(well, pw_info, time_step, param, rate_converter, pvtRegionIdx, num_components, num_phases, index_of_well, perf_data)
, MSWEval(static_cast<WellInterfaceIndices<FluidSystem,Indices>&>(*this))
, regularize_(false)
@@ -371,7 +371,7 @@ namespace Opm
const Scalar seg_pressure = segment_pressure[seg];
std::vector<Scalar> cq_s(this->num_components_, 0.);
Scalar perf_press = 0.0;
PerforationRates perf_rates;
PerforationRates<Scalar> perf_rates;
computePerfRate(intQuants, mob, Tw, seg, perf, seg_pressure,
allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
@@ -870,7 +870,7 @@ namespace Opm
const std::vector<Value>& cmix_s,
std::vector<Value>& cq_s,
Value& perf_press,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const
{
// pressure difference between the segment and the perforation
@@ -1006,7 +1006,7 @@ namespace Opm
const bool& allow_cf,
std::vector<Value>& cq_s,
Value& perf_press,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const
{
@@ -1874,7 +1874,7 @@ namespace Opm
const std::vector<Scalar> Tw = this->wellIndex(perf, int_quants, trans_mult, wellstate_nupcol);
std::vector<EvalWell> cq_s(this->num_components_, 0.0);
EvalWell perf_press;
PerforationRates perfRates;
PerforationRates<Scalar> perfRates;
computePerfRate(int_quants, mob, Tw, seg, perf, seg_pressure,
allow_cf, cq_s, perf_press, perfRates, deferred_logger);
@@ -2190,7 +2190,7 @@ namespace Opm
const std::vector<Scalar> Tw = this->wellIndex(perf, int_quants, trans_mult, wellstate_nupcol);
std::vector<Scalar> cq_s(this->num_components_, 0.0);
Scalar perf_press = 0.0;
PerforationRates perf_rates;
PerforationRates<Scalar> perf_rates;
computePerfRate(int_quants, mob, Tw, seg, perf, seg_pressure,
allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
for (int comp = 0; comp < this->num_components_; ++comp) {

2
opm/simulators/wells/ParallelWBPCalculation.cpp
View File

@@ -173,7 +173,7 @@ Opm::ParallelWBPCalculation::evalCellSource(Evaluator evalCellSrc)
std::size_t
Opm::ParallelWBPCalculation::
createCalculator(const Well& well,
const ParallelWellInfo& parallelWellInfo,
const ParallelWellInfo<double>& parallelWellInfo,
const std::vector<int>& localConnIdx,
EvaluatorFactory makeWellSourceEvaluator)
{

4
opm/simulators/wells/ParallelWBPCalculation.hpp
View File

@@ -36,7 +36,7 @@
namespace Opm {
class GridDims;
class ParallelWellInfo;
template<class Scalar> class ParallelWellInfo;
class PAvg;
class Well;
}
@@ -112,7 +112,7 @@ public:
/// well.
std::size_t
createCalculator(const Well& well,
const ParallelWellInfo& parallelWellInfo,
const ParallelWellInfo<double>& parallelWellInfo,
const std::vector<int>& localConnIdx,
EvaluatorFactory makeWellSourceEvaluator);

275
opm/simulators/wells/ParallelWellInfo.cpp
View File

@@ -45,15 +45,33 @@ struct CommPolicy<double*>
return 1;
}
};
template<>
struct CommPolicy<float*>
{
using Type = float*;
using IndexedType = float;
using IndexedTypeFlag = Dune::SizeOne;
static const void* getAddress(const float*& v, int index)
{
return v + index;
}
static int getSize(const float*&, int)
{
return 1;
}
};
#endif
}
namespace Opm
{
GlobalPerfContainerFactory::GlobalPerfContainerFactory(const IndexSet& local_indices,
const Parallel::Communication comm,
const int num_local_perfs)
template<class Scalar>
GlobalPerfContainerFactory<Scalar>::
GlobalPerfContainerFactory(const IndexSet& local_indices,
const Parallel::Communication comm,
const int num_local_perfs)
: local_indices_(local_indices), comm_(comm)
{
if ( comm_.size() > 1 )
@@ -101,13 +119,13 @@ GlobalPerfContainerFactory::GlobalPerfContainerFactory(const IndexSet& local_ind
}
}
std::vector<double> GlobalPerfContainerFactory::createGlobal(const std::vector<double>& local_perf_container,
std::size_t num_components) const
template<class Scalar>
std::vector<Scalar> GlobalPerfContainerFactory<Scalar>::
createGlobal(const std::vector<Scalar>& local_perf_container,
std::size_t num_components) const
{
// Could be become templated later.
using Value = double;
using Value = Scalar;
if (comm_.size() > 1)
{
@@ -152,8 +170,11 @@ std::vector<double> GlobalPerfContainerFactory::createGlobal(const std::vector<d
}
}
void GlobalPerfContainerFactory::copyGlobalToLocal(const std::vector<double>& global, std::vector<double>& local,
std::size_t num_components) const
template<class Scalar>
void GlobalPerfContainerFactory<Scalar>::
copyGlobalToLocal(const std::vector<Scalar>& global,
std::vector<Scalar>& local,
std::size_t num_components) const
{
if (global.empty())
{
@@ -182,19 +203,22 @@ void GlobalPerfContainerFactory::copyGlobalToLocal(const std::vector<double>& gl
}
}
int GlobalPerfContainerFactory::numGlobalPerfs() const
template<class Scalar>
int GlobalPerfContainerFactory<Scalar>::numGlobalPerfs() const
{
return num_global_perfs_;
}
CommunicateAboveBelow::CommunicateAboveBelow([[maybe_unused]] const Parallel::Communication& comm)
template<class Scalar>
CommunicateAboveBelow<Scalar>::
CommunicateAboveBelow([[maybe_unused]] const Parallel::Communication& comm)
#if HAVE_MPI
: comm_(comm), interface_(comm_)
#endif
{}
void CommunicateAboveBelow::clear()
template<class Scalar>
void CommunicateAboveBelow<Scalar>::clear()
{
#if HAVE_MPI
above_indices_ = {};
@@ -205,7 +229,8 @@ void CommunicateAboveBelow::clear()
num_local_perfs_ = 0;
}
void CommunicateAboveBelow::beginReset()
template<class Scalar>
void CommunicateAboveBelow<Scalar>::beginReset()
{
clear();
#if HAVE_MPI
@@ -217,7 +242,8 @@ void CommunicateAboveBelow::beginReset()
#endif
}
int CommunicateAboveBelow::endReset()
template<class Scalar>
int CommunicateAboveBelow<Scalar>::endReset()
{
#if HAVE_MPI
if (comm_.size() > 1)
@@ -226,7 +252,7 @@ int CommunicateAboveBelow::endReset()
current_indices_.endResize();
remote_indices_.setIndexSets(current_indices_, above_indices_, comm_);
// It is mandatory to not set includeSelf to true, as this will skip some entries.
remote_indices_.rebuild<true>();
remote_indices_.template rebuild<true>();
using FromSet = Dune::EnumItem<Attribute,owner>;
using ToSet = Dune::AllSet<Attribute>;
interface_.build(remote_indices_, FromSet(), ToSet());
@@ -236,8 +262,9 @@ int CommunicateAboveBelow::endReset()
return num_local_perfs_;
}
template<class Scalar>
template<class RAIterator>
void CommunicateAboveBelow::partialSumPerfValues(RAIterator begin, RAIterator end) const
void CommunicateAboveBelow<Scalar>::partialSumPerfValues(RAIterator begin, RAIterator end) const
{
if (this->comm_.size() < 2)
{
@@ -295,28 +322,27 @@ void CommunicateAboveBelow::partialSumPerfValues(RAIterator begin, RAIterator en
}
}
using dIter = typename std::vector<double>::iterator;
template void CommunicateAboveBelow::partialSumPerfValues<dIter>(dIter begin, dIter end) const;
template<class Scalar>
struct CopyGatherScatter
{
static const double& gather(const double* a, std::size_t i)
static const Scalar& gather(const Scalar* a, std::size_t i)
{
return a[i];
}
static void scatter(double* a, const double& v, std::size_t i)
static void scatter(Scalar* a, const Scalar& v, std::size_t i)
{
a[i] = v;
}
};
std::vector<double> CommunicateAboveBelow::communicateAbove(double first_above,
const double* current,
std::size_t size)
template<class Scalar>
std::vector<Scalar> CommunicateAboveBelow<Scalar>::
communicateAbove(Scalar first_above,
const Scalar* current,
std::size_t size)
{
std::vector<double> above(size, first_above);
std::vector<Scalar> above(size, first_above);
#if HAVE_MPI
if (comm_.size() > 1)
@@ -326,7 +352,7 @@ std::vector<double> CommunicateAboveBelow::communicateAbove(double first_above,
// passing const double*& and double* as parameter is
// incompatible with function decl template<Data> forward(const Data&, Data&))
// That would need the first argument to be double* const&
communicator_.forward<CopyGatherScatter>(const_cast<double*>(current), aboveData);
communicator_.forward<CopyGatherScatter<Scalar>>(const_cast<Scalar*>(current), aboveData);
}
else
#endif
@@ -339,11 +365,14 @@ std::vector<double> CommunicateAboveBelow::communicateAbove(double first_above,
}
return above;
}
std::vector<double> CommunicateAboveBelow::communicateBelow(double last_below,
const double* current,
std::size_t size)
template<class Scalar>
std::vector<Scalar> CommunicateAboveBelow<Scalar>::
communicateBelow(Scalar last_below,
const Scalar* current,
std::size_t size)
{
std::vector<double> below(size, last_below);
std::vector<Scalar> below(size, last_below);
#if HAVE_MPI
if (comm_.size() > 1)
@@ -353,7 +382,7 @@ std::vector<double> CommunicateAboveBelow::communicateBelow(double last_below,
// passing const double*& and double* as parameter is
// incompatible with function decl template<Data> backward(Data&, const Data&)
// That would need the first argument to be double* const&
communicator_.backward<CopyGatherScatter>(belowData, const_cast<double*>(current));
communicator_.backward<CopyGatherScatter<Scalar>>(belowData, const_cast<Scalar*>(current));
}
else
#endif
@@ -367,9 +396,11 @@ std::vector<double> CommunicateAboveBelow::communicateBelow(double last_below,
return below;
}
void CommunicateAboveBelow::pushBackEclIndex([[maybe_unused]] int above,
[[maybe_unused]] int current,
[[maybe_unused]] bool isOwner)
template<class Scalar>
void CommunicateAboveBelow<Scalar>::
pushBackEclIndex([[maybe_unused]] int above,
[[maybe_unused]] int current,
[[maybe_unused]] bool isOwner)
{
#if HAVE_MPI
if (comm_.size() > 1)
@@ -386,8 +417,8 @@ void CommunicateAboveBelow::pushBackEclIndex([[maybe_unused]] int above,
++num_local_perfs_;
}
void ParallelWellInfo::DestroyComm::operator()(Parallel::Communication* comm)
template<class Scalar>
void ParallelWellInfo<Scalar>::DestroyComm::operator()(Parallel::Communication* comm)
{
#if HAVE_MPI
// Only delete custom communicators.
@@ -406,28 +437,31 @@ void ParallelWellInfo::DestroyComm::operator()(Parallel::Communication* comm)
delete comm;
}
const CommunicateAboveBelow::IndexSet& CommunicateAboveBelow::getIndexSet() const
template<class Scalar>
const typename CommunicateAboveBelow<Scalar>::IndexSet&
CommunicateAboveBelow<Scalar>::getIndexSet() const
{
return current_indices_;
}
int CommunicateAboveBelow::numLocalPerfs() const
template<class Scalar>
int CommunicateAboveBelow<Scalar>::numLocalPerfs() const
{
return num_local_perfs_;
}
ParallelWellInfo::ParallelWellInfo(const std::string& name,
bool hasLocalCells)
template<class Scalar>
ParallelWellInfo<Scalar>::ParallelWellInfo(const std::string& name,
bool hasLocalCells)
: name_(name), hasLocalCells_ (hasLocalCells),
isOwner_(true), rankWithFirstPerf_(-1),
comm_(new Parallel::Communication(Dune::MPIHelper::getLocalCommunicator())),
commAboveBelow_(new CommunicateAboveBelow(*comm_))
commAboveBelow_(new CommunicateAboveBelow<Scalar>(*comm_))
{}
ParallelWellInfo::ParallelWellInfo(const std::pair<std::string, bool>& well_info,
[[maybe_unused]] Parallel::Communication allComm)
template<class Scalar>
ParallelWellInfo<Scalar>::ParallelWellInfo(const std::pair<std::string, bool>& well_info,
[[maybe_unused]] Parallel::Communication allComm)
: name_(well_info.first), hasLocalCells_(well_info.second),
rankWithFirstPerf_(-1)
{
@@ -439,12 +473,12 @@ ParallelWellInfo::ParallelWellInfo(const std::pair<std::string, bool>& well_info
#else
comm_.reset(new Parallel::Communication(Dune::MPIHelper::getLocalCommunicator()));
#endif
commAboveBelow_.reset(new CommunicateAboveBelow(*comm_));
commAboveBelow_.reset(new CommunicateAboveBelow<Scalar>(*comm_));
isOwner_ = (comm_->rank() == 0);
}
void ParallelWellInfo::communicateFirstPerforation(bool hasFirst)
template<class Scalar>
void ParallelWellInfo<Scalar>::communicateFirstPerforation(bool hasFirst)
{
int first = hasFirst;
std::vector<int> firstVec(comm_->size());
@@ -455,29 +489,32 @@ void ParallelWellInfo::communicateFirstPerforation(bool hasFirst)
rankWithFirstPerf_ = found - firstVec.begin();
}
void ParallelWellInfo::pushBackEclIndex(int above, int current)
template<class Scalar>
void ParallelWellInfo<Scalar>::pushBackEclIndex(int above, int current)
{
commAboveBelow_->pushBackEclIndex(above, current);
}
void ParallelWellInfo::beginReset()
template<class Scalar>
void ParallelWellInfo<Scalar>::beginReset()
{
commAboveBelow_->beginReset();
}
void ParallelWellInfo::endReset()
template<class Scalar>
void ParallelWellInfo<Scalar>::endReset()
{
int local_num_perfs = commAboveBelow_->endReset();
globalPerfCont_
.reset(new GlobalPerfContainerFactory(commAboveBelow_->getIndexSet(),
*comm_,
local_num_perfs));
.reset(new GlobalPerfContainerFactory<Scalar>(commAboveBelow_->getIndexSet(),
*comm_,
local_num_perfs));
}
template<class Scalar>
template<typename It>
typename It::value_type
ParallelWellInfo::sumPerfValues(It begin, It end) const
ParallelWellInfo<Scalar>::sumPerfValues(It begin, It end) const
{
using V = typename It::value_type;
/// \todo cater for overlap later. Currently only owner
@@ -485,19 +522,16 @@ ParallelWellInfo::sumPerfValues(It begin, It end) const
return communication().sum(local);
}
using cdIter = typename std::vector<double>::const_iterator;
template typename cdIter::value_type ParallelWellInfo::sumPerfValues<cdIter>(cdIter,cdIter) const;
template typename dIter::value_type ParallelWellInfo::sumPerfValues<dIter>(dIter,dIter) const;
void ParallelWellInfo::clear()
template<class Scalar>
void ParallelWellInfo<Scalar>::clear()
{
commAboveBelow_->clear();
globalPerfCont_.reset();
}
template<class Scalar>
template<class T>
T ParallelWellInfo::broadcastFirstPerforationValue(const T& t) const
T ParallelWellInfo<Scalar>::broadcastFirstPerforationValue(const T& t) const
{
T res = t;
if (rankWithFirstPerf_ >= 0) {
@@ -515,41 +549,50 @@ T ParallelWellInfo::broadcastFirstPerforationValue(const T& t) const
return res;
}
template int ParallelWellInfo::broadcastFirstPerforationValue<int>(const int&) const;
template double ParallelWellInfo::broadcastFirstPerforationValue<double>(const double&) const;
template int ParallelWellInfo<double>::broadcastFirstPerforationValue<int>(const int&) const;
template double ParallelWellInfo<double>::broadcastFirstPerforationValue<double>(const double&) const;
std::vector<double> ParallelWellInfo::communicateAboveValues(double zero_value,
const double* current_values,
std::size_t size) const
template<class Scalar>
std::vector<Scalar> ParallelWellInfo<Scalar>::
communicateAboveValues(Scalar zero_value,
const Scalar* current_values,
std::size_t size) const
{
return commAboveBelow_->communicateAbove(zero_value, current_values,
size);
}
std::vector<double> ParallelWellInfo::communicateAboveValues(double zero_value,
const std::vector<double>& current_values) const
template<class Scalar>
std::vector<Scalar> ParallelWellInfo<Scalar>::
communicateAboveValues(Scalar zero_value,
const std::vector<Scalar>& current_values) const
{
return commAboveBelow_->communicateAbove(zero_value, current_values.data(),
current_values.size());
}
std::vector<double> ParallelWellInfo::communicateBelowValues(double last_value,
const double* current_values,
std::size_t size) const
template<class Scalar>
std::vector<Scalar> ParallelWellInfo<Scalar>::
communicateBelowValues(Scalar last_value,
const Scalar* current_values,
std::size_t size) const
{
return commAboveBelow_->communicateBelow(last_value, current_values,
size);
}
std::vector<double> ParallelWellInfo::communicateBelowValues(double last_value,
const std::vector<double>& current_values) const
template<class Scalar>
std::vector<Scalar> ParallelWellInfo<Scalar>::
communicateBelowValues(Scalar last_value,
const std::vector<Scalar>& current_values) const
{
return commAboveBelow_->communicateBelow(last_value, current_values.data(),
current_values.size());
}
const GlobalPerfContainerFactory&
ParallelWellInfo::getGlobalPerfContainerFactory() const
template<class Scalar>
const GlobalPerfContainerFactory<Scalar>&
ParallelWellInfo<Scalar>::getGlobalPerfContainerFactory() const
{
if(globalPerfCont_)
return *globalPerfCont_;
@@ -558,12 +601,15 @@ ParallelWellInfo::getGlobalPerfContainerFactory() const
"No ecl indices have been added via beginReset, pushBackEclIndex, endReset");
}
bool operator<(const ParallelWellInfo& well1, const ParallelWellInfo& well2)
template<class Scalar>
bool operator<(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2)
{
return well1.name() < well2.name() || (! (well2.name() < well1.name()) && well1.hasLocalCells() < well2.hasLocalCells());
return well1.name() < well2.name() || (! (well2.name() < well1.name()) &&
well1.hasLocalCells() < well2.hasLocalCells());
}
bool operator==(const ParallelWellInfo& well1, const ParallelWellInfo& well2)
template<class Scalar>
bool operator==(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2)
{
bool ret = well1.name() == well2.name() && well1.hasLocalCells() == well2.hasLocalCells()
&& well1.isOwner() == well2.isOwner();
@@ -575,56 +621,67 @@ bool operator==(const ParallelWellInfo& well1, const ParallelWellInfo& well2)
return ret;
}
bool operator!=(const ParallelWellInfo& well1, const ParallelWellInfo& well2)
template<class Scalar>
bool operator!=(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2)
{
return ! (well1 == well2);
}
bool operator<(const std::pair<std::string, bool>& pair, const ParallelWellInfo& well)
template<class Scalar>
bool operator<(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well)
{
return pair.first < well.name() || ( !( well.name() < pair.first ) && pair.second < well.hasLocalCells() );
}
bool operator<( const ParallelWellInfo& well, const std::pair<std::string, bool>& pair)
template<class Scalar>
bool operator<( const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair)
{
return well.name() < pair.first || ( !( pair.first < well.name() ) && well.hasLocalCells() < pair.second );
}
bool operator==(const std::pair<std::string, bool>& pair, const ParallelWellInfo& well)
template<class Scalar>
bool operator==(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well)
{
return pair.first == well.name() && pair.second == well.hasLocalCells();
}
bool operator==(const ParallelWellInfo& well, const std::pair<std::string, bool>& pair)
template<class Scalar>
bool operator==(const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair)
{
return pair == well;
}
bool operator!=(const std::pair<std::string, bool>& pair, const ParallelWellInfo& well)
template<class Scalar>
bool operator!=(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well)
{
return pair.first != well.name() || pair.second != well.hasLocalCells();
}
bool operator!=(const ParallelWellInfo& well, const std::pair<std::string, bool>& pair)
template<class Scalar>
bool operator!=(const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair)
{
return pair != well;
}
CheckDistributedWellConnections::CheckDistributedWellConnections(const Well& well,
const ParallelWellInfo& info)
template<class Scalar>
CheckDistributedWellConnections<Scalar>::
CheckDistributedWellConnections(const Well& well,
const ParallelWellInfo<Scalar>& info)
: well_(well), pwinfo_(info)
{
foundConnections_.resize(well.getConnections().size(), 0);
}
void
CheckDistributedWellConnections::connectionFound(std::size_t i)
template<class Scalar>
void CheckDistributedWellConnections<Scalar>::
connectionFound(std::size_t i)
{
foundConnections_[i] = 1;
}
bool
CheckDistributedWellConnections::checkAllConnectionsFound()
template<class Scalar>
bool CheckDistributedWellConnections<Scalar>::
checkAllConnectionsFound()
{
// Ecl does not hold any information of remote connections.
assert(pwinfo_.communication().max(foundConnections_.size()) == foundConnections_.size());
@@ -652,4 +709,30 @@ CheckDistributedWellConnections::checkAllConnectionsFound()
}
return !missingCells;
}
template<class Scalar> using dIter = typename std::vector<Scalar>::iterator;
template<class Scalar> using cdIter = typename std::vector<Scalar>::const_iterator;
#define INSTANCE(T) \
template class CheckDistributedWellConnections<T>; \
template class CommunicateAboveBelow<T>; \
template class GlobalPerfContainerFactory<T>; \
template class ParallelWellInfo<T>; \
template typename cdIter<T>::value_type \
ParallelWellInfo<T>::sumPerfValues<cdIter<T>>(cdIter<T>,cdIter<T>) const; \
template typename dIter<T>::value_type \
ParallelWellInfo<T>::sumPerfValues<dIter<T>>(dIter<T>,dIter<T>) const; \
template void CommunicateAboveBelow<T>::partialSumPerfValues<dIter<T>>(dIter<T>, dIter<T>) const; \
template bool operator<(const ParallelWellInfo<T>&, const ParallelWellInfo<T>&); \
template bool operator<(const ParallelWellInfo<T>&, const std::pair<std::string, bool>&); \
template bool operator<(const std::pair<std::string, bool>&, const ParallelWellInfo<T>&); \
template bool operator==(const ParallelWellInfo<T>&, const ParallelWellInfo<T>&); \
template bool operator==(const ParallelWellInfo<T>& well, const std::pair<std::string, bool>&); \
template bool operator==(const std::pair<std::string, bool>&, const ParallelWellInfo<T>&); \
template bool operator!=(const ParallelWellInfo<T>&, const ParallelWellInfo<T>&); \
template bool operator!=(const std::pair<std::string, bool>&, const ParallelWellInfo<T>&); \
template bool operator!=(const ParallelWellInfo<T>&, const std::pair<std::string, bool>&);
INSTANCE(double)
} // end namespace Opm

77
opm/simulators/wells/ParallelWellInfo.hpp
View File

@@ -29,13 +29,13 @@
#include <memory>
namespace Opm
{
namespace Opm {
class Well;
/// \brief Class to facilitate getting values associated with the above/below perforation
///
template<class Scalar>
class CommunicateAboveBelow
{
public:
@@ -83,8 +83,8 @@ public:
/// \param current C-array of the values at the perforations
/// \param size The size of the C-array and the returned vector
/// \return a vector containing the values for the perforation above.
std::vector<double> communicateAbove(double first_value,
const double* current,
std::vector<Scalar> communicateAbove(Scalar first_value,
const Scalar* current,
std::size_t size);
/// \brief Creates an array of values for the perforation below.
@@ -92,8 +92,8 @@ public:
/// \param current C-array of the values at the perforations
/// \param size The size of the C-array and the returned vector
/// \return a vector containing the values for the perforation above.
std::vector<double> communicateBelow(double first_value,
const double* current,
std::vector<Scalar> communicateBelow(Scalar first_value,
const Scalar* current,
std::size_t size);
/// \brief Do a (in place) partial sum on values attached to all perforations.
@@ -131,11 +131,12 @@ private:
/// even if the data is distributed. This class is supposed to help with that by
/// computing the global data arrays for the well and copy computed values back to
/// the distributed representation.
template<class Scalar>
class GlobalPerfContainerFactory
{
public:
using IndexSet = CommunicateAboveBelow::IndexSet;
using Attribute = CommunicateAboveBelow::Attribute;
using IndexSet = typename CommunicateAboveBelow<Scalar>::IndexSet;
using Attribute = typename CommunicateAboveBelow<Scalar>::Attribute;
using GlobalIndex = typename IndexSet::IndexPair::GlobalIndex;
/// \brief Constructor
@@ -149,14 +150,14 @@ public:
/// \param num_components the number of components per perforation.
/// \return A container with values attached to all perforations of a well.
/// Values are ordered by the index of the perforation in the ECL schedule.
std::vector<double> createGlobal(const std::vector<double>& local_perf_container,
std::vector<Scalar> createGlobal(const std::vector<Scalar>& local_perf_container,
std::size_t num_components) const;
/// \brief Copies the values of the global perforation to the local representation
/// \param global values attached to all peforations of a well (as if the well would live on one process)
/// \param num_components the number of components per perforation.
/// \param[out] local The values attached to the local perforations only.
void copyGlobalToLocal(const std::vector<double>& global, std::vector<double>& local,
void copyGlobalToLocal(const std::vector<Scalar>& global, std::vector<Scalar>& local,
std::size_t num_components) const;
int numGlobalPerfs() const;
@@ -180,6 +181,7 @@ private:
/// \brief Class encapsulating some information about parallel wells
///
/// e.g. It provides a communicator for well information
template<class Scalar>
class ParallelWellInfo
{
public:
@@ -218,30 +220,30 @@ public:
/// \param current C-array of the values at the perforations
/// \param size The size of the C-array and the returned vector
/// \return a vector containing the values for the perforation above.
std::vector<double> communicateAboveValues(double first_value,
const double* current,
std::vector<Scalar> communicateAboveValues(Scalar first_value,
const Scalar* current,
std::size_t size) const;
/// \brief Creates an array of values for the perforation above.
/// \param first_value Value to use for above of the first perforation
/// \param current vector of current values
std::vector<double> communicateAboveValues(double first_value,
const std::vector<double>& current) const;
std::vector<Scalar> communicateAboveValues(Scalar first_value,
const std::vector<Scalar>& current) const;
/// \brief Creates an array of values for the perforation below.
/// \param last_value Value to use for below of the last perforation
/// \param current C-array of the values at the perforations
/// \param size The size of the C-array and the returned vector
/// \return a vector containing the values for the perforation above.
std::vector<double> communicateBelowValues(double last_value,
const double* current,
std::vector<Scalar> communicateBelowValues(Scalar last_value,
const Scalar* current,
std::size_t size) const;
/// \brief Creates an array of values for the perforation above.
/// \param last_value Value to use for below of the last perforation
/// \param current vector of current values
std::vector<double> communicateBelowValues(double last_value,
const std::vector<double>& current) const;
std::vector<Scalar> communicateBelowValues(Scalar last_value,
const std::vector<Scalar>& current) const;
/// \brief Adds information about the ecl indices of the perforations.
///
@@ -301,7 +303,7 @@ public:
/// That is a container that holds data for every perforation no matter where
/// it is stored. Container is ordered via ascendings index of the perforations
/// in the ECL schedule.
const GlobalPerfContainerFactory& getGlobalPerfContainerFactory() const;
const GlobalPerfContainerFactory<Scalar>& getGlobalPerfContainerFactory() const;
private:
@@ -326,19 +328,20 @@ private:
std::unique_ptr<Parallel::Communication, DestroyComm> comm_;
/// \brief used to communicate the values for the perforation above.
std::unique_ptr<CommunicateAboveBelow> commAboveBelow_;
std::unique_ptr<CommunicateAboveBelow<Scalar>> commAboveBelow_;
std::unique_ptr<GlobalPerfContainerFactory> globalPerfCont_;
std::unique_ptr<GlobalPerfContainerFactory<Scalar>> globalPerfCont_;
};
/// \brief Class checking that all connections are on active cells
///
/// Works for distributed wells, too
template<class Scalar>
class CheckDistributedWellConnections
{
public:
CheckDistributedWellConnections(const Well& well,
const ParallelWellInfo& info);
const ParallelWellInfo<Scalar>& info);
/// \brief Inidicate that the i-th completion was found
///
@@ -351,26 +354,36 @@ public:
private:
std::vector<std::size_t> foundConnections_;
const Well& well_;
const ParallelWellInfo& pwinfo_;
const ParallelWellInfo<Scalar>& pwinfo_;
};
bool operator<(const ParallelWellInfo& well1, const ParallelWellInfo& well2);
template<class Scalar>
bool operator<(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2);
bool operator==(const ParallelWellInfo& well1, const ParallelWellInfo& well2);
template<class Scalar>
bool operator==(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2);
bool operator!=(const ParallelWellInfo& well1, const ParallelWellInfo& well2);
template<class Scalar>
bool operator!=(const ParallelWellInfo<Scalar>& well1, const ParallelWellInfo<Scalar>& well2);
bool operator<(const std::pair<std::string, bool>& pair, const ParallelWellInfo& well);
template<class Scalar>
bool operator<(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well);
bool operator<( const ParallelWellInfo& well, const std::pair<std::string, bool>& pair);
template<class Scalar>
bool operator<( const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair);
bool operator==(const std::pair<std::string, bool>& pair, const ParallelWellInfo& well);
template<class Scalar>
bool operator==(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well);
bool operator==(const ParallelWellInfo& well, const std::pair<std::string, bool>& pair);
template<class Scalar>
bool operator==(const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair);
bool operator!=(const std::pair<std::string, bool>& pair, const ParallelWellInfo& well);
template<class Scalar>
bool operator!=(const std::pair<std::string, bool>& pair, const ParallelWellInfo<Scalar>& well);
bool operator!=(const ParallelWellInfo& well, const std::pair<std::string, bool>& pair);
template<class Scalar>
bool operator!=(const ParallelWellInfo<Scalar>& well, const std::pair<std::string, bool>& pair);
} // end namespace Opm
#endif // OPM_PARALLELWELLINFO_HEADER_INCLUDED

17
opm/simulators/wells/PerforationData.hpp
View File

@@ -22,26 +22,27 @@
#include <cstddef>
namespace Opm
{
namespace Opm {
/// Static data associated with a well perforation.
template<class Scalar>
struct PerforationData
{
int cell_index{};
double connection_transmissibility_factor{};
double connection_d_factor{};
Scalar connection_transmissibility_factor{};
Scalar connection_d_factor{};
int satnum_id{};
/// \brief The original index of the perforation in ECL Schedule
std::size_t ecl_index{};
};
template<class Scalar>
struct PerforationRates
{
double dis_gas = 0.0;
double dis_gas_in_water = 0.0;
double vap_oil = 0.0;
double vap_wat = 0.0;
Scalar dis_gas = 0.0;
Scalar dis_gas_in_water = 0.0;
Scalar vap_oil = 0.0;
Scalar vap_wat = 0.0;
};
} // namespace Opm

95
opm/simulators/wells/RateConverter.cpp
View File

@@ -31,12 +31,12 @@
namespace {
template <typename Rates>
std::pair<double, double>
template <typename Scalar, typename Rates>
std::pair<Scalar, Scalar>
dissolvedVaporisedRatio(const int io,
const int ig,
const double rs,
const double rv,
const Scalar rs,
const Scalar rv,
const Rates& surface_rates)
{
if ((io < 0) || (ig < 0)) {
@@ -49,8 +49,8 @@ dissolvedVaporisedRatio(const int io,
const auto Rv = surface_rates[io] / (surface_rates[ig] + eps);
return {
std::clamp(static_cast<double>(Rs), 0.0, rs),
std::clamp(static_cast<double>(Rv), 0.0, rv)
std::clamp(static_cast<Scalar>(Rs), Scalar{0.0}, rs),
std::clamp(static_cast<Scalar>(Rv), Scalar{0.0}, rv)
};
}
@@ -80,8 +80,8 @@ sumRates(std::unordered_map<RegionId,Attributes>& attributes_hpv,
comm.sum(ra.data.data(), ra.data.size());
assert(ra.pv > 0.);
}
const double pv_sum = ra.pv;
for (double& d : ra.data)
const Scalar pv_sum = ra.pv;
for (Scalar & d : ra.data)
d /= pv_sum;
ra.pv = pv_sum;
}
@@ -94,9 +94,9 @@ calcInjCoeff(const RegionId r, const int pvtRegionIdx, Coeff& coeff) const
{
const auto& pu = phaseUsage_;
const auto& ra = attr_.attributes(r);
const double p = ra.pressure;
const double T = ra.temperature;
const double saltConcentration = ra.saltConcentration;
const Scalar p = ra.pressure;
const Scalar T = ra.temperature;
const Scalar saltConcentration = ra.saltConcentration;
const int iw = RegionAttributeHelpers::PhasePos::water(pu);
const int io = RegionAttributeHelpers::PhasePos::oil (pu);
@@ -107,18 +107,29 @@ calcInjCoeff(const RegionId r, const int pvtRegionIdx, Coeff& coeff) const
if (RegionAttributeHelpers::PhaseUsed::water(pu)) {
// q[w]_r = q[w]_s / bw
const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0, saltConcentration);
const Scalar bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx,
T,
p,
Scalar{0.0},
saltConcentration);
coeff[iw] = 1.0 / bw;
}
if (RegionAttributeHelpers::PhaseUsed::oil(pu)) {
const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0);
const Scalar bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx,
T,
p,
Scalar{0.0});
coeff[io] += 1.0 / bo;
}
if (RegionAttributeHelpers::PhaseUsed::gas(pu)) {
const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0, 0.0);
const Scalar bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx,
T,
p,
Scalar{0.0},
Scalar{0.0});
coeff[ig] += 1.0 / bg;
}
}
@@ -154,15 +165,15 @@ calcCoeff(const RegionId r, const int pvtRegionIdx, const Rates& surface_rates,
template <class FluidSystem, class Region>
template <class Coeff>
void SurfaceToReservoirVoidage<FluidSystem,Region>::
calcCoeff( const int pvtRegionIdx,
const double p,
const double Rs,
const double Rv,
const double Rsw,
const double Rvw,
const double T,
const double saltConcentration,
Coeff& coeff) const
calcCoeff(const int pvtRegionIdx,
const Scalar p,
const Scalar Rs,
const Scalar Rv,
const Scalar Rsw,
const Scalar Rvw,
const Scalar T,
const Scalar saltConcentration,
Coeff& coeff) const
{
const auto& pu = phaseUsage_;
@@ -173,14 +184,14 @@ calcCoeff( const int pvtRegionIdx,
std::fill(& coeff[0], & coeff[0] + phaseUsage_.num_phases, 0.0);
// Determinant of 'R' matrix
const double detRw = 1.0 - (Rsw * Rvw);
const Scalar detRw = 1.0 - (Rsw * Rvw);
if (RegionAttributeHelpers::PhaseUsed::water(pu)) {
// q[w]_r = 1/(bw * (1 - rsw*rvw)) * (q[w]_s - rvw*q[g]_s)
const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rsw, saltConcentration);
const Scalar bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rsw, saltConcentration);
const double den = bw * detRw;
const Scalar den = bw * detRw;
coeff[iw] += 1.0 / den;
@@ -190,7 +201,7 @@ calcCoeff( const int pvtRegionIdx,
}
// Determinant of 'R' matrix
const double detR = 1.0 - (Rs * Rv);
const Scalar detR = 1.0 - (Rs * Rv);
// Currently we only support either gas in water or gas in oil
// not both
@@ -204,8 +215,8 @@ calcCoeff( const int pvtRegionIdx,
if (RegionAttributeHelpers::PhaseUsed::oil(pu)) {
// q[o]_r = 1/(bo * (1 - rs*rv)) * (q[o]_s - rv*q[g]_s)
const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
const double den = bo * detR;
const Scalar bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
const Scalar den = bo * detR;
coeff[io] += 1.0 / den;
@@ -216,15 +227,15 @@ calcCoeff( const int pvtRegionIdx,
if (RegionAttributeHelpers::PhaseUsed::gas(pu)) {
// q[g]_r = 1/(bg * (1 - rs*rv)) * (q[g]_s - rs*q[o]_s)
const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, Rvw);
const Scalar bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, Rvw);
if (FluidSystem::enableDissolvedGasInWater()) {
const double denw = bg * detRw;
const Scalar denw = bg * detRw;
coeff[ig] += 1.0 / denw;
if (RegionAttributeHelpers::PhaseUsed::water(pu)) {
coeff[iw] -= Rsw / denw;
}
} else {
const double den = bg * detR;
const Scalar den = bg * detR;
coeff[ig] += 1.0 / den;
if (RegionAttributeHelpers::PhaseUsed::oil(pu)) {
coeff[io] -= Rs / den;
@@ -238,13 +249,13 @@ template <class FluidSystem, class Region>
template <typename SurfaceRates, typename VoidageRates>
void SurfaceToReservoirVoidage<FluidSystem,Region>::
calcReservoirVoidageRates(const int pvtRegionIdx,
const double p,
const double rs,
const double rv,
const double rsw,
const double rvw,
const double T,
const double saltConcentration,
const Scalar p,
const Scalar rs,
const Scalar rv,
const Scalar rsw,
const Scalar rvw,
const Scalar T,
const Scalar saltConcentration,
const SurfaceRates& surface_rates,
VoidageRates& voidage_rates) const
{
@@ -347,10 +358,10 @@ calcReservoirVoidageRates(const RegionId r,
template <class FluidSystem, class Region>
template <class Rates>
std::pair<double, double>
std::pair<typename FluidSystem::Scalar, typename FluidSystem::Scalar>
SurfaceToReservoirVoidage<FluidSystem,Region>::
inferDissolvedVaporisedRatio(const double rsMax,
const double rvMax,
inferDissolvedVaporisedRatio(const Scalar rsMax,
const Scalar rvMax,
const Rates& surface_rates) const
{
const auto io = RegionAttributeHelpers::PhasePos::oil(this->phaseUsage_);

72
opm/simulators/wells/RateConverter.hpp
View File

@@ -69,6 +69,8 @@ namespace Opm {
template <class FluidSystem, class Region>
class SurfaceToReservoirVoidage {
public:
using Scalar = typename FluidSystem::Scalar;
/**
* Constructor.
*
@@ -130,12 +132,12 @@ namespace Opm {
const auto& intQuants = elemCtx.intensiveQuantities(/*spaceIdx=*/0, /*timeIdx=*/0);
const auto& fs = intQuants.fluidState();
// use pore volume weighted averages.
const double pv_cell =
const Scalar pv_cell =
simulator.model().dofTotalVolume(cellIdx)
* intQuants.porosity().value();
// only count oil and gas filled parts of the domain
double hydrocarbon = 1.0;
Scalar hydrocarbon = 1.0;
const auto& pu = phaseUsage_;
if (RegionAttributeHelpers::PhaseUsed::water(pu)) {
hydrocarbon -= fs.saturation(FluidSystem::waterPhaseIdx).value();
@@ -145,7 +147,7 @@ namespace Opm {
assert(reg >= 0);
// sum p, rs, rv, and T.
const double hydrocarbonPV = pv_cell*hydrocarbon;
const Scalar hydrocarbonPV = pv_cell*hydrocarbon;
if (hydrocarbonPV > 0.) {
auto& attr = attributes_hpv[reg];
attr.pv += hydrocarbonPV;
@@ -249,15 +251,15 @@ namespace Opm {
template <class Coeff>
void
calcCoeff( const int pvtRegionIdx,
const double p,
const double rs,
const double rv,
const double rsw,
const double rvw,
const double T,
const double saltConcentration,
Coeff& coeff) const;
calcCoeff( const int pvtRegionIdx,
const Scalar p,
const Scalar rs,
const Scalar rv,
const Scalar rsw,
const Scalar rvw,
const Scalar T,
const Scalar saltConcentration,
Coeff& coeff) const;
template <class Coeff>
void
@@ -327,20 +329,20 @@ namespace Opm {
*/
template <typename SurfaceRates, typename VoidageRates>
void calcReservoirVoidageRates(const int pvtRegionIdx,
const double p,
const double rs,
const double rv,
const double rsw,
const double rvw,
const double T,
const double saltConcentration,
const Scalar p,
const Scalar rs,
const Scalar rv,
const Scalar rsw,
const Scalar rvw,
const Scalar T,
const Scalar saltConcentration,
const SurfaceRates& surface_rates,
VoidageRates& voidage_rates) const;
template <class Rates>
std::pair<double, double>
inferDissolvedVaporisedRatio(const double rsMax,
const double rvMax,
std::pair<Scalar, Scalar>
inferDissolvedVaporisedRatio(const Scalar rsMax,
const Scalar rvMax,
const Rates& surface_rates) const;
/**
@@ -357,13 +359,13 @@ namespace Opm {
*/
template <class SolventModule>
void
calcCoeffSolvent(const RegionId r, const int pvtRegionIdx, double& coeff) const
calcCoeffSolvent(const RegionId r, const int pvtRegionIdx, Scalar& coeff) const
{
const auto& ra = attr_.attributes(r);
const double p = ra.pressure;
const double T = ra.temperature;
const Scalar p = ra.pressure;
const Scalar T = ra.temperature;
const auto& solventPvt = SolventModule::solventPvt();
const double bs = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
const Scalar bs = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
coeff = 1.0 / bs;
}
@@ -406,15 +408,15 @@ namespace Opm {
return *this;
}
std::array<double,8> data;
double& pressure;
double& temperature;
double& rs;
double& rv;
double& rsw;
double& rvw;
double& pv;
double& saltConcentration;
std::array<Scalar,8> data;
Scalar& pressure;
Scalar& temperature;
Scalar& rs;
Scalar& rv;
Scalar& rsw;
Scalar& rvw;
Scalar& pv;
Scalar& saltConcentration;
};
void sumRates(std::unordered_map<RegionId,Attributes>& attributes_hpv,

31
opm/simulators/wells/RegionAttributeHelpers.hpp
View File

@@ -28,7 +28,6 @@
#include <dune/grid/common/gridenums.hh>
#include <algorithm>
#include <cmath>
#include <memory>
#include <stdexcept>
#include <type_traits>
@@ -63,7 +62,7 @@ namespace Opm {
* In a parallel run only cells owned contribute to the cell average.
* \tparam is_parallel Whether this is a parallel run.
*/
template<bool is_parallel>
template<class Scalar, bool is_parallel>
struct AverageIncrementCalculator
{
/**
@@ -77,12 +76,12 @@ namespace Opm {
* by this process (value 1), or not (value 0).
* \param cell The cell index.
*/
std::tuple<double, double, double, double, int>
operator()(const std::vector<double>& pressure,
const std::vector<double>& temperature,
const std::vector<double>& rs,
const std::vector<double>& rv,
const std::vector<double>& ownership,
std::tuple<Scalar, Scalar, Scalar, Scalar, int>
operator()(const std::vector<Scalar>& pressure,
const std::vector<Scalar>& temperature,
const std::vector<Scalar>& rs,
const std::vector<Scalar>& rv,
const std::vector<Scalar>& ownership,
std::size_t cell){
if ( ownership[cell] )
{
@@ -98,15 +97,15 @@ namespace Opm {
}
}
};
template<>
struct AverageIncrementCalculator<false>
template<class Scalar>
struct AverageIncrementCalculator<Scalar, false>
{
std::tuple<double, double, double, double, int>
operator()(const std::vector<double>& pressure,
const std::vector<double>& temperature,
const std::vector<double>& rs,
const std::vector<double>& rv,
const std::vector<double>&,
std::tuple<Scalar, Scalar, Scalar, Scalar, int>
operator()(const std::vector<Scalar>& pressure,
const std::vector<Scalar>& temperature,
const std::vector<Scalar>& rs,
const std::vector<Scalar>& rv,
const std::vector<Scalar>&,
std::size_t cell){
return std::make_tuple(pressure[cell],
temperature[cell],

30
opm/simulators/wells/RegionAverageCalculator.hpp
View File

@@ -29,13 +29,7 @@
#include <algorithm>
#include <cassert>
#include <cmath>
#include <memory>
#include <stdexcept>
#include <type_traits>
#include <unordered_map>
#include <utility>
#include <vector>
/**
* \file
@@ -64,6 +58,7 @@ namespace Opm {
template <class FluidSystem, class Region>
class AverageRegionalPressure {
public:
using Scalar = typename FluidSystem::Scalar;
/**
* Constructor.
*
@@ -127,12 +122,12 @@ namespace Opm {
const auto& intQuants = elemCtx.intensiveQuantities(/*spaceIdx=*/0, /*timeIdx=*/0);
const auto& fs = intQuants.fluidState();
// use pore volume weighted averages.
const double pv_cell =
const Scalar pv_cell =
simulator.model().dofTotalVolume(cellIdx)
* intQuants.porosity().value();
// only count oil and gas filled parts of the domain
double hydrocarbon = 1.0;
Scalar hydrocarbon = 1.0;
const auto& pu = phaseUsage_;
if (RegionAttributeHelpers::PhaseUsed::water(pu)) {
hydrocarbon -= fs.saturation(FluidSystem::waterPhaseIdx).value();
@@ -142,7 +137,7 @@ namespace Opm {
assert(reg >= 0);
// sum p, rs, rv, and T.
const double hydrocarbonPV = pv_cell*hydrocarbon;
const Scalar hydrocarbonPV = pv_cell*hydrocarbon;
if (hydrocarbonPV > 0.) {
auto& attr = attributes_hpv[reg];
attr.pv += hydrocarbonPV;
@@ -171,19 +166,19 @@ namespace Opm {
for (int reg = 1; reg <= numRegions ; ++ reg) {
auto& ra = attr_.attributes(reg);
const double hpv_sum = comm.sum(attributes_hpv[reg].pv);
const Scalar hpv_sum = comm.sum(attributes_hpv[reg].pv);
// TODO: should we have some epsilon here instead of zero?
if (hpv_sum > 0.) {
const auto& attri_hpv = attributes_hpv[reg];
const double p_hpv_sum = comm.sum(attri_hpv.pressure);
const Scalar p_hpv_sum = comm.sum(attri_hpv.pressure);
ra.pressure = p_hpv_sum / hpv_sum;
} else {
// using the pore volume to do the averaging
const auto& attri_pv = attributes_pv[reg];
const double pv_sum = comm.sum(attri_pv.pv);
const Scalar pv_sum = comm.sum(attri_pv.pv);
// pore volums can be zero if a fipnum region is empty
if (pv_sum > 0) {
const double p_pv_sum = comm.sum(attri_pv.pressure);
const Scalar p_pv_sum = comm.sum(attri_pv.pressure);
ra.pressure = p_pv_sum / pv_sum;
}
}
@@ -201,12 +196,12 @@ namespace Opm {
* Average pressure
*
*/
double
Scalar
pressure(const RegionId r) const
{
if (r == 0 ) // region 0 is the whole field
{
double pressure = 0.0;
Scalar pressure = 0.0;
int num_active_regions = 0;
for (const auto& attr : attr_.attributes()) {
const auto& value = *attr.second;
@@ -243,9 +238,8 @@ namespace Opm {
{}
double pressure;
double pv;
Scalar pressure;
Scalar pv;
};
RegionAttributeHelpers::RegionAttributes<RegionId, Attributes> attr_;

8
opm/simulators/wells/SingleWellState.cpp
View File

@@ -31,10 +31,10 @@ namespace Opm {
template<class Scalar>
SingleWellState<Scalar>::
SingleWellState(const std::string& name_,
const ParallelWellInfo& pinfo,
const ParallelWellInfo<Scalar>& pinfo,
bool is_producer,
Scalar pressure_first_connection,
const std::vector<PerforationData>& perf_input,
const std::vector<PerforationData<Scalar>>& perf_input,
const PhaseUsage& pu_,
Scalar temp)
: name(name_)
@@ -63,7 +63,7 @@ SingleWellState(const std::string& name_,
template<class Scalar>
SingleWellState<Scalar> SingleWellState<Scalar>::
serializationTestObject(const ParallelWellInfo& pinfo)
serializationTestObject(const ParallelWellInfo<Scalar>& pinfo)
{
SingleWellState result("testing", pinfo, true, 1.0, {}, PhaseUsage{}, 2.0);
result.perf_data = PerfData<Scalar>::serializationTestObject();
@@ -138,7 +138,7 @@ void SingleWellState<Scalar>::updateStatus(Well::Status new_status)
template<class Scalar>
void SingleWellState<Scalar>::
reset_connection_factors(const std::vector<PerforationData>& new_perf_data)
reset_connection_factors(const std::vector<PerforationData<Scalar>>& new_perf_data)
{
if (this->perf_data.size() != new_perf_data.size()) {
throw std::invalid_argument {

12
opm/simulators/wells/SingleWellState.hpp
View File

@@ -33,7 +33,7 @@
namespace Opm {
struct PerforationData;
template<class Scalar> struct PerforationData;
class SummaryState;
class Well;
@@ -41,14 +41,14 @@ template<class Scalar>
class SingleWellState {
public:
SingleWellState(const std::string& name,
const ParallelWellInfo& pinfo,
const ParallelWellInfo<Scalar>& pinfo,
bool is_producer,
Scalar presssure_first_connection,
const std::vector<PerforationData>& perf_input,
const std::vector<PerforationData<Scalar>>& perf_input,
const PhaseUsage& pu,
Scalar temp);
static SingleWellState serializationTestObject(const ParallelWellInfo& pinfo);
static SingleWellState serializationTestObject(const ParallelWellInfo<Scalar>& pinfo);
template<class Serializer>
void serializeOp(Serializer& serializer)
@@ -79,7 +79,7 @@ public:
bool operator==(const SingleWellState&) const;
std::string name;
std::reference_wrapper<const ParallelWellInfo> parallel_info;
std::reference_wrapper<const ParallelWellInfo<Scalar>> parallel_info;
WellStatus status{WellStatus::OPEN};
bool producer;
@@ -120,7 +120,7 @@ public:
/// \param[in] new_perf_data New perforation data. Only
/// PerforationData::connection_transmissibility_factor actually
/// used (overwrites existing internal values).
void reset_connection_factors(const std::vector<PerforationData>& new_perf_data);
void reset_connection_factors(const std::vector<PerforationData<Scalar>>& new_perf_data);
void update_producer_targets(const Well& ecl_well, const SummaryState& st);
void update_injector_targets(const Well& ecl_well, const SummaryState& st);
void update_targets(const Well& ecl_well, const SummaryState& st);

16
opm/simulators/wells/StandardWell.hpp
View File

@@ -120,7 +120,7 @@ namespace Opm
using BVectorWell = typename StdWellEval::BVectorWell;
StandardWell(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const ModelParameters& param,
const RateConverterType& rate_converter,
@@ -128,7 +128,7 @@ namespace Opm
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data);
const std::vector<PerforationData<Scalar>>& perf_data);
virtual void init(const PhaseUsage* phase_usage_arg,
const std::vector<Scalar>& depth_arg,
@@ -291,7 +291,7 @@ namespace Opm
const int perf,
const bool allow_cf,
std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const;
template<class Value>
@@ -309,7 +309,7 @@ namespace Opm
const Value& skin_pressure,
const std::vector<Value>& cmix_s,
std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const;
void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator,
@@ -469,7 +469,7 @@ namespace Opm
template<class Value>
void gasOilPerfRateInj(const std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rv,
const Value& rs,
const Value& pressure,
@@ -478,20 +478,20 @@ namespace Opm
template<class Value>
void gasOilPerfRateProd(std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rv,
const Value& rs,
const Value& rvw) const;
template<class Value>
void gasWaterPerfRateProd(std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rvw,
const Value& rsw) const;
template<class Value>
void gasWaterPerfRateInj(const std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rvw,
const Value& rsw,
const Value& pressure,

2
opm/simulators/wells/StandardWellEquations.cpp
View File

@@ -42,7 +42,7 @@ namespace Opm
template<class Scalar, int numEq>
StandardWellEquations<Scalar,numEq>::
StandardWellEquations(const ParallelWellInfo& parallel_well_info)
StandardWellEquations(const ParallelWellInfo<Scalar>& parallel_well_info)
: parallelB_(duneB_, parallel_well_info)
{
duneB_.setBuildMode(OffDiagMatWell::row_wise);

4
opm/simulators/wells/StandardWellEquations.hpp
View File

@@ -34,7 +34,7 @@
namespace Opm
{
class ParallelWellInfo;
template<class Scalar> class ParallelWellInfo;
template<class Scalar, int numEq> class StandardWellEquationAccess;
#if COMPILE_BDA_BRIDGE
class WellContributions;
@@ -65,7 +65,7 @@ public:
// block vector type
using BVector = Dune::BlockVector<Dune::FieldVector<Scalar,numEq>>;
StandardWellEquations(const ParallelWellInfo& parallel_well_info);
StandardWellEquations(const ParallelWellInfo<Scalar>& parallel_well_info);
//! \brief Setup sparsity pattern for the matrices.
//! \param num_cells Total number of cells

26
opm/simulators/wells/StandardWell_impl.hpp
View File

@@ -70,7 +70,7 @@ namespace Opm
template<typename TypeTag>
StandardWell<TypeTag>::
StandardWell(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const ModelParameters& param,
const RateConverterType& rate_converter,
@@ -78,7 +78,7 @@ namespace Opm
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data)
const std::vector<PerforationData<Scalar>>& perf_data)
: Base(well, pw_info, time_step, param, rate_converter, pvtRegionIdx, num_components, num_phases, index_of_well, perf_data)
, StdWellEval(static_cast<const WellInterfaceIndices<FluidSystem,Indices>&>(*this))
, regularize_(false)
@@ -130,7 +130,7 @@ namespace Opm
const int perf,
const bool allow_cf,
std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const
{
auto obtain = [this](const Eval& value)
@@ -237,7 +237,7 @@ namespace Opm
const Value& skin_pressure,
const std::vector<Value>& cmix_s,
std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
DeferredLogger& deferred_logger) const
{
// Pressure drawdown (also used to determine direction of flow)
@@ -502,8 +502,8 @@ namespace Opm
std::vector<EvalWell> mob(this->num_components_, {this->primary_variables_.numWellEq() + Indices::numEq, 0.});
getMobility(simulator, perf, mob, deferred_logger);
PerforationRates perf_rates;
Scalar trans_mult = simulator.problem().template wellTransMultiplier<double>(intQuants, cell_idx);
PerforationRates<Scalar> perf_rates;
Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(intQuants, cell_idx);
const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
const std::vector<Scalar> Tw = this->wellIndex(perf, intQuants, trans_mult, wellstate_nupcol);
computePerfRate(intQuants, mob, bhp, Tw, perf, allow_cf,
@@ -1469,7 +1469,7 @@ namespace Opm
const std::vector<Scalar> Tw = this->wellIndex(perf, intQuants, trans_mult, wellstate_nupcol);
std::vector<Scalar> cq_s(this->num_components_, 0.);
PerforationRates perf_rates;
PerforationRates<Scalar> perf_rates;
computePerfRate(intQuants, mob, bhp, Tw, perf, allow_cf,
cq_s, perf_rates, deferred_logger);
@@ -1838,7 +1838,7 @@ namespace Opm
const EvalWell& bhp = this->primary_variables_.eval(Bhp);
std::vector<EvalWell> cq_s(this->num_components_, {this->primary_variables_.numWellEq() + Indices::numEq, 0.});
PerforationRates perf_rates;
PerforationRates<Scalar> perf_rates;
Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(int_quant, cell_idx);
const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
const std::vector<Scalar> Tw = this->wellIndex(perf, int_quant, trans_mult, wellstate_nupcol);
@@ -2463,7 +2463,7 @@ namespace Opm
Scalar trans_mult = simulator.problem().template wellTransMultiplier<Scalar>(intQuants, cell_idx);
const auto& wellstate_nupcol = simulator.problem().wellModel().nupcolWellState().well(this->index_of_well_);
const std::vector<Scalar> Tw = this->wellIndex(perf, intQuants, trans_mult, wellstate_nupcol);
PerforationRates perf_rates;
PerforationRates<Scalar> perf_rates;
computePerfRate(intQuants, mob, bhp.value(), Tw, perf, allow_cf,
cq_s, perf_rates, deferred_logger);
for (int comp = 0; comp < this->num_components_; ++comp) {
@@ -2596,7 +2596,7 @@ namespace Opm
void
StandardWell<TypeTag>::
gasOilPerfRateInj(const std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rv,
const Value& rs,
const Value& pressure,
@@ -2644,7 +2644,7 @@ namespace Opm
void
StandardWell<TypeTag>::
gasOilPerfRateProd(std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rv,
const Value& rs,
const Value& rvw) const
@@ -2680,7 +2680,7 @@ namespace Opm
void
StandardWell<TypeTag>::
gasWaterPerfRateProd(std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rvw,
const Value& rsw) const
{
@@ -2704,7 +2704,7 @@ namespace Opm
void
StandardWell<TypeTag>::
gasWaterPerfRateInj(const std::vector<Value>& cq_s,
PerforationRates& perf_rates,
PerforationRates<Scalar>& perf_rates,
const Value& rvw,
const Value& rsw,
const Value& pressure,

2
opm/simulators/wells/WGState.cpp
View File

@@ -35,7 +35,7 @@ WGState<Scalar>::WGState(const PhaseUsage& pu) :
template<class Scalar>
WGState<Scalar> WGState<Scalar>::
serializationTestObject(const ParallelWellInfo& pinfo)
serializationTestObject(const ParallelWellInfo<Scalar>& pinfo)
{
WGState result(PhaseUsage{});
result.well_state = WellState<Scalar>::serializationTestObject(pinfo);

4
opm/simulators/wells/WGState.hpp
View File

@@ -26,7 +26,7 @@
namespace Opm {
class ParallelWellInfo;
template<class Scalar> class ParallelWellInfo;
/*
Microscopic class to handle well, group and well test state.
@@ -39,7 +39,7 @@ struct WGState
{
WGState(const PhaseUsage& pu);
static WGState serializationTestObject(const ParallelWellInfo& pinfo);
static WGState serializationTestObject(const ParallelWellInfo<Scalar>& pinfo);
void wtest_state(WellTestState wtest_state);

1
opm/simulators/wells/WellFilterCake.cpp
View File

@@ -25,6 +25,7 @@
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/PerforationData.hpp>
#include <opm/simulators/wells/WellInterfaceGeneric.hpp>
#include <opm/simulators/wells/WellState.hpp>

3
opm/simulators/wells/WellHelpers.cpp
View File

@@ -41,7 +41,8 @@ namespace wellhelpers {
template<typename Scalar>
ParallelStandardWellB<Scalar>::
ParallelStandardWellB(const Matrix& B, const ParallelWellInfo& parallel_well_info)
ParallelStandardWellB(const Matrix& B,
const ParallelWellInfo<Scalar>& parallel_well_info)
: B_(B), parallel_well_info_(parallel_well_info)
{}

7
opm/simulators/wells/WellHelpers.hpp
View File

@@ -28,7 +28,7 @@
namespace Opm {
class ParallelWellInfo;
template<class Scalar> class ParallelWellInfo;
struct WellProductionControls;
struct WellInjectionControls;
enum class WellProducerCMode;
@@ -53,7 +53,8 @@ public:
using Block = Dune::DynamicMatrix<Scalar>;
using Matrix = Dune::BCRSMatrix<Block>;
ParallelStandardWellB(const Matrix& B, const ParallelWellInfo& parallel_well_info);
ParallelStandardWellB(const Matrix& B,
const ParallelWellInfo<Scalar>& parallel_well_info);
//! y = A x
template<class X, class Y>
@@ -65,7 +66,7 @@ public:
private:
const Matrix& B_;
const ParallelWellInfo& parallel_well_info_;
const ParallelWellInfo<Scalar>& parallel_well_info_;
};
template<class Scalar>

5
opm/simulators/wells/WellInterface.hpp
View File

@@ -60,7 +60,6 @@ namespace Opm {
#include <opm/material/densead/Evaluation.hpp>
#include <cassert>
#include <vector>
namespace Opm
@@ -134,7 +133,7 @@ public:
Indices::numPhases >;
/// Constructor
WellInterface(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const ModelParameters& param,
const RateConverterType& rate_converter,
@@ -142,7 +141,7 @@ public:
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data);
const std::vector<PerforationData<Scalar>>& perf_data);
/// Virtual destructor
virtual ~WellInterface() = default;

4
opm/simulators/wells/WellInterfaceFluidSystem.cpp
View File

@@ -46,14 +46,14 @@ namespace Opm
template<class FluidSystem>
WellInterfaceFluidSystem<FluidSystem>::
WellInterfaceFluidSystem(const Well& well,
const ParallelWellInfo& parallel_well_info,
const ParallelWellInfo<Scalar>& parallel_well_info,
const int time_step,
const RateConverterType& rate_converter,
const int pvtRegionIdx,
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data)
const std::vector<PerforationData<Scalar>>& perf_data)
: WellInterfaceGeneric<Scalar>(well, parallel_well_info, time_step,
pvtRegionIdx, num_components, num_phases,
index_of_well, perf_data)

4
opm/simulators/wells/WellInterfaceFluidSystem.hpp
View File

@@ -68,14 +68,14 @@ public:
protected:
WellInterfaceFluidSystem(const Well& well,
const ParallelWellInfo& parallel_well_info,
const ParallelWellInfo<Scalar>& parallel_well_info,
const int time_step,
const RateConverterType& rate_converter,
const int pvtRegionIdx,
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data);
const std::vector<PerforationData<Scalar>>& perf_data);
// updating the voidage rates in well_state when requested
void calculateReservoirRates(SingleWellState<Scalar>& ws) const;

6
opm/simulators/wells/WellInterfaceGeneric.cpp
View File

@@ -54,13 +54,13 @@ namespace Opm {
template<class Scalar>
WellInterfaceGeneric<Scalar>::
WellInterfaceGeneric(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const int pvtRegionIdx,
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data)
const std::vector<PerforationData<Scalar>>& perf_data)
: well_ecl_(well)
, parallel_well_info_(pw_info)
, current_step_(time_step)
@@ -143,7 +143,7 @@ adaptRatesForVFP(std::vector<Scalar>& rates) const
}
template<class Scalar>
const std::vector<PerforationData>&
const std::vector<PerforationData<Scalar>>&
WellInterfaceGeneric<Scalar>::perforationData() const
{
return *perf_data_;

16
opm/simulators/wells/WellInterfaceGeneric.hpp
View File

@@ -36,8 +36,8 @@ namespace Opm
class DeferredLogger;
class GuideRate;
class ParallelWellInfo;
struct PerforationData;
template<class Scalar> class ParallelWellInfo;
template<class Scalar> struct PerforationData;
struct PhaseUsage;
class SummaryState;
template<class Scalar> class VFPProperties;
@@ -51,16 +51,16 @@ template<class Scalar>
class WellInterfaceGeneric {
public:
WellInterfaceGeneric(const Well& well,
const ParallelWellInfo& parallel_well_info,
const ParallelWellInfo<Scalar>& parallel_well_info,
const int time_step,
const int pvtRegionIdx,
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data);
const std::vector<PerforationData<Scalar>>& perf_data);
/// \brief Get the perforations of the well
const std::vector<PerforationData>& perforationData() const;
const std::vector<PerforationData<Scalar>>& perforationData() const;
/// Well name.
const std::string& name() const;
@@ -131,7 +131,7 @@ public:
const VFPProperties<Scalar>* vfpProperties() const { return vfp_properties_; }
const ParallelWellInfo& parallelWellInfo() const { return parallel_well_info_; }
const ParallelWellInfo<Scalar>& parallelWellInfo() const { return parallel_well_info_; }
const std::vector<Scalar>& perfDepth() const { return perf_depth_; }
@@ -271,7 +271,7 @@ protected:
Well well_ecl_;
const ParallelWellInfo& parallel_well_info_;
const ParallelWellInfo<Scalar>& parallel_well_info_;
const int current_step_;
// The pvt region of the well. We assume
@@ -286,7 +286,7 @@ protected:
// the index of well in Wells struct
int index_of_well_;
const std::vector<PerforationData>* perf_data_;
const std::vector<PerforationData<Scalar>>* perf_data_;
// the vectors used to describe the inflow performance relationship (IPR)
// Q = IPR_A - BHP * IPR_B

4
opm/simulators/wells/WellInterfaceIndices.cpp
View File

@@ -36,14 +36,14 @@ namespace Opm
template<class FluidSystem, class Indices>
WellInterfaceIndices<FluidSystem,Indices>::
WellInterfaceIndices(const Well& well,
const ParallelWellInfo& parallel_well_info,
const ParallelWellInfo<Scalar>& parallel_well_info,
const int time_step,
const typename WellInterfaceFluidSystem<FluidSystem>::RateConverterType& rate_converter,
const int pvtRegionIdx,
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data)
const std::vector<PerforationData<Scalar>>& perf_data)
: WellInterfaceFluidSystem<FluidSystem>(well,
parallel_well_info,
time_step,

4
opm/simulators/wells/WellInterfaceIndices.hpp
View File

@@ -56,14 +56,14 @@ public:
protected:
WellInterfaceIndices(const Well& well,
const ParallelWellInfo& parallel_well_info,
const ParallelWellInfo<Scalar>& parallel_well_info,
const int time_step,
const typename WellInterfaceFluidSystem<FluidSystem>::RateConverterType& rate_converter,
const int pvtRegionIdx,
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data);
const std::vector<PerforationData<Scalar>>& perf_data);
};
}

5
opm/simulators/wells/WellInterface_impl.hpp
View File

@@ -40,6 +40,7 @@
#include <dune/common/version.hh>
#include <cassert>
#include <cstddef>
#include <utility>
@@ -52,7 +53,7 @@ namespace Opm
template<typename TypeTag>
WellInterface<TypeTag>::
WellInterface(const Well& well,
const ParallelWellInfo& pw_info,
const ParallelWellInfo<Scalar>& pw_info,
const int time_step,
const ModelParameters& param,
const RateConverterType& rate_converter,
@@ -60,7 +61,7 @@ namespace Opm
const int num_components,
const int num_phases,
const int index_of_well,
const std::vector<PerforationData>& perf_data)
const std::vector<PerforationData<Scalar>>& perf_data)
: WellInterfaceIndices<FluidSystem,Indices>(well,
pw_info,
time_step,

31
opm/simulators/wells/WellState.cpp
View File

@@ -30,6 +30,7 @@
#include <opm/simulators/utils/ParallelCommunication.hpp>
#include <opm/simulators/wells/ParallelWellInfo.hpp>
#include <opm/simulators/wells/PerforationData.hpp>
#include <opm/grid/common/p2pcommunicator.hh>
#include <opm/output/data/Wells.hpp>
@@ -127,7 +128,7 @@ void PackUnpackXConn::unpack([[maybe_unused]] const int link,
namespace Opm {
template<class Scalar>
WellState<Scalar>::WellState(const ParallelWellInfo& pinfo)
WellState<Scalar>::WellState(const ParallelWellInfo<Scalar>& pinfo)
: phase_usage_{}
{
wells_.add("test4",
@@ -136,7 +137,7 @@ WellState<Scalar>::WellState(const ParallelWellInfo& pinfo)
template<class Scalar>
WellState<Scalar> WellState<Scalar>::
serializationTestObject(const ParallelWellInfo& pinfo)
serializationTestObject(const ParallelWellInfo<Scalar>& pinfo)
{
WellState result(PhaseUsage{});
result.alq_state = ALQState<Scalar>::serializationTestObject();
@@ -149,8 +150,8 @@ serializationTestObject(const ParallelWellInfo& pinfo)
template<class Scalar>
void WellState<Scalar>::base_init(const std::vector<Scalar>& cellPressures,
const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
const SummaryState& summary_state)
{
// clear old name mapping
@@ -170,9 +171,9 @@ void WellState<Scalar>::base_init(const std::vector<Scalar>& cellPressures,
template<class Scalar>
void WellState<Scalar>::initSingleProducer(const Well& well,
const ParallelWellInfo& well_info,
const ParallelWellInfo<Scalar>& well_info,
Scalar pressure_first_connection,
const std::vector<PerforationData>& well_perf_data,
const std::vector<PerforationData<Scalar>>& well_perf_data,
const SummaryState& summary_state)
{
const auto& pu = this->phase_usage_;
@@ -199,9 +200,9 @@ void WellState<Scalar>::initSingleProducer(const Well& well,
template<class Scalar>
void WellState<Scalar>::initSingleInjector(const Well& well,
const ParallelWellInfo& well_info,
const ParallelWellInfo<Scalar>& well_info,
Scalar pressure_first_connection,
const std::vector<PerforationData>& well_perf_data,
const std::vector<PerforationData<Scalar>>& well_perf_data,
const SummaryState& summary_state)
{
const auto& pu = this->phase_usage_;
@@ -228,8 +229,8 @@ void WellState<Scalar>::initSingleInjector(const Well& well,
template<class Scalar>
void WellState<Scalar>::initSingleWell(const std::vector<Scalar>& cellPressures,
const Well& well,
const std::vector<PerforationData>& well_perf_data,
const ParallelWellInfo& well_info,
const std::vector<PerforationData<Scalar>>& well_perf_data,
const ParallelWellInfo<Scalar>& well_info,
const SummaryState& summary_state)
{
Scalar pressure_first_connection = -1;
@@ -250,10 +251,10 @@ template<class Scalar>
void WellState<Scalar>::init(const std::vector<Scalar>& cellPressures,
const Schedule& schedule,
const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
const int report_step,
const WellState* prevState,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
const SummaryState& summary_state)
{
// call init on base class
@@ -419,11 +420,11 @@ void WellState<Scalar>::init(const std::vector<Scalar>& cellPressures,
template<class Scalar>
void WellState<Scalar>::resize(const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
const Schedule& schedule,
const bool handle_ms_well,
const std::size_t numCells,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
const SummaryState& summary_state)
{
const std::vector<Scalar> tmp(numCells, 0.0); // <- UGLY HACK to pass the size
@@ -1043,7 +1044,7 @@ bool WellState<Scalar>::operator==(const WellState& rhs) const
}
template<class Scalar>
const ParallelWellInfo&
const ParallelWellInfo<Scalar>&
WellState<Scalar>::parallelWellInfo(std::size_t well_index) const
{
const auto& ws = this->well(well_index);

34
opm/simulators/wells/WellState.hpp
View File

@@ -36,7 +36,6 @@
#include <opm/simulators/wells/ALQState.hpp>
#include <opm/simulators/wells/GlobalWellInfo.hpp>
#include <opm/simulators/wells/PerfData.hpp>
#include <opm/simulators/wells/PerforationData.hpp>
#include <opm/simulators/wells/SegmentState.hpp>
#include <opm/simulators/wells/SingleWellState.hpp>
#include <opm/simulators/wells/WellContainer.hpp>
@@ -51,7 +50,8 @@
namespace Opm
{
class ParallelWellInfo;
template<class Scalar> class ParallelWellInfo;
template<class Scalar> struct PerforationData;
class Schedule;
enum class WellStatus;
@@ -68,13 +68,13 @@ public:
static const int Gas = BlackoilPhases::Vapour;
// Only usable for testing purposes
explicit WellState(const ParallelWellInfo& pinfo);
explicit WellState(const ParallelWellInfo<Scalar>& pinfo);
explicit WellState(const PhaseUsage& pu)
: phase_usage_(pu)
{}
static WellState serializationTestObject(const ParallelWellInfo& pinfo);
static WellState serializationTestObject(const ParallelWellInfo<Scalar>& pinfo);
std::size_t size() const
{
@@ -91,7 +91,7 @@ public:
return this->size();
}
const ParallelWellInfo& parallelWellInfo(std::size_t well_index) const;
const ParallelWellInfo<Scalar>& parallelWellInfo(std::size_t well_index) const;
/// Allocate and initialize if wells is non-null. Also tries
/// to give useful initial values to the bhp(), wellRates()
@@ -99,18 +99,18 @@ public:
void init(const std::vector<Scalar>& cellPressures,
const Schedule& schedule,
const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
const int report_step,
const WellState* prevState,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
const SummaryState& summary_state);
void resize(const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
const Schedule& schedule,
const bool handle_ms_well,
const std::size_t numCells,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
const SummaryState& summary_state);
void setCurrentWellRates(const std::string& wellName,
@@ -378,26 +378,26 @@ private:
/// with -1e100.
void base_init(const std::vector<Scalar>& cellPressures,
const std::vector<Well>& wells_ecl,
const std::vector<std::reference_wrapper<ParallelWellInfo>>& parallel_well_info,
const std::vector<std::vector<PerforationData>>& well_perf_data,
const std::vector<std::reference_wrapper<ParallelWellInfo<Scalar>>>& parallel_well_info,
const std::vector<std::vector<PerforationData<Scalar>>>& well_perf_data,
const SummaryState& summary_state);
void initSingleWell(const std::vector<Scalar>& cellPressures,
const Well& well,
const std::vector<PerforationData>& well_perf_data,
const ParallelWellInfo& well_info,
const std::vector<PerforationData<Scalar>>& well_perf_data,
const ParallelWellInfo<Scalar>& well_info,
const SummaryState& summary_state);
void initSingleProducer(const Well& well,
const ParallelWellInfo& well_info,
const ParallelWellInfo<Scalar>& well_info,
Scalar pressure_first_connection,
const std::vector<PerforationData>& well_perf_data,
const std::vector<PerforationData<Scalar>>& well_perf_data,
const SummaryState& summary_state);
void initSingleInjector(const Well& well,
const ParallelWellInfo& well_info,
const ParallelWellInfo<Scalar>& well_info,
Scalar pressure_first_connection,
const std::vector<PerforationData>& well_perf_data,
const std::vector<PerforationData<Scalar>>& well_perf_data,
const SummaryState& summary_state);
};

8
tests/test_RestartSerialization.cpp
View File

@@ -149,7 +149,7 @@ TEST_FOR_TYPE_NAMED(BVec, BlockVectorWrapper)
BOOST_AUTO_TEST_CASE(SingleWellState)
{
Opm::ParallelWellInfo dummy;
Opm::ParallelWellInfo<double> dummy;
auto data_out = Opm::SingleWellState<double>::serializationTestObject(dummy);
Opm::Serialization::MemPacker packer;
Opm::Serializer ser(packer);
@@ -178,7 +178,7 @@ BOOST_AUTO_TEST_CASE(WellContainer)
BOOST_AUTO_TEST_CASE(WellState)
{
Opm::ParallelWellInfo dummy;
Opm::ParallelWellInfo<double> dummy;
auto data_out = Opm::WellState<double>::serializationTestObject(dummy);
Opm::Serialization::MemPacker packer;
Opm::Serializer ser(packer);
@@ -193,7 +193,7 @@ BOOST_AUTO_TEST_CASE(WellState)
BOOST_AUTO_TEST_CASE(WGState)
{
Opm::ParallelWellInfo dummy;
Opm::ParallelWellInfo<double> dummy;
auto data_out = Opm::WGState<double>::serializationTestObject(dummy);
Opm::Serialization::MemPacker packer;
Opm::Serializer ser(packer);
@@ -354,7 +354,7 @@ public:
}
private:
ParallelWellInfo dummy;
ParallelWellInfo<double> dummy;
};
}

8
tests/test_parallel_wbp_calculation.cpp
View File

@@ -457,10 +457,10 @@ namespace {
return w;
}
Opm::ParallelWellInfo
Opm::ParallelWellInfo<double>
parallelWellInfo(const Opm::Parallel::Communication& comm)
{
auto pwi = Opm::ParallelWellInfo {
auto pwi = Opm::ParallelWellInfo<double> {
std::pair { std::string{ "P" }, true }, comm
};
@@ -469,7 +469,7 @@ namespace {
const auto numLocalPerf = 3;
const auto perfOffset = comm.rank() * numLocalPerf;
auto prev = Opm::ParallelWellInfo::INVALID_ECL_INDEX;
auto prev = Opm::ParallelWellInfo<double>::INVALID_ECL_INDEX;
for (auto perf = 0; perf < numLocalPerf; ++perf) {
const auto curr = perfOffset + perf;
pwi.pushBackEclIndex(prev, curr);
@@ -552,7 +552,7 @@ namespace {
Opm::Parallel::Communication comm;
Opm::GridDims cellIndexMap;
Opm::ParallelWBPCalculation wbpCalcService;
Opm::ParallelWellInfo pwi;
Opm::ParallelWellInfo<double> pwi;
};
} // Anonymous namespace

26
tests/test_parallelwellinfo.cpp
View File

@@ -97,7 +97,7 @@ std::ostream& operator<<(std::ostream& os, const std::pair<std::string, bool>& p
}
namespace Opm
{
std::ostream& operator<<(std::ostream& os, const Opm::ParallelWellInfo& w)
std::ostream& operator<<(std::ostream& os, const Opm::ParallelWellInfo<double>& w)
{
return os << "{" << w.name() << " "<< w.hasLocalCells() << " "<<
w.isOwner() << "}";
@@ -117,7 +117,7 @@ BOOST_AUTO_TEST_CASE(ParallelWellComparison)
else
pairs = {{"Test1", false},{"Test2", true}, {"Test1", true} };
std::vector<Opm::ParallelWellInfo> well_info;
std::vector<Opm::ParallelWellInfo<double>> well_info;
for (const auto& wellinfo : pairs) {
well_info.emplace_back(wellinfo, Opm::Parallel::Communication());
@@ -138,16 +138,16 @@ BOOST_AUTO_TEST_CASE(ParallelWellComparison)
BOOST_CHECK(well_info[0] != well_info[1]);
Opm::ParallelWellInfo well0, well1;
Opm::ParallelWellInfo<double> well0, well1;
BOOST_CHECK(well0 == well1);
#if HAVE_MPI
BOOST_CHECK(well0.communication()==helper.getLocalCommunicator());
#endif
Opm::ParallelWellInfo well2("Test", false);
Opm::ParallelWellInfo<double> well2("Test", false);
std::pair<std::string, bool> pwell={"Test", true};
BOOST_CHECK(well2 < pwell);
Opm::ParallelWellInfo well3("Test", true);
Opm::ParallelWellInfo<double> well3("Test", true);
BOOST_CHECK(! (well3 < pwell));
pwell.second = false;
BOOST_CHECK(! (well3 < pwell));
@@ -171,7 +171,7 @@ BOOST_AUTO_TEST_CASE(ParallelWellComparison)
BOOST_AUTO_TEST_CASE(CommunicateAboveBelowSelf)
{
auto comm = Dune::MPIHelper::getLocalCommunicator();
Opm::CommunicateAboveBelow commAboveBelow{ comm };
Opm::CommunicateAboveBelow<double> commAboveBelow{ comm };
for(std::size_t count=0; count < 2; ++count)
{
std::vector<int> eclIndex = {0, 1, 2, 3, 7 , 8, 10, 11};
@@ -206,7 +206,7 @@ BOOST_AUTO_TEST_CASE(CommunicateAboveBelowSelf)
BOOST_AUTO_TEST_CASE(CommunicateAboveBelowSelf1)
{
auto comm = Dune::MPIHelper::getLocalCommunicator();
Opm::CommunicateAboveBelow commAboveBelow{ comm };
Opm::CommunicateAboveBelow<double> commAboveBelow{ comm };
for(std::size_t count=0; count < 2; ++count)
{
std::vector<int> eclIndex = {0};
@@ -291,7 +291,7 @@ BOOST_AUTO_TEST_CASE(CommunicateAboveBelowParallel)
{
auto comm = Opm::Parallel::Communication(Dune::MPIHelper::getCommunicator());
Opm::CommunicateAboveBelow commAboveBelow{ comm };
Opm::CommunicateAboveBelow<double> commAboveBelow{ comm };
for(std::size_t count=0; count < 2; ++count)
{
auto globalEclIndex = createGlobalEclIndex(comm);
@@ -351,7 +351,7 @@ BOOST_AUTO_TEST_CASE(PartialSumself)
{
auto comm = Dune::MPIHelper::getLocalCommunicator();
Opm::CommunicateAboveBelow commAboveBelow{ comm };
Opm::CommunicateAboveBelow<double> commAboveBelow{ comm };
std::vector<int> eclIndex = {0, 1, 2, 3, 7 , 8, 10, 11};
std::vector<double> current(eclIndex.size());
std::transform(eclIndex.begin(), eclIndex.end(), current.begin(),
@@ -383,7 +383,7 @@ BOOST_AUTO_TEST_CASE(PartialSumParallel)
auto comm = Opm::Parallel::Communication(Dune::MPIHelper::getCommunicator());
Opm::CommunicateAboveBelow commAboveBelow{ comm };
Opm::CommunicateAboveBelow<double> commAboveBelow{ comm };
auto globalEclIndex = createGlobalEclIndex(comm);
std::vector<double> globalCurrent(globalEclIndex.size());
initRandomNumbers(std::begin(globalCurrent), std::end(globalCurrent),
@@ -411,7 +411,7 @@ void testGlobalPerfFactoryParallel(int num_component, bool local_consecutive = f
{
auto comm = Opm::Parallel::Communication(Dune::MPIHelper::getCommunicator());
Opm::ParallelWellInfo wellInfo{ {"Test", true }, comm };
Opm::ParallelWellInfo<double> wellInfo{ {"Test", true }, comm };
auto globalEclIndex = createGlobalEclIndex(comm);
std::vector<double> globalCurrent(globalEclIndex.size() * num_component);
std::vector<double> globalAdd(globalEclIndex.size() * num_component);
@@ -425,7 +425,7 @@ void testGlobalPerfFactoryParallel(int num_component, bool local_consecutive = f
local_consecutive);
// A hack to get local values to add.
Opm::ParallelWellInfo dummy{ {"Test", true }, comm };
Opm::ParallelWellInfo<double> dummy{ {"Test", true }, comm };
auto localAdd = populateCommAbove(dummy, comm, globalEclIndex,
globalAdd, num_component,
local_consecutive);
@@ -476,7 +476,7 @@ BOOST_AUTO_TEST_CASE(GlobalPerfFactoryParallel1)
BOOST_AUTO_TEST_CASE(EmptyWell) {
auto comm = Opm::Parallel::Communication(Dune::MPIHelper::getCommunicator());
Opm::ParallelWellInfo pw({"WELL1", true}, comm);
Opm::ParallelWellInfo<double> pw({"WELL1", true}, comm);
pw.communicateFirstPerforation(false);
double local_p = 1;
auto global_p = pw.broadcastFirstPerforationValue(local_p);

14
tests/test_wellmodel.cpp
View File

@@ -83,7 +83,7 @@ struct SetupTest {
std::shared_ptr<Opm::Python> python;
std::unique_ptr<const Opm::Schedule> schedule;
std::unique_ptr<Opm::SummaryState> summaryState;
std::vector<std::vector<Opm::PerforationData>> well_perf_data;
std::vector<std::vector<Opm::PerforationData<double>>> well_perf_data;
int current_timestep;
};
@@ -125,13 +125,13 @@ BOOST_AUTO_TEST_CASE(TestStandardWellInput) {
rateConverter.reset(new RateConverterType (phaseUsage,
std::vector<int>(10, 0)));
Opm::PerforationData dummy;
std::vector<Opm::PerforationData> pdata(well.getConnections().size(), dummy);
Opm::PerforationData<double> dummy;
std::vector<Opm::PerforationData<double>> pdata(well.getConnections().size(), dummy);
for (auto c = 0*pdata.size(); c < pdata.size(); ++c) {
pdata[c].ecl_index = c;
}
Opm::ParallelWellInfo pinfo{well.name()};
Opm::ParallelWellInfo<double> pinfo{well.name()};
BOOST_CHECK_THROW( StandardWell( well, pinfo, -1, param, *rateConverter, 0, 3, 3, 0, pdata), std::invalid_argument);
}
@@ -160,13 +160,13 @@ BOOST_AUTO_TEST_CASE(TestBehavoir) {
// Compute reservoir volumes for RESV controls.
rateConverter.reset(new RateConverterType (phaseUsage,
std::vector<int>(10, 0)));
Opm::PerforationData dummy;
std::vector<Opm::PerforationData> pdata(wells_ecl[w].getConnections().size(), dummy);
Opm::PerforationData<double> dummy;
std::vector<Opm::PerforationData<double>> pdata(wells_ecl[w].getConnections().size(), dummy);
for (auto c = 0*pdata.size(); c < pdata.size(); ++c) {
pdata[c].ecl_index = c;
}
Opm::ParallelWellInfo pinfo{wells_ecl[w].name()};
Opm::ParallelWellInfo<double> pinfo{wells_ecl[w].name()};
wells.emplace_back(new StandardWell(wells_ecl[w], pinfo, current_timestep, param, *rateConverter, 0, 3, 3, w, pdata) );
}
}

23
tests/test_wellstate.cpp
View File

@@ -28,6 +28,7 @@
#include <opm/common/ErrorMacros.hpp>
#include <opm/simulators/wells/GlobalWellInfo.hpp>
#include <opm/simulators/wells/ParallelWellInfo.hpp>
#include <opm/simulators/wells/PerforationData.hpp>
#include <opm/simulators/wells/WellState.hpp>
#include <opm/simulators/wells/SingleWellState.hpp>
#include <opm/simulators/wells/SegmentState.hpp>
@@ -106,7 +107,7 @@ struct Setup
+ std::to_string(k) + " not found in grid (well = " + well.name() + ").");
OPM_THROW(std::runtime_error, msg);
} else {
Opm::PerforationData pd;
Opm::PerforationData<double> pd;
pd.cell_index = active_index;
pd.connection_transmissibility_factor = completion.CF();
pd.connection_d_factor = completion.dFactor();
@@ -132,13 +133,13 @@ struct Setup
std::shared_ptr<Opm::Python> python;
Opm::Schedule sched;
Opm::SummaryState st;
std::vector<std::vector<Opm::PerforationData>> well_perf_data;
std::vector<std::vector<Opm::PerforationData<double>>> well_perf_data;
};
namespace {
Opm::WellState<double>
buildWellState(const Setup& setup, const std::size_t timeStep,
std::vector<Opm::ParallelWellInfo>& pinfos)
std::vector<Opm::ParallelWellInfo<double>>& pinfos)
{
auto state = Opm::WellState<double>{setup.pu};
@@ -148,7 +149,7 @@ namespace {
auto wells = setup.sched.getWells(timeStep);
pinfos.resize(wells.size());
std::vector<std::reference_wrapper<Opm::ParallelWellInfo>> ppinfos;
std::vector<std::reference_wrapper<Opm::ParallelWellInfo<double>>> ppinfos;
auto pw = pinfos.begin();
for (const auto& well : wells)
@@ -258,7 +259,7 @@ BOOST_AUTO_TEST_CASE(Linearisation)
const Setup setup{ "msw.data" };
const auto tstep = std::size_t{0};
std::vector<Opm::ParallelWellInfo> pinfos;
std::vector<Opm::ParallelWellInfo<double>> pinfos;
const auto wstate = buildWellState(setup, tstep, pinfos);
const auto& ws = wstate.well("PROD01");
@@ -275,7 +276,7 @@ BOOST_AUTO_TEST_CASE(Pressure)
const Setup setup{ "msw.data" };
const auto tstep = std::size_t{0};
std::vector<Opm::ParallelWellInfo> pinfos;
std::vector<Opm::ParallelWellInfo<double>> pinfos;
auto wstate = buildWellState(setup, tstep, pinfos);
const auto& wells = setup.sched.getWells(tstep);
@@ -314,7 +315,7 @@ BOOST_AUTO_TEST_CASE(Rates)
const Setup setup{ "msw.data" };
const auto tstep = std::size_t{0};
std::vector<Opm::ParallelWellInfo> pinfos;
std::vector<Opm::ParallelWellInfo<double>> pinfos;
auto wstate = buildWellState(setup, tstep, pinfos);
const auto wells = setup.sched.getWells(tstep);
@@ -367,7 +368,7 @@ BOOST_AUTO_TEST_CASE(STOP_well)
*/
const Setup setup{ "wells_manager_data_wellSTOP.data" };
std::vector<Opm::ParallelWellInfo> pinfos;
std::vector<Opm::ParallelWellInfo<double>> pinfos;
auto wstate = buildWellState(setup, 0, pinfos);
for (std::size_t well_index = 0; well_index < setup.sched.numWells(0); well_index++) {
const auto& ws = wstate.well(well_index);
@@ -524,7 +525,7 @@ BOOST_AUTO_TEST_CASE(TESTSegmentState) {
BOOST_AUTO_TEST_CASE(TESTSegmentState2) {
const Setup setup{ "msw.data" };
std::vector<Opm::ParallelWellInfo> pinfo;
std::vector<Opm::ParallelWellInfo<double>> pinfo;
const auto wstate = buildWellState(setup, 0, pinfo);
const auto& well = setup.sched.getWell("PROD01", 0);
const auto& ws = wstate.well("PROD01");
@@ -580,8 +581,8 @@ BOOST_AUTO_TEST_CASE(TESTPerfData) {
BOOST_AUTO_TEST_CASE(TestSingleWellState) {
Opm::ParallelWellInfo pinfo;
std::vector<Opm::PerforationData> connections = {{0,1,1,0,0},{1,1,1,0,1},{2,1,1,0,2}};
Opm::ParallelWellInfo<double> pinfo;
std::vector<Opm::PerforationData<double>> connections = {{0,1,1,0,0},{1,1,1,0,1},{2,1,1,0,2}};
Opm::PhaseUsage pu;
// This is totally bonkers, but the pu needs a complete deck to initialize properly