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 <opm/grid/common/GridAdapter.hpp>
19#include <opm/grid/UnstructuredGrid.h>
20#include <opm/core/transport/implicit/NormSupport.hpp>
21#include <opm/core/transport/implicit/ImplicitTransport.hpp>
22#include <opm/common/utility/parameters/ParameterGroup.hpp>
24#include <opm/porsol/common/ReservoirPropertyCapillary.hpp>
33 template <
class Vector>
37 norm(
const Vector& v) {
38 return ImplicitTransportDefault::AccumulationNorm <Vector, ImplicitTransportDefault::MaxAbs>::norm(v);
42 template <
class Vector>
46 norm(
const Vector& v) {
47 return v.infinity_norm();
55 : press_ (g->number_of_cells),
56 fpress_(g->number_of_faces),
57 flux_ (g->number_of_faces),
58 sat_ (np * g->number_of_cells)
61 ::std::vector<double>& pressure () {
return press_ ; }
62 ::std::vector<double>& facepressure() {
return fpress_; }
63 ::std::vector<double>& faceflux () {
return flux_ ; }
64 ::std::vector<double>& saturation () {
return sat_ ; }
66 const ::std::vector<double>& faceflux ()
const {
return flux_; }
67 const ::std::vector<double>& saturation ()
const {
return sat_ ; }
70 ::std::vector<double> press_ ;
71 ::std::vector<double> fpress_;
72 ::std::vector<double> flux_ ;
73 ::std::vector<double> sat_ ;
78 Rock(::std::size_t nc, ::std::size_t dim)
80 perm_(nc * dim * dim),
83 const ::std::vector<double>& perm()
const {
return perm_; }
84 const ::std::vector<double>& poro()
const {
return poro_; }
87 perm_homogeneous(
double k) {
88 setVector(0.0, perm_);
90 const ::std::size_t d2 = dim_ * dim_;
92 for (::std::size_t c = 0, nc = poro_.size(); c < nc; ++c) {
93 for (::std::size_t i = 0; i < dim_; ++i) {
94 perm_[c*d2 + i*(dim_ + 1)] = k;
100 poro_homogeneous(
double phi) {
101 setVector(phi, poro_);
106 setVector(
double x, ::std::vector<double>& v) {
107 ::std::fill(v.begin(), v.end(), x);
111 ::std::vector<double> perm_;
112 ::std::vector<double> poro_;
121 template <
class Sat ,
125 mobility(
int c,
const Sat& s, Mob& mob, DMob& dmob)
const {
126 const double s1 = s[0];
129 r_.phaseMobilitiesDeriv(c, s1, dmob);
132 template <
class Sat ,
136 pc(
int c,
const Sat& s, PC& pc_arg, DPC& dpc)
const {
137 const double s1 = s[0];
142 double density(
int p)
const {
150 double s_min(
int c)
const {
154 double s_max(
int c)
const {
164 typedef Dune::FieldVector<double, 1> ScalarVectorBlockType;
165 typedef Dune::FieldMatrix<double, 1, 1> ScalarMatrixBlockType;
167 typedef Dune::BlockVector<ScalarVectorBlockType> ScalarBlockVector;
168 typedef Dune::BCRSMatrix <ScalarMatrixBlockType> ScalarBCRSMatrix;
171 solve(
const ScalarBCRSMatrix& A,
172 const ScalarBlockVector& b,
173 ScalarBlockVector& x)
175 Dune::MatrixAdapter<ScalarBCRSMatrix,
177 ScalarBlockVector> opA(A);
179 Dune::SeqILU<ScalarBCRSMatrix,
181 ScalarBlockVector> precond(A, 1.0);
187 Dune::BiCGSTABSolver<ScalarBlockVector>
188 solver(opA, precond, tol, maxit, verb);
190 ScalarBlockVector bcpy(b);
191 Dune::InverseOperatorResult res;
192 solver.apply(x, bcpy, res);
202 ::std::vector< int > cell;
203 ::std::vector<double> pressure;
204 ::std::vector<double> flux;
205 ::std::vector<double> saturation;
206 ::std::vector<double> periodic_cells;
207 ::std::vector<double> periodic_faces;
212 append_transport_source(
int c,
double p,
double v,
const Arr& s,
215 src.cell .push_back(c);
216 src.pressure .push_back(p);
217 src.flux .push_back(v);
218 src.saturation.insert(src.saturation.end(),
224 compute_porevolume(
const UnstructuredGrid* g,
226 std::vector<double>& porevol);
Definition ImplicitTransportDefs.hpp:162
Definition ImplicitTransportDefs.hpp:43
Definition ImplicitTransportDefs.hpp:34
A property class for incompressible two-phase flow.
Definition ReservoirPropertyCapillary.hpp:80
void phaseMobilities(int cell_index, double saturation, Vector &mobility) const
Mobilities for both phases.
Definition ReservoirPropertyCapillary_impl.hpp:93
double densitySecondPhase() const
Density of second (oil) phase.
Definition ReservoirPropertyCommon_impl.hpp:373
double densityFirstPhase() const
Density of first (water) phase.
Definition ReservoirPropertyCommon_impl.hpp:366
double capillaryPressure(int cell_index, double saturation) const
Capillary pressure.
Definition ReservoirPropertyCommon_impl.hpp:467
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 s_max(int c) const
Definition ReservoirPropertyCommon_impl.hpp:506
Definition ImplicitTransportDefs.hpp:52
A property class for porous media rock.
Definition ImplicitTransportDefs.hpp:76
Rock()
Default constructor.
Definition Rock_impl.hpp:37
Definition ImplicitTransportDefs.hpp:196
Definition ImplicitTransportDefs.hpp:115
Inverting small matrices.
Definition ImplicitAssembly.hpp:43