2017-08-17 08:49:54 -05:00
/*
Copyright 2017 SINTEF Digital , Mathematics and Cybernetics .
Copyright 2017 Statoil ASA .
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_MULTISEGMENTWELL_HEADER_INCLUDED
# define OPM_MULTISEGMENTWELL_HEADER_INCLUDED
# include <opm/autodiff/WellInterface.hpp>
namespace Opm
{
template < typename TypeTag >
class MultisegmentWell : public WellInterface < TypeTag >
{
public :
typedef WellInterface < TypeTag > Base ;
using typename Base : : WellState ;
2017-09-01 10:37:51 -05:00
using typename Base : : Simulator ;
2017-09-08 04:09:54 -05:00
using typename Base : : IntensiveQuantities ;
using typename Base : : FluidSystem ;
2017-09-01 10:37:51 -05:00
using typename Base : : ModelParameters ;
2017-09-08 07:58:57 -05:00
using typename Base : : MaterialLaw ;
2017-12-04 02:14:08 -06:00
using typename Base : : Indices ;
2017-11-23 01:37:30 -06:00
using typename Base : : RateConverterType ;
2017-09-26 03:07:06 -05:00
/// the number of reservior equations
using Base : : numEq ;
2017-08-17 08:49:54 -05:00
2017-10-13 03:16:44 -05:00
using Base : : has_solvent ;
using Base : : has_polymer ;
2017-11-08 06:57:36 -06:00
using Base : : Water ;
using Base : : Oil ;
using Base : : Gas ;
2017-10-13 03:16:44 -05:00
2017-08-17 08:49:54 -05:00
// TODO: for now, not considering the polymer, solvent and so on to simplify the development process.
2017-10-16 09:46:28 -05:00
2017-09-12 06:16:57 -05:00
// TODO: we need to have order for the primary variables and also the order for the well equations.
// sometimes, they are similar, while sometimes, they can have very different forms.
2017-09-26 03:07:06 -05:00
// TODO: the following system looks not rather flexible. Looking into all kinds of possibilities
// TODO: gas is always there? how about oil water case?
// Is it gas oil two phase case?
2017-12-04 02:14:08 -06:00
static const bool gasoil = numEq = = 2 & & ( Indices : : compositionSwitchIdx > = 0 ) ;
2017-09-26 03:07:06 -05:00
static const int GTotal = 0 ;
static const int WFrac = gasoil ? - 1000 : 1 ;
static const int GFrac = gasoil ? 1 : 2 ;
static const int SPres = gasoil ? 2 : 3 ;
2017-08-17 08:49:54 -05:00
/// the number of well equations // TODO: it should have a more general strategy for it
2017-09-26 03:07:06 -05:00
static const int numWellEq = GET_PROP_VALUE ( TypeTag , EnablePolymer ) ? numEq : numEq + 1 ;
2017-08-17 08:49:54 -05:00
using typename Base : : Scalar ;
/// the matrix and vector types for the reservoir
using typename Base : : Mat ;
using typename Base : : BVector ;
using typename Base : : Eval ;
// sparsity pattern for the matrices
// [A C^T [x = [ res
// B D ] x_well] res_well]
// the vector type for the res_well and x_well
typedef Dune : : FieldVector < Scalar , numWellEq > VectorBlockWellType ;
typedef Dune : : BlockVector < VectorBlockWellType > BVectorWell ;
// the matrix type for the diagonal matrix D
typedef Dune : : FieldMatrix < Scalar , numWellEq , numWellEq > DiagMatrixBlockWellType ;
typedef Dune : : BCRSMatrix < DiagMatrixBlockWellType > DiagMatWell ;
// the matrix type for the non-diagonal matrix B and C^T
typedef Dune : : FieldMatrix < Scalar , numWellEq , numEq > OffDiagMatrixBlockWellType ;
typedef Dune : : BCRSMatrix < OffDiagMatrixBlockWellType > OffDiagMatWell ;
// TODO: for more efficient implementation, we should have EvalReservoir, EvalWell, and EvalRerservoirAndWell
// EvalR (Eval), EvalW, EvalRW
// TODO: for now, we only use one type to save some implementation efforts, while improve later.
typedef DenseAd : : Evaluation < double , /*size=*/ numEq + numWellEq > EvalWell ;
2017-11-23 01:37:30 -06:00
MultisegmentWell ( const Well * well , const int time_step , const Wells * wells ,
const ModelParameters & param ,
const RateConverterType & rate_converter ,
2017-11-30 09:31:48 -06:00
const int pvtRegionIdx ,
const int num_components ) ;
2017-08-17 08:49:54 -05:00
virtual void init ( const PhaseUsage * phase_usage_arg ,
const std : : vector < double > & depth_arg ,
const double gravity_arg ,
const int num_cells ) ;
virtual void initPrimaryVariablesEvaluation ( ) const ;
virtual void assembleWellEq ( Simulator & ebosSimulator ,
const double dt ,
WellState & well_state ,
bool only_wells ) ;
/// updating the well state based the control mode specified with current
// TODO: later will check wheter we need current
2017-11-30 10:14:29 -06:00
virtual void updateWellStateWithTarget ( WellState & well_state ) const ;
2017-08-17 08:49:54 -05:00
/// check whether the well equations get converged for this well
2018-10-25 04:57:47 -05:00
virtual ConvergenceReport getWellConvergence ( const std : : vector < double > & B_avg ) const ;
2017-08-17 08:49:54 -05:00
/// Ax = Ax - C D^-1 B x
virtual void apply ( const BVector & x , BVector & Ax ) const ;
/// r = r - C D^-1 Rw
virtual void apply ( BVector & r ) const ;
/// using the solution x to recover the solution xw for wells and applying
/// xw to update Well State
2017-10-16 09:46:28 -05:00
virtual void recoverWellSolutionAndUpdateWellState ( const BVector & x ,
2017-08-17 08:49:54 -05:00
WellState & well_state ) const ;
/// computing the well potentials for group control
virtual void computeWellPotentials ( const Simulator & ebosSimulator ,
const WellState & well_state ,
2017-09-05 04:09:58 -05:00
std : : vector < double > & well_potentials ) ;
2017-08-17 08:49:54 -05:00
virtual void updatePrimaryVariables ( const WellState & well_state ) const ;
2017-10-16 09:46:28 -05:00
virtual void solveEqAndUpdateWellState ( WellState & well_state ) ; // const?
2017-08-29 10:26:36 -05:00
2017-09-05 04:20:43 -05:00
virtual void calculateExplicitQuantities ( const Simulator & ebosSimulator ,
const WellState & well_state ) ; // should be const?
2017-09-05 04:09:58 -05:00
2017-08-28 10:30:22 -05:00
/// number of segments for this well
/// int number_of_segments_;
int numberOfSegments ( ) const ;
int numberOfPerforations ( ) const ;
2017-08-17 08:49:54 -05:00
protected :
int number_segments_ ;
// components of the pressure drop to be included
WellSegment : : CompPressureDropEnum compPressureDrop ( ) const ;
// multi-phase flow model
WellSegment : : MultiPhaseModelEnum multiphaseModel ( ) const ;
2018-06-17 07:42:14 -05:00
// get the WellSegments from the well_ecl_
const WellSegments & segmentSet ( ) const ;
2017-08-17 08:49:54 -05:00
2017-09-04 05:42:23 -05:00
// protected member variables from the Base class
2017-08-17 08:49:54 -05:00
using Base : : well_ecl_ ;
using Base : : number_of_perforations_ ; // TODO: can use well_ecl_?
2017-09-01 10:37:51 -05:00
using Base : : current_step_ ;
2017-09-04 05:42:23 -05:00
using Base : : index_of_well_ ;
using Base : : number_of_phases_ ;
2017-08-17 08:49:54 -05:00
2017-10-03 09:25:51 -05:00
// TODO: the current implementation really relies on the order of the
// perforation does not change from the parser to Wells structure.
using Base : : well_cells_ ;
2017-10-16 09:46:28 -05:00
using Base : : param_ ;
2017-08-17 08:49:54 -05:00
using Base : : well_index_ ;
2017-09-04 05:42:23 -05:00
using Base : : well_type_ ;
2017-09-04 07:47:27 -05:00
using Base : : first_perf_ ;
2017-09-08 07:58:57 -05:00
using Base : : saturation_table_number_ ;
using Base : : well_efficiency_factor_ ;
2017-09-12 07:26:32 -05:00
using Base : : gravity_ ;
2017-08-17 08:49:54 -05:00
using Base : : well_controls_ ;
2017-09-28 09:01:11 -05:00
using Base : : perf_depth_ ;
2017-11-30 09:31:48 -06:00
using Base : : num_components_ ;
2017-08-17 08:49:54 -05:00
2017-09-04 05:42:23 -05:00
// protected functions from the Base class
2017-09-04 07:47:27 -05:00
using Base : : phaseUsage ;
2017-09-04 10:58:37 -05:00
using Base : : name ;
2017-09-08 07:58:57 -05:00
using Base : : flowPhaseToEbosCompIdx ;
2017-12-04 02:14:08 -06:00
using Base : : ebosCompIdxToFlowCompIdx ;
2017-09-08 07:58:57 -05:00
using Base : : getAllowCrossFlow ;
2017-11-23 01:37:30 -06:00
using Base : : scalingFactor ;
2017-09-04 05:42:23 -05:00
2017-08-17 08:49:54 -05:00
// TODO: trying to use the information from the Well opm-parser as much
// as possible, it will possibly be re-implemented later for efficiency reason.
2017-08-31 11:22:48 -05:00
// the completions that is related to each segment
2017-10-16 07:39:07 -05:00
// the completions's ids are their index in the vector well_index_, well_cell_
2017-08-31 11:22:48 -05:00
// This is also assuming the order of the completions in Well is the same with
// the order of the completions in wells.
// it is for convinience reason. we can just calcuate the inforation for segment once then using it for all the perofrations
// belonging to this segment
2017-09-01 10:37:51 -05:00
std : : vector < std : : vector < int > > segment_perforations_ ;
// the inlet segments for each segment. It is for convinience and efficiency reason
std : : vector < std : : vector < int > > segment_inlets_ ;
2017-08-31 11:22:48 -05:00
2017-08-17 08:49:54 -05:00
// segment number is an ID of the segment, it is specified in the deck
// get the loation of the segment with a segment number in the segmentSet
2017-10-16 07:39:07 -05:00
int segmentNumberToIndex ( const int segment_number ) const ;
2017-08-17 08:49:54 -05:00
// TODO, the following should go to a class for computing purpose
// two off-diagonal matrices
mutable OffDiagMatWell duneB_ ;
mutable OffDiagMatWell duneC_ ;
// diagonal matrix for the well
2017-09-13 10:48:08 -05:00
mutable DiagMatWell duneD_ ;
2017-08-17 08:49:54 -05:00
// residuals of the well equations
2017-09-01 10:37:51 -05:00
mutable BVectorWell resWell_ ;
2017-08-17 08:49:54 -05:00
// the values for the primary varibles
// based on different solutioin strategies, the wells can have different primary variables
2017-08-29 10:26:36 -05:00
mutable std : : vector < std : : array < double , numWellEq > > primary_variables_ ;
2017-08-17 08:49:54 -05:00
// the Evaluation for the well primary variables, which contain derivativles and are used in AD calculation
2017-08-29 10:26:36 -05:00
mutable std : : vector < std : : array < EvalWell , numWellEq > > primary_variables_evaluation_ ;
2017-09-28 09:01:11 -05:00
// depth difference between perforations and the perforated grid cells
2017-10-02 11:29:03 -05:00
std : : vector < double > cell_perforation_depth_diffs_ ;
2017-09-05 08:16:04 -05:00
// pressure correction due to the different depth of the perforation and
// center depth of the grid block
2017-10-02 11:29:03 -05:00
std : : vector < double > cell_perforation_pressure_diffs_ ;
2017-09-05 08:16:04 -05:00
2017-09-08 04:09:54 -05:00
// depth difference between the segment and the peforation
// or in another way, the depth difference between the perforation and
// the segment the perforation belongs to
2017-10-02 11:29:03 -05:00
std : : vector < double > perforation_segment_depth_diffs_ ;
2017-09-08 04:09:54 -05:00
2017-09-05 08:16:04 -05:00
// the intial component compistion of segments
2017-09-05 08:34:31 -05:00
std : : vector < std : : vector < double > > segment_comp_initial_ ;
2017-09-05 08:16:04 -05:00
2017-09-08 04:09:54 -05:00
// the densities of segment fluids
// we should not have this member variable
std : : vector < EvalWell > segment_densities_ ;
2017-09-28 04:00:10 -05:00
// the viscosity of the segments
2017-09-28 03:04:04 -05:00
std : : vector < EvalWell > segment_viscosities_ ;
2017-09-28 06:11:35 -05:00
// the mass rate of the segments
std : : vector < EvalWell > segment_mass_rates_ ;
2017-09-12 07:26:32 -05:00
std : : vector < double > segment_depth_diffs_ ;
2017-09-01 10:37:51 -05:00
void initMatrixAndVectors ( const int num_cells ) const ;
2017-08-29 10:26:36 -05:00
// protected functions
// EvalWell getBhp(); this one should be something similar to getSegmentPressure();
// EvalWell getQs(); this one should be something similar to getSegmentRates()
// EValWell wellVolumeFractionScaled, wellVolumeFraction, wellSurfaceVolumeFraction ... these should have different names, and probably will be needed.
// bool crossFlowAllowed(const Simulator& ebosSimulator) const; probably will be needed
// xw = inv(D)*(rw - C*x)
void recoverSolutionWell ( const BVector & x , BVectorWell & xw ) const ;
// updating the well_state based on well solution dwells
void updateWellState ( const BVectorWell & dwells ,
2017-10-12 09:46:25 -05:00
const bool inner_iteration ,
2017-08-29 10:26:36 -05:00
WellState & well_state ) const ;
2017-10-12 05:01:48 -05:00
// initialize the segment rates with well rates
// when there is no more accurate way to initialize the segment rates, we initialize
// the segment rates based on well rates with a simple strategy
void initSegmentRatesWithWellRates ( WellState & well_state ) const ;
2017-09-05 04:20:43 -05:00
// computing the accumulation term for later use in well mass equations
2017-09-05 08:16:04 -05:00
void computeInitialComposition ( ) ;
// compute the pressure difference between the perforation and cell center
void computePerfCellPressDiffs ( const Simulator & ebosSimulator ) ;
// fraction value of the primary variables
// should we just use member variables to store them instead of calculating them again and again
2017-12-04 02:14:08 -06:00
EvalWell volumeFraction ( const int seg , const unsigned comp_idx ) const ;
2017-09-05 08:16:04 -05:00
// F_p / g_p, the basic usage of this value is because Q_p = G_t * F_p / G_p
EvalWell volumeFractionScaled ( const int seg , const int comp_idx ) const ;
// basically Q_p / \sigma_p Q_p
EvalWell surfaceVolumeFraction ( const int seg , const int comp_idx ) const ;
2017-09-05 04:20:43 -05:00
2017-09-08 04:09:54 -05:00
void computePerfRate ( const IntensiveQuantities & int_quants ,
const std : : vector < EvalWell > & mob_perfcells ,
const int seg ,
2017-10-02 11:29:03 -05:00
const int perf ,
2017-09-08 04:09:54 -05:00
const EvalWell & segment_pressure ,
const bool & allow_cf ,
std : : vector < EvalWell > & cq_s ) const ;
// convert a Eval from reservoir to contain the derivative related to wells
EvalWell extendEval ( const Eval & in ) const ;
2017-10-03 09:25:51 -05:00
// compute the fluid properties, such as densities, viscosities, and so on, in the segments
2017-09-08 04:09:54 -05:00
// They will be treated implicitly, so they need to be of Evaluation type
2017-09-21 08:58:01 -05:00
void computeSegmentFluidProperties ( const Simulator & ebosSimulator ) ;
2017-09-08 04:09:54 -05:00
EvalWell getSegmentPressure ( const int seg ) const ;
2017-09-08 07:58:57 -05:00
EvalWell getSegmentRate ( const int seg , const int comp_idx ) const ;
2017-09-11 10:31:42 -05:00
EvalWell getSegmentGTotal ( const int seg ) const ;
2017-09-08 07:58:57 -05:00
// get the mobility for specific perforation
void getMobility ( const Simulator & ebosSimulator ,
const int perf ,
std : : vector < EvalWell > & mob ) const ;
2017-09-11 10:31:42 -05:00
2017-09-15 06:42:15 -05:00
void assembleControlEq ( ) const ;
2017-09-12 06:16:57 -05:00
2017-09-15 06:42:15 -05:00
void assemblePressureEq ( const int seg ) const ;
2017-09-12 07:26:32 -05:00
// hytrostatic pressure loss
2017-09-18 07:10:09 -05:00
EvalWell getHydroPressureLoss ( const int seg ) const ;
2017-09-12 10:13:30 -05:00
2017-09-28 04:00:10 -05:00
// frictinal pressure loss
EvalWell getFrictionPressureLoss ( const int seg ) const ;
2017-09-28 08:12:09 -05:00
void handleAccelerationPressureLoss ( const int seg ) const ;
2017-09-12 10:13:30 -05:00
// handling the overshooting and undershooting of the fractions
2017-09-12 10:32:27 -05:00
void processFractions ( const int seg ) const ;
2017-09-15 09:28:10 -05:00
void updateWellStateFromPrimaryVariables ( WellState & well_state ) const ;
2017-09-26 03:07:06 -05:00
2017-09-28 08:12:09 -05:00
bool frictionalPressureLossConsidered ( ) const ;
bool accelerationalPressureLossConsidered ( ) const ;
2017-10-10 07:30:23 -05:00
// TODO: try to make ebosSimulator const, as it should be
void iterateWellEquations ( Simulator & ebosSimulator ,
const double dt ,
WellState & well_state ) ;
void assembleWellEqWithoutIteration ( Simulator & ebosSimulator ,
const double dt ,
WellState & well_state ,
bool only_wells ) ;
2017-09-27 11:16:45 -05:00
} ;
2017-09-13 10:10:34 -05:00
2017-08-17 08:49:54 -05:00
}
# include "MultisegmentWell_impl.hpp"
# endif // OPM_MULTISEGMENTWELL_HEADER_INCLUDED