All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
BlackoilMultiSegmentModel.hpp
1 /*
2  Copyright 2013, 2015 SINTEF ICT, Applied Mathematics.
3 
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #ifndef OPM_BLACKOILMULTISEGMENTMODEL_HEADER_INCLUDED
21 #define OPM_BLACKOILMULTISEGMENTMODEL_HEADER_INCLUDED
22 
23 #include <opm/core/simulator/BlackoilState.hpp>
24 #include <opm/autodiff/BlackoilModelBase.hpp>
25 #include <opm/autodiff/BlackoilModelParameters.hpp>
26 #include <opm/autodiff/WellStateMultiSegment.hpp>
27 #include <opm/autodiff/WellMultiSegment.hpp>
28 #include <opm/autodiff/StandardWells.hpp>
29 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
30 
31 #include <opm/autodiff/MultisegmentWells.hpp>
32 
33 namespace Opm {
34 
36  {
37  explicit BlackoilMultiSegmentSolutionState(const int np)
39  , segp ( ADB::null())
40  , segqs ( ADB::null())
41  {
42  }
43  ADB segp; // the segment pressures
44  ADB segqs; // the segment phase rate in surface volume
45  };
46 
54  template<class Grid>
55  class BlackoilMultiSegmentModel : public BlackoilModelBase<Grid, MultisegmentWells, BlackoilMultiSegmentModel<Grid> >
56  {
57  public:
58 
60  typedef typename Base::ReservoirState ReservoirState;
61  typedef typename Base::WellState WellState;
63 
64  friend Base;
65 
66  // --------- Public methods ---------
67 
84  BlackoilMultiSegmentModel(const typename Base::ModelParameters& param,
85  const Grid& grid ,
86  const BlackoilPropsAdFromDeck& fluid,
87  const DerivedGeology& geo ,
88  const RockCompressibility* rock_comp_props,
89  const MultisegmentWells& well_model,
90  const NewtonIterationBlackoilInterface& linsolver,
91  std::shared_ptr< const EclipseState > eclState,
92  const bool has_disgas,
93  const bool has_vapoil,
94  const bool terminal_output);
95 
100  void prepareStep(const SimulatorTimerInterface& timer,
101  const ReservoirState& reservoir_state,
102  const WellState& well_state);
103 
104 
109  SimulatorReport
110  assemble(const ReservoirState& reservoir_state,
111  WellState& well_state,
112  const bool initial_assembly);
113 
114  using Base::numPhases;
115  using Base::numMaterials;
116  using Base::materialName;
117 
118  protected:
119  // --------- Data members ---------
120 
121  // For non-segmented wells, it should be the density calculated with AVG or SEG way.
122  // while usually SEG way by default.
123  using Base::pvdt_;
124  using Base::geo_;
125  using Base::active_;
126  using Base::sd_;
127  using Base::fluid_;
129  using Base::grid_;
130  using Base::canph_;
131  using Base::residual_;
132  using Base::isSg_;
133  using Base::isRs_;
134  using Base::isRv_;
135  using Base::has_disgas_;
136  using Base::has_vapoil_;
137  using Base::cells_;
138  using Base::param_;
139  using Base::linsolver_;
140  using Base::phaseCondition_;
141  using Base::vfp_properties_;
142  using Base::well_model_;
143 
144  using Base::wellModel;
145  // using Base::wells;
146  using Base::wellsActive;
148  using Base::phaseCondition;
149  using Base::fluidRvSat;
150  using Base::fluidRsSat;
151  using Base::fluidDensity;
153  using Base::computeGasPressure;
154  using Base::dpMaxRel;
155  using Base::dsMax;
156  using Base::drMaxRel;
158  using Base::maxResidualAllowed;
159  using Base::variableState;
160  // using Base::variableWellStateIndices;
161  using Base::asImpl;
162  using Base::variableReservoirStateInitials;
163 
164 
165  const std::vector<WellMultiSegmentConstPtr>& wellsMultiSegment() const { return well_model_.msWells(); }
166 
167  const MultisegmentWells::MultisegmentWellOps& msWellOps() const { return well_model_.wellOps(); }
168 
169  SimulatorReport
170  solveWellEq(const std::vector<ADB>& mob_perfcells,
171  const std::vector<ADB>& b_perfcells,
172  const ReservoirState& reservoir_state,
173  SolutionState& state,
174  WellState& well_state);
175 
176  void
177  makeConstantState(SolutionState& state) const;
178 
179  // TODO: added since the interfaces of the function are different
180  // TODO: for StandardWells and MultisegmentWells
181  void
182  computeWellConnectionPressures(const SolutionState& state,
183  const WellState& well_state);
184 
185  };
186 
188  template <class GridT>
190  {
191  typedef BlackoilState ReservoirState;
195  };
196 
197 
198 
199 
200 } // namespace Opm
201 
202 #include "BlackoilMultiSegmentModel_impl.hpp"
203 
204 #endif // OPM_BLACKOILMULTISEGMENTMODEL_HEADER_INCLUDED
BlackoilMultiSegmentModel(const typename Base::ModelParameters &param, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const MultisegmentWells &well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclState, const bool has_disgas, const bool has_vapoil, const bool terminal_output)
Construct the model.
Definition: BlackoilMultiSegmentModel_impl.hpp:60
static AutoDiffBlock null()
Construct an empty AutoDiffBlock.
Definition: AutoDiffBlock.hpp:104
double convergenceReduction(const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic > &B, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic > &tempV, const Eigen::Array< double, Eigen::Dynamic, Eigen::Dynamic > &R, std::vector< double > &R_sum, std::vector< double > &maxCoeff, std::vector< double > &B_avg, std::vector< double > &maxNormWell, int nc) const
Compute the reduction within the convergence check.
Definition: BlackoilMultiSegmentModel.hpp:35
A model implementation for three-phase black oil.
Definition: BlackoilModelBase.hpp:76
void prepareStep(const SimulatorTimerInterface &timer, const ReservoirState &reservoir_state, const WellState &well_state)
Called once before each time step.
Definition: BlackoilMultiSegmentModel_impl.hpp:83
Definition: MultisegmentWells.hpp:62
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
SimulatorReport assemble(const ReservoirState &reservoir_state, WellState &well_state, const bool initial_assembly)
Assemble the residual and Jacobian of the nonlinear system.
Definition: BlackoilMultiSegmentModel_impl.hpp:132
void updatePrimalVariableFromState(const ReservoirState &state)
update the primal variable for Sg, Rv or Rs.
int numPhases() const
The number of active fluid phases in the model.
Traits to encapsulate the types used by classes using or extending this model.
Definition: BlackoilModelBase.hpp:60
A model implementation for three-phase black oil with support for multi-segment wells.
Definition: BlackoilMultiSegmentModel.hpp:55
MultisegmentWells & wellModel()
return the WellModel object
Definition: BlackoilModelBase.hpp:268
bool terminal_output_
Whether we print something to std::cout.
Definition: BlackoilModelBase.hpp:353
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
int numMaterials() const
The number of active materials in the model.
Struct for containing iteration variables.
Definition: DefaultBlackoilSolutionState.hpp:29
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
BlackoilMultiSegmentModel< Grid > & asImpl()
Access the most-derived class used for static polymorphism (CRTP).
Definition: BlackoilModelBase.hpp:370
const std::string & materialName(int material_index) const
The name of an active material in the model.
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
The state of a set of multi-sgemnet wells.
Definition: WellStateMultiSegment.hpp:45
void updatePhaseCondFromPrimalVariable(const ReservoirState &state)
Update the phaseCondition_ member based on the primalVariable_ member.
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilModelBase.hpp:386
Class for handling the multi-segment well model.
Definition: MultisegmentWells.hpp:55