more or less getting the old implementation in

and compile now.
This commit is contained in:
Kai Bao 2017-06-15 17:19:49 +02:00
parent 0cf6699591
commit 2d02503091
8 changed files with 116 additions and 19 deletions

View File

@ -50,7 +50,7 @@ list (APPEND MAIN_SOURCE_FILES
opm/autodiff/WellMultiSegment.cpp
opm/autodiff/MultisegmentWells.cpp
opm/autodiff/WellInterface.cpp
opm/autodiff/StandardWell.hpp
opm/autodiff/StandardWell.cpp
opm/autodiff/MissingFeatures.cpp
opm/polymer/PolymerState.cpp
opm/polymer/PolymerBlackoilState.cpp

View File

@ -27,11 +27,11 @@
namespace Opm
{
StandardWell::
StandardWell(const Well* well, const size_t time_step, const Wells* wells)
StandardWell(const Well* well, const int time_step, const Wells* wells)
: WellInterface(well, time_step, wells)
, perf_densities_(number_of_perforations_)
, perf_pressure_diffs_(number_of_perforations_)
, well_variables_(blocksize) // the number of the primary variables
, well_variables_(numWellEq) // the number of the primary variables
{
dune_B_.setBuildMode( Mat::row_wise );
dune_C_.setBuildMode( Mat::row_wise );
@ -111,7 +111,7 @@ namespace Opm
for (int phase = 0; phase < np; ++phase) {
well_variables_[phase] = 0.0;
well_variables_[phase].setValue(well_state.wellSolutions()[index_of_well_ + nw * phase]);
well_variables_[phase].setDerivative(blocksize + phase, 1.0);
well_variables_[phase].setDerivative(numEq + phase, 1.0);
}
}
@ -182,7 +182,9 @@ namespace Opm
{
EvalWell qs = 0.0;
const WellControls* wc = well_controls_;
return qs; // temporary
/* const WellControls* wc = well_controls_;
const int np = number_of_phases_;
// the target from the well controls
@ -199,7 +201,7 @@ namespace Opm
}
if (well_controls_get_current_type(wc) == BHP || well_controls_get_current_type(wc) == THP) {
return wellVariables_[XvarWell];
return well_variables_[XvarWell];
}
// rate control
@ -211,8 +213,8 @@ namespace Opm
// Producers
if (well_controls_get_current_type(wc) == BHP || well_controls_get_current_type(wc) == THP ) {
return wellVariables_[XvarWell] * wellVolumeFractionScaled(phase);
}
return well_variables_[XvarWell] * wellVolumeFractionScaled(phase);
} */
}
@ -251,21 +253,21 @@ namespace Opm
wellVolumeFraction(const int phase) const
{
if (phase == Water) {
return wellVariables_[WFrac];
return well_variables_[WFrac];
}
if (phase == Gas) {
return wellVariables_[GFrac];
return well_variables_[GFrac];
}
// Oil fraction
EvalWell well_fraction = 1.0;
if (active_[Water]) {
well_fraction -= wellVariables_[WFrac];
if (active()[Water]) {
well_fraction -= well_variables_[WFrac];
}
if (active_[Gas]) {
well_fraction -= wellVariables_[GFrac];
if (active()[Gas]) {
well_fraction -= well_variables_[GFrac];
}
return well_fraction;
}
@ -287,6 +289,6 @@ namespace Opm
assert(sum_volume_fraction_scaled.value() != 0.);
return wellVolumeFractionScaled(phase) / sum_volume_fraction_scaled;
}
}
}

View File

@ -53,6 +53,19 @@ namespace Opm
GFrac = 2
};
typedef double Scalar;
// static const int numEq = BlackoilIndices::numEq;
static const int numEq = 3;
static const int numWellEq = numEq; //number of wellEq is the same as numEq in the model
static const int solventCompIdx = 3; //TODO get this from ebos
typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
typedef Dune::BlockVector<VectorBlockType> BVector;
typedef DenseAd::Evaluation<double, /*size=*/numEq + numWellEq> EvalWell;
typedef DenseAd::Evaluation<double, /*size=*/numEq> Eval;
// for now, using the matrix and block version in StandardWellsDense.
// TODO: for bettern generality, it should contain blocksize_field and blocksize_well.
// They are allowed to be different and it will create four types of matrix blocks and two types of
@ -71,7 +84,7 @@ namespace Opm
using WellInterface::MatrixBlockType;
using WellInterface::VectorBlockType; */
StandardWell(const Well* well, const size_t time_step, const Wells* wells);
StandardWell(const Well* well, const int time_step, const Wells* wells);
/// the densities of the fluid in each perforation
virtual const std::vector<double>& perfDensities() const;

