36#ifndef OPENRS_EULERUPSTREAMRESIDUAL_HEADER
37#define OPENRS_EULERUPSTREAMRESIDUAL_HEADER
41#include <opm/common/utility/parameters/ParameterGroup.hpp>
42#include <opm/common/utility/numeric/SparseVector.hpp>
47 namespace EulerUpstreamResidualDetails {
49 template <
class UpstreamSolver,
class PressureSolution>
56 template <
class Gr
idInterface,
class ReservoirProperties,
class BoundaryConditions>
60 template <
class S,
class P>
62 typedef typename GridInterface::CellIterator
CIt;
63 typedef typename CIt::FaceIterator
FIt;
65 typedef ReservoirProperties
RP;
74 const ReservoirProperties& resprop,
75 const BoundaryConditions& boundary);
78 const ReservoirProperties& resprop,
79 const BoundaryConditions& boundary);
82 template <
class FlowSolution>
85 const FlowSolution& flow_sol,
86 const Opm::SparseVector<double>& injection_rates,
87 const bool method_viscous,
88 const bool method_gravity,
89 const bool method_capillary,
90 std::vector<double>& sat_delta)
const;
94 const GridInterface&
grid()
const;
102 estimateCapPressureGradient(
const FIt& f,
const FIt& nbf)
const;
104 const GridInterface* pgrid_;
105 const ReservoirProperties* preservoir_properties_;
106 const BoundaryConditions* pboundary_;
110 std::vector<FIt> bid_to_face_;
113 std::vector<CIt> cell_iters_;
116 mutable std::vector<double> cap_pressures_;
117 mutable const Opm::SparseVector<double>* pinjection_rates_;
118 mutable bool method_viscous_;
119 mutable bool method_gravity_;
120 mutable bool method_capillary_;
Definition: EulerUpstreamResidual.hpp:58
EulerUpstreamResidual()
Definition: EulerUpstreamResidual_impl.hpp:356
void initObj(const GridInterface &grid, const ReservoirProperties &resprop, const BoundaryConditions &boundary)
Definition: EulerUpstreamResidual_impl.hpp:376
void computeResidual(const std::vector< double > &saturation, const typename GridInterface::Vector &gravity, const FlowSolution &flow_sol, const Opm::SparseVector< double > &injection_rates, const bool method_viscous, const bool method_gravity, const bool method_capillary, std::vector< double > &sat_delta) const
const ReservoirProperties & reservoirProperties() const
Definition: EulerUpstreamResidual_impl.hpp:430
FIt::Vector Vector
Definition: EulerUpstreamResidual.hpp:64
GridInterface::CellIterator CIt
Definition: EulerUpstreamResidual.hpp:62
EulerUpstreamResidual(const GridInterface &grid, const ReservoirProperties &resprop, const BoundaryConditions &boundary)
const GridInterface & grid() const
Definition: EulerUpstreamResidual_impl.hpp:422
const BoundaryConditions & boundaryConditions() const
Definition: EulerUpstreamResidual_impl.hpp:437
void computeCapPressures(const std::vector< double > &saturation) const
Definition: EulerUpstreamResidual_impl.hpp:445
ReservoirProperties RP
Definition: EulerUpstreamResidual.hpp:65
CIt::FaceIterator FIt
Definition: EulerUpstreamResidual.hpp:63
Dune::BlockVector< Dune::FieldVector< double, 1 > > Vector
A vector holding our RHS.
Definition: matrixops.hpp:33
Definition: ImplicitAssembly.hpp:43
Definition: EulerUpstreamResidual_impl.hpp:70