New WellStateFullyImplicitBlackoilDense

The WellStateFullyImplicitBlackoil is reverted to master and all
extentions moved to WellStateFullyImplicitBlackoilDense
This commit is contained in:
Tor Harald Sandve 2016-11-01 10:07:28 +01:00
parent ae0ebf1200
commit 362968c315
8 changed files with 195 additions and 116 deletions

View File

@ -202,6 +202,7 @@ list (APPEND PUBLIC_HEADER_FILES
opm/autodiff/TransportSolverTwophaseAd.hpp
opm/autodiff/WellDensitySegmented.hpp
opm/autodiff/WellStateFullyImplicitBlackoil.hpp
opm/autodiff/WellStateFullyImplicitBlackoilDense.hpp
opm/autodiff/WellStateFullyImplicitBlackoilSolvent.hpp
opm/autodiff/SimulatorFullyImplicitBlackoilOutput.hpp
opm/autodiff/VFPProperties.hpp

View File

@ -105,7 +105,7 @@ namespace Opm {
public:
// --------- Types and enums ---------
typedef BlackoilState ReservoirState;
typedef WellStateFullyImplicitBlackoil WellState;
typedef WellStateFullyImplicitBlackoilDense WellState;
typedef BlackoilModelParameters ModelParameters;
typedef DefaultBlackoilSolutionState SolutionState;

View File

@ -27,7 +27,7 @@
#include <opm/autodiff/NonlinearSolver.hpp>
#include <opm/autodiff/BlackoilModelEbos.hpp>
#include <opm/autodiff/BlackoilModelParameters.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoilDense.hpp>
#include <opm/autodiff/StandardWellsDense.hpp>
#include <opm/autodiff/RateConverter.hpp>
#include <opm/autodiff/SimFIBODetails.hpp>
@ -54,7 +54,7 @@ public:
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, Indices) BlackoilIndices;
typedef WellStateFullyImplicitBlackoil WellState;
typedef WellStateFullyImplicitBlackoilDense WellState;
typedef BlackoilState ReservoirState;
typedef BlackoilOutputWriterEbos OutputWriter;
typedef BlackoilModelEbos Model;

View File

@ -204,7 +204,7 @@ namespace Opm
BlackoilOutputWriterEbos::
restore(SimulatorTimerInterface& timer,
BlackoilState& state,
WellStateFullyImplicitBlackoil& wellState,
WellStateFullyImplicitBlackoilDense& wellState,
const std::string& filename,
const int desiredResportStep )
{

View File

@ -36,7 +36,7 @@
#include <opm/autodiff/GridHelpers.hpp>
#include <opm/autodiff/ParallelDebugOutput.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoilDense.hpp>
#include <opm/autodiff/ThreadHandle.hpp>
#include <opm/autodiff/AutoDiffBlock.hpp>
@ -128,7 +128,7 @@ namespace Opm
void restore(SimulatorTimerInterface& timer,
BlackoilState& state,
WellStateFullyImplicitBlackoil& wellState,
WellStateFullyImplicitBlackoilDense& wellState,
const std::string& filename,
const int desiredReportStep);
@ -138,7 +138,7 @@ namespace Opm
const double* permeability,
const Grid& grid,
SimulationDataContainer& simulatorstate,
WellStateFullyImplicitBlackoil& wellstate);
WellStateFullyImplicitBlackoilDense& wellstate);
bool isRestart() const;
@ -230,7 +230,7 @@ namespace Opm
const double* permeability,
const Grid& grid,
SimulationDataContainer& simulatorstate,
WellStateFullyImplicitBlackoil& wellstate)
WellStateFullyImplicitBlackoilDense& wellstate)
{
// gives a dummy dynamic_list_econ_limited
DynamicListEconLimited dummy_list_econ_limited;

View File

@ -49,7 +49,7 @@
#include <opm/autodiff/BlackoilDetails.hpp>
#include <opm/autodiff/BlackoilModelParameters.hpp>
#include <opm/autodiff/LinearisedBlackoilResidual.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoilDense.hpp>
#include<dune/common/fmatrix.hh>
#include<dune/istl/bcrsmatrix.hh>
#include<dune/istl/matrixmatrix.hh>
@ -69,7 +69,7 @@ namespace Opm {
using ADB = AutoDiffBlock<double>;
typedef DenseAd::Evaluation<double, /*size=*/6> EvalWell;
typedef WellStateFullyImplicitBlackoil WellState;
typedef WellStateFullyImplicitBlackoilDense WellState;
typedef BlackoilModelParameters ModelParameters;
//typedef AutoDiffBlock<double> ADB;

View File

@ -24,7 +24,6 @@
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/core/simulator/WellState.hpp>
#include <opm/autodiff/BlackoilModelEnums.hpp>
#include <opm/core/props/BlackoilPhases.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <vector>
@ -107,68 +106,12 @@ namespace Opm
well_potentials_.clear();
well_potentials_.resize(nperf * np, 0.0);
well_solutions_.clear();
well_solutions_.resize(nw * np, 0.0);
std::vector<double> g = {1.0,1.0,0.01};
for (int w = 0; w < nw; ++w) {
const WellControls* wc = wells->ctrls[w];
// The current control in the well state overrides
// the current control set in the Wells struct, which
// is instead treated as a default.
const int current = current_controls_[w];
const WellType& well_type = wells->type[w];
switch (well_controls_iget_type(wc, current)) {
case THP: // Intentional fall-through
case BHP:
{
if (well_type == INJECTOR) {
for (int p = 0; p < np; ++p) {
well_solutions_[w] += wellRates()[np*w + p] * wells->comp_frac[np*w + p];
}
} else {
for (int p = 0; p < np; ++p) {
well_solutions_[w] += g[p] * wellRates()[np*w + p];
}
}
}
break;
case RESERVOIR_RATE: // Intentional fall-through
case SURFACE_RATE:
{
wellSolutions()[w] = bhp()[w];
}
break;
}
assert(np == 3);
double total_rates = 0.0;
for (int p = 0; p < np; ++p) {
total_rates += g[p] * wellRates()[np*w + p];
}
if(std::abs(total_rates) > 0) {
wellSolutions()[nw + w] = g[Water] * wellRates()[np*w + Water] / total_rates; //wells->comp_frac[np*w + Water]; // Water;
wellSolutions()[2*nw + w] = g[Gas] * wellRates()[np*w + Gas] / total_rates ; //wells->comp_frac[np*w + Gas]; //Gas
} else {
wellSolutions()[nw + w] = wells->comp_frac[np*w + Water];
wellSolutions()[2*nw + w] = wells->comp_frac[np*w + Gas];
}
}
// intialize wells that have been there before
// order may change so the mapping is based on the well name
if( ! prevState.wellMap().empty() )
{
typedef typename WellMapType :: const_iterator const_iterator;
const_iterator end = prevState.wellMap().end();
int nw_old = prevState.bhp().size();
for (int w = 0; w < nw; ++w) {
std::string name( wells->name[ w ] );
const_iterator it = prevState.wellMap().find( name );
@ -181,31 +124,11 @@ namespace Opm
bhp()[ newIndex ] = prevState.bhp()[ oldIndex ];
// wellrates
double total_well_rates = 0.0;
for( int i=0, idx=newIndex*np, oldidx=oldIndex*np; i<np; ++i, ++idx, ++oldidx )
{
total_well_rates += prevState.wellRates()[ oldidx ];
}
//if (std::abs(total_well_rates) > 0) {
for( int i=0, idx=newIndex*np, oldidx=oldIndex*np; i<np; ++i, ++idx, ++oldidx )
{
wellRates()[ idx ] = prevState.wellRates()[ oldidx ];
}
// wellSolutions
//if (std::abs(total_well_rates) > 0.0) {
//wellSolutions()[ 0*nw + newIndex ] = prevState.wellSolutions()[0 * nw_old + oldIndex ];
//if (wells->type[w] == PRODUCER) {
for( int i = 0; i < np; ++i)
{
wellSolutions()[ i*nw + newIndex ] = prevState.wellSolutions()[i * nw_old + oldIndex ];
}
//}
//}
// perfPhaseRates
int oldPerf_idx = (*it).second[ 1 ];
const int num_perf_old_well = (*it).second[ 2 ];
@ -243,9 +166,6 @@ namespace Opm
// If the set of controls have changed, this may not be identical
// to the last control, but it must be a valid control.
currentControls()[ newIndex ] = old_control_index;
WellControls* wc = wells->ctrls[newIndex];
well_controls_set_current( wc, old_control_index);
}
}
@ -272,11 +192,7 @@ namespace Opm
std::vector<double>& wellPotentials() { return well_potentials_; }
const std::vector<double>& wellPotentials() const { return well_potentials_; }
/// One rate per phase and well connection.
std::vector<double>& wellSolutions() { return well_solutions_; }
const std::vector<double>& wellSolutions() const { return well_solutions_; }
data::Wells report(const PhaseUsage& pu) const override {
data::Wells report(const PhaseUsage &pu) const override {
data::Wells res = WellState::report(pu);
const int nw = this->numWells();
@ -327,31 +243,10 @@ namespace Opm
return res;
}
WellStateFullyImplicitBlackoil() = default;
WellStateFullyImplicitBlackoil( const WellStateFullyImplicitBlackoil& rhs ) :
BaseType(rhs),
perfphaserates_( rhs.perfphaserates_ ),
current_controls_( rhs.current_controls_ ),
well_potentials_( rhs.well_potentials_ ),
well_solutions_( rhs.well_solutions_ )
{}
WellStateFullyImplicitBlackoil& operator=( const WellStateFullyImplicitBlackoil& rhs ) {
BaseType::operator =(rhs);
this->perfPhaseRates() = rhs.perfPhaseRates();
this->currentControls() = rhs.currentControls();
this->wellPotentials() = rhs.wellPotentials();
this->wellSolutions() = rhs.wellSolutions();
return *this;
}
private:
std::vector<double> perfphaserates_;
std::vector<int> current_controls_;
std::vector<double> well_potentials_;
std::vector<double> well_solutions_;
};
} // namespace Opm

View File

@ -0,0 +1,183 @@
/*
Copyright 2016 IRIS
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_WELLSTATEFULLYIMPLICITBLACKOILDENSE_HEADER_INCLUDED
#define OPM_WELLSTATEFULLYIMPLICITBLACKOILDENSE_HEADER_INCLUDED
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/core/simulator/WellState.hpp>
#include <opm/autodiff/BlackoilModelEnums.hpp>
#include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
#include <opm/core/props/BlackoilPhases.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <vector>
#include <cassert>
#include <string>
#include <utility>
#include <map>
#include <algorithm>
#include <array>
namespace Opm
{
/// The state of a set of wells, tailored for use by the fully
/// implicit blackoil simulator.
class WellStateFullyImplicitBlackoilDense
: public WellStateFullyImplicitBlackoil
{
typedef WellStateFullyImplicitBlackoil BaseType;
public:
typedef BaseType :: WellMapType WellMapType;
using BaseType :: wellRates;
using BaseType :: bhp;
using BaseType :: perfPress;
using BaseType :: wellMap;
using BaseType :: numWells;
using BaseType :: numPhases;
using BaseType :: perfPhaseRates;
using BaseType :: currentControls;
using BaseType :: wellPotentials;
/// Allocate and initialize if wells is non-null. Also tries
/// to give useful initial values to the bhp(), wellRates()
/// and perfPhaseRates() fields, depending on controls
template <class State, class PrevState>
void init(const Wells* wells, const State& state, const PrevState& prevState)
{
// call init on base class
BaseType :: init(wells, state, prevState);
// if there are no well, do nothing in init
if (wells == 0) {
return;
}
const int nw = wells->number_of_wells;
if( nw == 0 ) return ;
const int np = wells->number_of_phases;
well_solutions_.clear();
well_solutions_.resize(nw * np, 0.0);
std::vector<double> g = {1.0,1.0,0.01};
for (int w = 0; w < nw; ++w) {
WellControls* wc = wells->ctrls[w];
// The current control in the well state overrides
// the current control set in the Wells struct, which
// is instead treated as a default.
const int current = currentControls()[w];
well_controls_set_current( wc, current);
const WellType& well_type = wells->type[w];
switch (well_controls_iget_type(wc, current)) {
case THP: // Intentional fall-through
case BHP:
{
if (well_type == INJECTOR) {
for (int p = 0; p < np; ++p) {
well_solutions_[w] += wellRates()[np*w + p] * wells->comp_frac[np*w + p];
}
} else {
for (int p = 0; p < np; ++p) {
well_solutions_[w] += g[p] * wellRates()[np*w + p];
}
}
}
break;
case RESERVOIR_RATE: // Intentional fall-through
case SURFACE_RATE:
{
wellSolutions()[w] = bhp()[w];
}
break;
}
assert(np == 3);
double total_rates = 0.0;
for (int p = 0; p < np; ++p) {
total_rates += g[p] * wellRates()[np*w + p];
}
if(std::abs(total_rates) > 0) {
wellSolutions()[nw + w] = g[Water] * wellRates()[np*w + Water] / total_rates; //wells->comp_frac[np*w + Water]; // Water;
wellSolutions()[2*nw + w] = g[Gas] * wellRates()[np*w + Gas] / total_rates ; //wells->comp_frac[np*w + Gas]; //Gas
} else {
wellSolutions()[nw + w] = wells->comp_frac[np*w + Water];
wellSolutions()[2*nw + w] = wells->comp_frac[np*w + Gas];
}
}
// intialize wells that have been there before
// order may change so the mapping is based on the well name
if( ! prevState.wellMap().empty() )
{
typedef typename WellMapType :: const_iterator const_iterator;
const_iterator end = prevState.wellMap().end();
int nw_old = prevState.bhp().size();
for (int w = 0; w < nw; ++w) {
std::string name( wells->name[ w ] );
const_iterator it = prevState.wellMap().find( name );
if( it != end )
{
const int oldIndex = (*it).second[ 0 ];
const int newIndex = w;
// wellSolutions
for( int i = 0; i < np; ++i)
{
wellSolutions()[ i*nw + newIndex ] = prevState.wellSolutions()[i * nw_old + oldIndex ];
}
}
}
}
}
template <class State>
void resize(const Wells* wells, const State& state) {
const WellStateFullyImplicitBlackoilDense dummy_state{}; // Init with an empty previous state only resizes
init(wells, state, dummy_state) ;
}
/// One rate per phase and well connection.
std::vector<double>& wellSolutions() { return well_solutions_; }
const std::vector<double>& wellSolutions() const { return well_solutions_; }
data::Wells report(const PhaseUsage& pu) const override {
data::Wells res = BaseType::report(pu);
return res;
}
private:
std::vector<double> well_solutions_;
};
} // namespace Opm
#endif // OPM_WELLSTATEFULLYIMPLICITBLACKOILDENSE_HEADER_INCLUDED