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>
34 struct 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);
96 typedef AutoDiffBlock<double> ADB;
99 typedef Eigen::Array<double,
102 Eigen::RowMajor> DataBlock;
104 struct ReservoirResidualQuant {
105 ReservoirResidualQuant();
106 std::vector<ADB> accum;
111 std::vector<ADB> ads;
114 struct SolutionState {
115 SolutionState(
const int np);
118 std::vector<ADB> saturation;
125 WellOps(
const Wells& wells);
126 Eigen::SparseMatrix<double> w2p;
127 Eigen::SparseMatrix<double> p2w;
130 enum { Water = BlackoilPropsAdInterface::Water,
131 Oil = BlackoilPropsAdInterface::Oil };
134 const UnstructuredGrid& grid_;
135 const BlackoilPropsAdInterface& fluid_;
136 const DerivedGeology& geo_;
137 const RockCompressibility* rock_comp_props_;
138 const PolymerPropsAd& polymer_props_ad_;
140 const NewtonIterationBlackoilInterface& linsolver_;
141 const std::vector<int> cells_;
146 std::vector<PhasePresence> phaseCondition_;
147 std::vector<ReservoirResidualQuant> rq_;
151 LinearisedBlackoilResidual residual_;
153 unsigned int newtonIterations_;
154 unsigned int linearIterations_;
158 constantState(
const PolymerBlackoilState& x,
159 const WellStateFullyImplicitBlackoil& xw);
162 variableState(
const PolymerBlackoilState& x,
163 const WellStateFullyImplicitBlackoil& xw);
166 computeAccum(
const SolutionState& state,
170 assemble(
const double dt,
171 const PolymerBlackoilState& x,
172 const WellStateFullyImplicitBlackoil& xw,
173 const std::vector<double>& polymer_inflow);
175 V solveJacobianSystem()
const;
177 void updateState(
const V& dx,
178 PolymerBlackoilState& state,
179 WellStateFullyImplicitBlackoil& well_state)
const;
182 computeRelPerm(
const SolutionState& state)
const;
185 computePressures(
const SolutionState& state)
const;
188 computeMassFlux(
const int actph ,
190 const std::vector<ADB>& kr ,
191 const SolutionState& state );
194 computeMassFlux(
const V& trans,
198 const SolutionState& state);
201 computeFracFlow(
const ADB& kro,
206 computeCmax(PolymerBlackoilState& state);
209 computeMc(
const SolutionState& state)
const;
212 rockPorosity(
const ADB& p)
const;
215 rockPermeability(
const ADB& p)
const;
218 residualNorm()
const;
221 fluidViscosity(
const int phase,
224 const std::vector<PhasePresence>& cond,
225 const std::vector<int>& cells)
const;
228 fluidReciprocFVF(
const int phase,
231 const std::vector<PhasePresence>& cond,
232 const std::vector<int>& cells)
const;
235 fluidDensity(
const int phase,
238 const std::vector<PhasePresence>& cond,
239 const std::vector<int>& cells)
const;
242 poroMult(
const ADB& p)
const;
245 transMult(
const ADB& p)
const;
247 const std::vector<PhasePresence>
248 phaseCondition()
const {
return phaseCondition_; }
251 classifyCondition(
const PolymerBlackoilState& state);
256 #endif // OPM_FULLYIMPLICITCOMPRESSIBLEPOLYMERSOLVER_HEADER_INCLUDED
int newtonIterations() const
Definition: CompressibleTpfaPolymer.hpp:32
Definition: WellStateFullyImplicitBlackoilPolymer.hpp:28
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)
Definition: PolymerBlackoilState.hpp:33
Definition: PolymerPropsAd.hpp:32
Definition: FullyImplicitCompressiblePolymerSolver.hpp:53
int step(const double dt, PolymerBlackoilState &state, WellStateFullyImplicitBlackoilPolymer &wstate)
int linearIterations() const
parameter::ParameterGroup SolverParameters
Not used by this class except to satisfy interface requirements.
Definition: FullyImplicitCompressiblePolymerSolver.hpp:93