1 #ifndef OPENRS_IMPLICITTRANSPORTDEFS_HEADER 2 #define OPENRS_IMPLICITTRANSPORTDEFS_HEADER 4 #include <opm/common/utility/platform_dependent/disable_warnings.h> 6 #include <dune/common/fvector.hh> 7 #include <dune/common/fmatrix.hh> 9 #include <dune/istl/operators.hh> 10 #include <dune/istl/solvers.hh> 11 #include <dune/istl/bvector.hh> 12 #include <dune/istl/bcrsmatrix.hh> 13 #include <dune/istl/preconditioners.hh> 15 #include <opm/common/utility/platform_dependent/reenable_warnings.h> 17 #include <dune/grid/common/GridAdapter.hpp> 19 #include <opm/core/grid.h> 20 #include <opm/core/linalg/LinearSolverIstl.hpp> 21 #include <opm/core/linalg/sparse_sys.h> 22 #include <opm/core/pressure/tpfa/ifs_tpfa.h> 23 #include <opm/core/pressure/tpfa/trans_tpfa.h> 24 #include <opm/core/transport/implicit/NormSupport.hpp> 25 #include <opm/core/transport/implicit/ImplicitTransport.hpp> 26 #include <opm/core/utility/parameters/ParameterGroup.hpp> 28 #include <opm/porsol/common/ReservoirPropertyCapillary.hpp> 36 #include <boost/lexical_cast.hpp> 39 template <
class Vector>
43 norm(
const Vector& v) {
44 return ImplicitTransportDefault::AccumulationNorm <Vector, ImplicitTransportDefault::MaxAbs>::norm(v);
48 template <
class Vector>
52 norm(
const Vector& v) {
53 return v.infinity_norm();
61 : press_ (g->number_of_cells),
62 fpress_(g->number_of_faces),
63 flux_ (g->number_of_faces),
64 sat_ (np * g->number_of_cells)
67 ::std::vector<double>& pressure () {
return press_ ; }
68 ::std::vector<double>& facepressure() {
return fpress_; }
69 ::std::vector<double>& faceflux () {
return flux_ ; }
70 ::std::vector<double>& saturation () {
return sat_ ; }
72 const ::std::vector<double>& faceflux ()
const {
return flux_; }
73 const ::std::vector<double>& saturation ()
const {
return sat_ ; }
76 ::std::vector<double> press_ ;
77 ::std::vector<double> fpress_;
78 ::std::vector<double> flux_ ;
79 ::std::vector<double> sat_ ;
84 Rock(::std::size_t nc, ::std::size_t dim)
86 perm_(nc * dim * dim),
89 const ::std::vector<double>& perm()
const {
return perm_; }
90 const ::std::vector<double>& poro()
const {
return poro_; }
93 perm_homogeneous(
double k) {
94 setVector(0.0, perm_);
96 const ::std::size_t d2 = dim_ * dim_;
98 for (::std::size_t c = 0, nc = poro_.size(); c < nc; ++c) {
99 for (::std::size_t i = 0; i < dim_; ++i) {
100 perm_[c*d2 + i*(dim_ + 1)] = k;
106 poro_homogeneous(
double phi) {
107 setVector(phi, poro_);
112 setVector(
double x, ::std::vector<double>& v) {
113 ::std::fill(v.begin(), v.end(), x);
117 ::std::vector<double> perm_;
118 ::std::vector<double> poro_;
127 template <
class Sat ,
131 mobility(
int c,
const Sat& s, Mob& mob, DMob& dmob)
const {
132 const double s1 = s[0];
135 r_.phaseMobilitiesDeriv(c, s1, dmob);
138 template <
class Sat ,
142 pc(
int c,
const Sat& s, PC& pc, DPC& dpc)
const {
143 const double s1 = s[0];
148 double density(
int p)
const {
156 double s_min(
int c)
const {
160 double s_max(
int c)
const {
170 typedef Dune::FieldVector<double, 1> ScalarVectorBlockType;
171 typedef Dune::FieldMatrix<double, 1, 1> ScalarMatrixBlockType;
173 typedef Dune::BlockVector<ScalarVectorBlockType> ScalarBlockVector;
174 typedef Dune::BCRSMatrix <ScalarMatrixBlockType> ScalarBCRSMatrix;
177 solve(
const ScalarBCRSMatrix& A,
178 const ScalarBlockVector& b,
179 ScalarBlockVector& x)
181 Dune::MatrixAdapter<ScalarBCRSMatrix,
183 ScalarBlockVector> opA(A);
185 Dune::SeqILU0<ScalarBCRSMatrix,
187 ScalarBlockVector> precond(A, 1.0);
193 Dune::BiCGSTABSolver<ScalarBlockVector>
194 solver(opA, precond, tol, maxit, verb);
196 ScalarBlockVector bcpy(b);
197 Dune::InverseOperatorResult res;
198 solver.apply(x, bcpy, res);
206 Opm::ParameterGroup params;
208 params.insertParameter(
"linsolver_tolerance",
209 boost::lexical_cast<std::string>(5.0e-9));
211 params.insertParameter(
"linsolver_verbosity",
212 boost::lexical_cast<std::string>(1));
214 params.insertParameter(
"linsolver_type",
215 boost::lexical_cast<std::string>(1));
217 ls_.reset(
new LinearSolverIstl(params));
221 solve(
struct CSRMatrix* A,
225 ls_->solve(A->m, A->nnz, A->ia, A->ja, A->sa, b, x);
229 std::unique_ptr<Opm::LinearSolverIstl> ls_;
238 ::std::vector< int > cell;
239 ::std::vector<double> pressure;
240 ::std::vector<double> flux;
241 ::std::vector<double> saturation;
242 ::std::vector<double> periodic_cells;
243 ::std::vector<double> periodic_faces;
248 append_transport_source(
int c,
double p,
double v,
const Arr& s,
251 src.cell .push_back(c);
252 src.pressure .push_back(p);
253 src.flux .push_back(v);
254 src.saturation.insert(src.saturation.end(),
260 compute_porevolume(
const UnstructuredGrid* g,
262 std::vector<double>& porevol);
265 #endif // OPENRS_IMPLICITTRANSPORTDEFS_HEADER Definition: ImplicitTransportDefs.hpp:202
Definition: ImplicitTransportDefs.hpp:168
double capillaryPressureDeriv(int c, double s) const
Derivative of Capillary pressure.
Definition: ReservoirPropertyCommon_impl.hpp:480
double s_min(int c) const
Definition: ReservoirPropertyCommon_impl.hpp:493
double densitySecondPhase() const
Density of second (oil) phase.
Definition: ReservoirPropertyCommon_impl.hpp:373
Definition: ImplicitTransportDefs.hpp:40
Definition: ImplicitTransportDefs.hpp:121
Class for immiscible dead oil and dry gas.
Definition: applier.hpp:18
Rock()
Default constructor.
Definition: Rock_impl.hpp:37
A property class for porous media rock.
Definition: ImplicitTransportDefs.hpp:82
double s_max(int c) const
Definition: ReservoirPropertyCommon_impl.hpp:506
Definition: ImplicitTransportDefs.hpp:58
Definition: ImplicitTransportDefs.hpp:232
Definition: ImplicitTransportDefs.hpp:49
double capillaryPressure(int cell_index, double saturation) const
Capillary pressure.
Definition: ReservoirPropertyCommon_impl.hpp:467
void phaseMobilities(int cell_index, double saturation, Vector &mobility) const
Mobilities for both phases.
Definition: ReservoirPropertyCapillary_impl.hpp:93
double densityFirstPhase() const
Density of first (water) phase.
Definition: ReservoirPropertyCommon_impl.hpp:366