My Project
EulerUpstreamImplicit.hpp
1 //===========================================================================
2 //
3 // File: EulerUpstreamImplicit.hpp
4 //
5 // Created: Tue Jun 16 14:07:53 2009
6 //
7 // Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8 // B�rd Skaflestad <bard.skaflestad@sintef.no>
9 //
10 // $Date$
11 //
12 // $Revision$
13 //
14 //===========================================================================
15 
16 /*
17 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
18 Copyright 2009, 2010 Statoil ASA.
19 
20 This file is part of The Open Reservoir Simulator Project (OpenRS).
21 
22 OpenRS is free software: you can redistribute it and/or modify
23 it under the terms of the GNU General Public License as published by
24 the Free Software Foundation, either version 3 of the License, or
25 (at your option) any later version.
26 
27 OpenRS is distributed in the hope that it will be useful,
28 but WITHOUT ANY WARRANTY; without even the implied warranty of
29 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 GNU General Public License for more details.
31 
32 You should have received a copy of the GNU General Public License
33 along with OpenRS. If not, see <http://www.gnu.org/licenses/>.
34 */
35 
36 #ifndef OPENRS_EULERUPSTREAMIMPLICIT_HEADER
37 #define OPENRS_EULERUPSTREAMIMPLICIT_HEADER
38 
39 #include <opm/common/utility/numeric/SparseVector.hpp>
40 #include <opm/common/utility/parameters/ParameterGroup.hpp>
41 
42 #include <opm/porsol/common/ImplicitTransportDefs.hpp>
43 #include <opm/porsol/common/BCRSMatrixBlockAssembler.hpp>
44 
45 #include <opm/grid/common/GridAdapter.hpp>
46 
47 #include <opm/core/transport/implicit/ImplicitAssembly.hpp>
48 #include <opm/core/transport/implicit/ImplicitTransport.hpp>
49 #include <opm/core/transport/implicit/JacobianSystem.hpp>
51 
52 
53 namespace Opm {
54 
57  template <class GridInterface, class ReservoirProperties, class BoundaryConditions>
59  {
60  public:
67  EulerUpstreamImplicit (const GridInterface& grid,
68  const ReservoirProperties& resprop,
69  const BoundaryConditions& boundary);
73  void init(const Opm::ParameterGroup& param);
77  void init(const Opm::ParameterGroup& param,
78  const GridInterface& grid,
79  const ReservoirProperties& resprop,
80  const BoundaryConditions& boundary);
84  void initObj(const GridInterface& grid,
85  const ReservoirProperties& resprop,
86  const BoundaryConditions& boundary);
90  void display();
91 
96  template <class PressureSolution>
97  bool transportSolve(std::vector<double>& saturation,
98  const double time,
99  const typename GridInterface::Vector& gravity,
100  const PressureSolution& pressure_sol,
101  const Opm::SparseVector<double>& injection_rates) const;
102 
103  protected:
104  // typedef typename GridInterface::CellIterator CIt;
105  // typedef typename CIt::FaceIterator FIt;
106  // typedef typename FIt::Vector Vector;
107 
108  template <class PressureSolution>
109  void smallTimeStep(std::vector<double>& saturation,
110  const double time,
111  const typename GridInterface::Vector& gravity,
112  const PressureSolution& pressure_sol,
113  const Opm::SparseVector<double>& injection_rates) const;
114 
115  void checkAndPossiblyClampSat(std::vector<double>& s) const;
116 
117 
118 
119  // Interface to implicit solver
122  typedef Dune::FieldVector<double, 1> ScalarVectorBlockType;
123  typedef Dune::FieldMatrix<double, 1, 1> ScalarMatrixBlockType;
124 
125  typedef Dune::BlockVector<ScalarVectorBlockType> ScalarBlockVector;
126  typedef Dune::BCRSMatrix <ScalarMatrixBlockType> ScalarBCRSMatrix;
127 
130 
133  JacSys ,
139 
140  // should be initialized by param
141  GridAdapter mygrid_;
142  ReservoirProperties myrp_;
143 
144  bool check_sat_;
145  bool clamp_sat_;
146  int max_repeats_;
147  std::vector<double> porevol_;
148 
149  std::vector<int> periodic_cells_;
150  std::vector<int> periodic_faces_;
151  std::vector<int> periodic_nbfaces_;
152  std::vector<int> periodic_hfaces_;
153  std::vector<int> direclet_cells_;
154  std::vector<double> direclet_sat_;
155  std::vector<double> direclet_hfaces_;
156  std::vector<double> htrans_;
158  };
159 
160 } // namespace Opm
161 
162 #include "EulerUpstreamImplicit_impl.hpp"
163 
164 #endif // OPENRS_EULERUPSTREAM_HEADER
Numerical model and support classes needed to model transport of two incompressible fluid phases.
Class for doing simple transport by implicit Euler upstream method for general grid.
Definition: EulerUpstreamImplicit.hpp:59
EulerUpstreamImplicit()
Definition: EulerUpstreamImplicit_impl.hpp:60
void init(const Opm::ParameterGroup &param, const GridInterface &grid, const ReservoirProperties &resprop, const BoundaryConditions &boundary)
void initObj(const GridInterface &grid, const ReservoirProperties &resprop, const BoundaryConditions &boundary)
Definition: EulerUpstreamImplicit_impl.hpp:98
bool transportSolve(std::vector< double > &saturation, const double time, const typename GridInterface::Vector &gravity, const PressureSolution &pressure_sol, const Opm::SparseVector< double > &injection_rates) const
Solve transport equation, evolving.
void display()
Definition: EulerUpstreamImplicit_impl.hpp:217
void init(const Opm::ParameterGroup &param)
Definition: EulerUpstreamImplicit_impl.hpp:74
EulerUpstreamImplicit(const GridInterface &grid, const ReservoirProperties &resprop, const BoundaryConditions &boundary)
Definition: JacobianSystem.hpp:222
Definition: JacobianSystem.hpp:108
Definition: JacobianSystem.hpp:89
Definition: JacobianSystem.hpp:65
Definition: JacobianSystem.hpp:78
Definition: ImplicitTransport.hpp:105
Definition: ImplicitTransportDefs.hpp:162
Definition: ImplicitTransportDefs.hpp:43
Definition: SinglePointUpwindTwoPhase.hpp:278
Definition: ImplicitTransportDefs.hpp:115
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43
Definition: ImplicitTransport.hpp:45