21#ifndef OPM_FULLYIMPLICITCOMPRESSIBLEPOLYMERSOLVER_HEADER_INCLUDED
22#define OPM_FULLYIMPLICITCOMPRESSIBLEPOLYMERSOLVER_HEADER_INCLUDED
24#include <opm/autodiff/AutoDiffBlock.hpp>
25#include <opm/autodiff/AutoDiffHelpers.hpp>
26#include <opm/autodiff/BlackoilPropsAdInterface.hpp>
27#include <opm/autodiff/NewtonIterationBlackoilInterface.hpp>
28#include <opm/autodiff/LinearisedBlackoilResidual.hpp>
32#include <opm/core/utility/parameters/ParameterGroup.hpp>
34struct UnstructuredGrid;
40 class RockCompressibility;
41 class NewtonIterationBlackoilInterface;
42 class PolymerBlackoilState;
43 class WellStateFullyImplicitBlackoil;
67 const BlackoilPropsAdInterface& fluid,
68 const DerivedGeology& geo ,
69 const RockCompressibility* rock_comp_props,
72 const NewtonIterationBlackoilInterface& linsolver);
103 typedef AutoDiffBlock<double> ADB;
106 typedef Eigen::Array<double,
109 Eigen::RowMajor> DataBlock;
111 struct ReservoirResidualQuant {
112 ReservoirResidualQuant();
113 std::vector<ADB> accum;
118 std::vector<ADB> ads;
121 struct SolutionState {
122 SolutionState(
const int np);
125 std::vector<ADB> saturation;
132 WellOps(
const Wells& wells);
133 Eigen::SparseMatrix<double> w2p;
134 Eigen::SparseMatrix<double> p2w;
137 enum { Water = BlackoilPropsAdInterface::Water,
138 Oil = BlackoilPropsAdInterface::Oil };
141 const UnstructuredGrid& grid_;
142 const BlackoilPropsAdInterface& fluid_;
143 const DerivedGeology& geo_;
144 const RockCompressibility* rock_comp_props_;
145 const PolymerPropsAd& polymer_props_ad_;
147 const NewtonIterationBlackoilInterface& linsolver_;
148 const std::vector<int> cells_;
153 std::vector<PhasePresence> phaseCondition_;
154 std::vector<ReservoirResidualQuant> rq_;
158 LinearisedBlackoilResidual residual_;
160 unsigned int newtonIterations_;
161 unsigned int linearIterations_;
165 constantState(
const PolymerBlackoilState& x,
166 const WellStateFullyImplicitBlackoil& xw);
169 variableState(
const PolymerBlackoilState& x,
170 const WellStateFullyImplicitBlackoil& xw);
173 computeAccum(
const SolutionState& state,
177 assemble(
const double dt,
178 const PolymerBlackoilState& x,
179 const WellStateFullyImplicitBlackoil& xw,
180 const std::vector<double>& polymer_inflow);
182 V solveJacobianSystem()
const;
184 void updateState(
const V& dx,
185 PolymerBlackoilState& state,
186 WellStateFullyImplicitBlackoil& well_state)
const;
189 computeRelPerm(
const SolutionState& state)
const;
192 computePressures(
const SolutionState& state)
const;
195 computeMassFlux(
const int actph ,
197 const std::vector<ADB>& kr ,
198 const SolutionState& state );
201 computeMassFlux(
const V& trans,
205 const SolutionState& state);
208 computeFracFlow(
const ADB& kro,
213 computeCmax(PolymerBlackoilState& state);
216 computeMc(
const SolutionState& state)
const;
219 rockPorosity(
const ADB& p)
const;
222 rockPermeability(
const ADB& p)
const;
225 residualNorm()
const;
228 fluidViscosity(
const int phase,
231 const std::vector<PhasePresence>& cond,
232 const std::vector<int>& cells)
const;
235 fluidReciprocFVF(
const int phase,
238 const std::vector<PhasePresence>& cond,
239 const std::vector<int>& cells)
const;
242 fluidDensity(
const int phase,
245 const std::vector<PhasePresence>& cond,
246 const std::vector<int>& cells)
const;
249 poroMult(
const ADB& p)
const;
252 transMult(
const ADB& p)
const;
254 const std::vector<PhasePresence>
255 phaseCondition()
const {
return phaseCondition_; }
258 classifyCondition(
const PolymerBlackoilState& state);
Definition: FullyImplicitCompressiblePolymerSolver.hpp:54
int nonlinearIterations() const
FullyImplicitCompressiblePolymerSolver(const UnstructuredGrid &grid, const BlackoilPropsAdInterface &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const PolymerPropsAd &polymer_props_ad, const Wells &wells, const NewtonIterationBlackoilInterface &linsolver)
double relativeChange(const PolymerBlackoilState &previous, const PolymerBlackoilState ¤t) const
Evaluate the relative changes in the physical variables.
int step(const double dt, PolymerBlackoilState &state, WellStateFullyImplicitBlackoilPolymer &wstate)
parameter::ParameterGroup SolverParameters
Not used by this class except to satisfy interface requirements.
Definition: FullyImplicitCompressiblePolymerSolver.hpp:93
int linearIterations() const
const FullyImplicitCompressiblePolymerSolver & model() const
There is no separate model class for this solver, return itself.
Definition: PolymerBlackoilState.hpp:34
Definition: PolymerPropsAd.hpp:33
Definition: WellStateFullyImplicitBlackoilPolymer.hpp:29
Definition: CompressibleTpfaPolymer.hpp:33