35#ifndef OPM_UPSCALERBASE_HEADER
36#define OPM_UPSCALERBASE_HEADER
38#include <opm/input/eclipse/Deck/Deck.hpp>
40#include <opm/common/utility/parameters/ParameterGroup.hpp>
42#include <opm/grid/CpGrid.hpp>
44#include <opm/porsol/common/GridInterfaceEuler.hpp>
45#include <opm/porsol/common/BoundaryConditions.hpp>
54 template <
class Traits>
60 typedef Dune::CpGrid GridType;
61 enum { Dimension = GridType::dimension };
63 typedef typename Traits::template ResProp<Dimension>::Type ResProp;
68 enum BoundaryConditionType { Fixed = 0, Linear = 1, Periodic = 2 };
78 void init(
const Opm::ParameterGroup& param);
81 void init(
const Opm::Deck& deck,
82 BoundaryConditionType bctype,
83 double perm_threshold,
84 double residual_tolerance = 1e-8,
85 int linsolver_verbosity = 0,
86 int linsolver_type = 3,
87 bool twodim_hack =
false,
88 int linsolver_maxit = 0,
89 double linsolver_prolongate_factor = 1.0,
90 int linsolver_smooth_steps = 1,
91 const double gravity = 0.0);
94 const GridType&
grid()
const;
131 typedef GridInterface::CellIterator CellIter;
132 typedef CellIter::FaceIterator FaceIter;
133 typedef BasicBoundaryConditions<true, true> BCs;
134 typedef typename Traits::template FlowSolver<GridInterface, BCs>::Type FlowSolver;
137 template <
class FlowSol>
138 double computeAverageVelocity(
const FlowSol& flow_solution,
140 const int pdrop_dir)
const;
142 double computeDelta(
const int flow_dir)
const;
144 template <
class Flu
idInterface>
145 permtensor_t upscaleEffectivePerm(
const FluidInterface& fluid);
147 virtual void initImpl(
const Opm::ParameterGroup& param);
149 virtual void initFinal(
const Opm::ParameterGroup& param);
152 BoundaryConditionType bctype_;
154 double residual_tolerance_;
155 int linsolver_maxit_;
156 double linsolver_prolongate_factor_;
157 int linsolver_verbosity_;
159 int linsolver_smooth_steps_;
163 GridInterface ginterf_;
166 FlowSolver flow_solver_;
171#include "UpscalerBase_impl.hpp"
Definition GridInterfaceEuler.hpp:403
A base class for upscaling.
Definition UpscalerBase.hpp:56
permtensor_t upscaleSinglePhase()
Does a single-phase upscaling.
Definition UpscalerBase_impl.hpp:219
void setPermeability(const int cell_index, const permtensor_t &k)
Set the permeability of a cell directly.
Definition UpscalerBase_impl.hpp:209
double upscalePorosity() const
Compute upscaled porosity.
Definition UpscalerBase_impl.hpp:374
ResProp::MutablePermTensor permtensor_t
A type for the upscaled permeability.
Definition UpscalerBase.hpp:66
const GridType & grid() const
Access the grid.
Definition UpscalerBase_impl.hpp:178
double upscaleNTG() const
Compute upscaled NTG.
Definition UpscalerBase_impl.hpp:400
double upscaleSOWCR(const bool NTG) const
Compute upscaled SOWCR.
Definition UpscalerBase_impl.hpp:432
void init(const Opm::ParameterGroup ¶m)
Initializes the upscaler from parameters.
Definition UpscalerBase_impl.hpp:65
void setBoundaryConditionType(BoundaryConditionType type)
Set boundary condition type.
Definition UpscalerBase_impl.hpp:188
double upscaleNetPorosity() const
Compute upscaled net porosity.
Definition UpscalerBase_impl.hpp:387
double upscaleSWCR(const bool NTG) const
Compute upscaled SWCR.
Definition UpscalerBase_impl.hpp:412
UpscalerBase()
Default constructor.
Definition UpscalerBase_impl.hpp:48
Inverting small matrices.
Definition ImplicitAssembly.hpp:43