View File

@ -49,6 +49,7 @@
#include <opm/autodiff/BlackoilModelParameters.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoilDense.hpp>
#include <opm/autodiff/RateConverter.hpp>
#include <opm/autodiff/WellInterface.hpp>
#include<dune/common/fmatrix.hh>
#include<dune/istl/bcrsmatrix.hh>
#include<dune/istl/matrixmatrix.hh>
@ -60,6 +61,7 @@
#include <math.h>
namespace Opm {
enum WellVariablePositions {
@ -296,6 +298,18 @@ enum WellVariablePositions {
const Wells* wells_;
const std::vector< const Well* > wells_ecl_;
// a vector of all the wells.
// eventually, the wells_ above should be gone.
// the name is just temporary
// later, might make share_ptr const later.
// TODO: forget why make it share_ptr instead of unique_ptr
std::vector<std::shared_ptr<WellInterface> > well_container_;
std::vector<std::shared_ptr<WellInterface> >
createWellContainer(const std::vector<const Well*>& wells_ecl,
const Wells* wells_arg,
const int time_step);
// Well collection is used to enforce the group control
WellCollection* well_collection_;

View File

@ -1,4 +1,5 @@
#include <opm/autodiff/StandardWell.hpp>
namespace Opm {
@ -15,6 +16,7 @@ namespace Opm {
: wells_active_(wells_arg!=nullptr)
, wells_(wells_arg)
, wells_ecl_(wells_ecl)
, well_container_(createWellContainer(wells_ecl, wells_arg, current_timeIdx) )
, well_collection_(well_collection)
, param_(param)
, terminal_output_(terminal_output)
@ -144,6 +146,63 @@ namespace Opm {
template<typename TypeTag>
std::vector<std::shared_ptr<WellInterface> >
StandardWellsDense<TypeTag>::
createWellContainer(const std::vector<const Well*>& wells_ecl,
const Wells* wells_arg,
const int time_step)
{
std::vector<std::shared_ptr<WellInterface> > wells_container;
// There might be no wells in the process
if (localWellsActive()) {
const int nw = wells_arg->number_of_wells;
wells_container.reserve(nw);
// With the following way, it will have the same order with wells struct
// Hopefully, it can generate the same residual history with master branch
for (int w = 0; w < nw; ++w) {
const std::string well_name = std::string(wells_arg->name[w]);
// finding the location of the well in wells_ecl
const int nw_wells_ecl = wells_ecl.size();
int index_well = 0;
for (; index_well < nw_wells_ecl; ++index_well) {
if (well_name == wells_ecl[index_well]->name()) {
break;
}
}
// It should be able to find in wells_ecl.
if (index_well == nw_wells_ecl) {
OPM_THROW(std::logic_error, "Could not find well " << well_name << " in wells_ecl ");
}
const Well* well_ecl = wells_ecl[index_well];
if (well_ecl->getStatus(time_step) == WellCommon::SHUT) {
continue;
}
if (well_ecl->isMultiSegment(time_step)) {
OPM_THROW(Opm::NumericalProblem, "Not handling Multisegment Wells for now");
}
// Basically, we are handling all the wells as StandardWell for the moment
// TODO: to be changed when we begin introducing MultisegmentWell
wells_container.push_back(std::make_shared<StandardWell>(well_ecl, time_step, wells_arg) );
}
}
return wells_container;
}
template<typename TypeTag>
SimulatorReport
StandardWellsDense<TypeTag>::

View File

@ -35,7 +35,7 @@ namespace Opm {
namespace wellhelpers
{
inline
double rateToCompare(const std::vector<double>& well_phase_flow_rate,
const int well,
@ -147,6 +147,15 @@ namespace Opm {
return dp;
}
inline
double computeHydrostaticCorrection(const double well_ref_depth, const double vfp_ref_depth,
const double rho, const double gravity) {
const double dh = vfp_ref_depth - well_ref_depth;
const double dp = rho * gravity * dh;
return dp;
}
template <class Vector>
inline
Vector computeHydrostaticCorrection(const Wells& wells, const Vector vfp_ref_depth,

View File

@ -29,7 +29,7 @@ namespace Opm
WellInterface::
WellInterface(const Well* well, const size_t time_step, const Wells* wells)
WellInterface(const Well* well, const int time_step, const Wells* wells)
{
// TODO: trying to use wells struct as little as possible here, be prepared to

View File

@ -54,7 +54,7 @@ namespace Opm
using WellState = WellStateFullyImplicitBlackoilDense;
/// Constructor
WellInterface(const Well* well, const size_t time_step, const Wells* wells);
WellInterface(const Well* well, const int time_step, const Wells* wells);
/// Well name.
const std::string& name() const;