My Project
EulerUpstreamResidual.hpp
1 //===========================================================================
2 //
3 // File: EulerUpstreamResidual.hpp
4 //
5 // Created: Thu May 6 11:14:23 2010
6 //
7 // Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8 // Jostein R Natvig <jostein.r.natvig@sintef.no>
9 //
10 // $Date$
11 //
12 // $Revision$
13 //
14 //===========================================================================
15 
16 /*
17  Copyright 2010 SINTEF ICT, Applied Mathematics.
18  Copyright 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_EULERUPSTREAMRESIDUAL_HEADER
37 #define OPENRS_EULERUPSTREAMRESIDUAL_HEADER
38 
39 
40 
41 #include <opm/common/utility/parameters/ParameterGroup.hpp>
42 #include <opm/common/utility/numeric/SparseVector.hpp>
43 
44 
45 namespace Opm {
46 
47  namespace EulerUpstreamResidualDetails {
48  // Forward declaration for friendship purposes.
49  template <class UpstreamSolver, class PressureSolution>
50  struct UpdateForCell;
51  }
52 
53 
56  template <class GridInterface, class ReservoirProperties, class BoundaryConditions>
58  {
59  public:
60  template <class S, class P>
62  typedef typename GridInterface::CellIterator CIt;
63  typedef typename CIt::FaceIterator FIt;
64  typedef typename FIt::Vector Vector;
65  typedef ReservoirProperties RP;
66 
73  EulerUpstreamResidual(const GridInterface& grid,
74  const ReservoirProperties& resprop,
75  const BoundaryConditions& boundary);
76 
77  void initObj(const GridInterface& grid,
78  const ReservoirProperties& resprop,
79  const BoundaryConditions& boundary);
80 
81 
82  template <class FlowSolution>
83  void computeResidual(const std::vector<double>& saturation,
84  const typename GridInterface::Vector& gravity,
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;
91 
92  void computeCapPressures(const std::vector<double>& saturation) const;
93 
94  const GridInterface& grid() const;
95  const ReservoirProperties& reservoirProperties() const;
96  const BoundaryConditions& boundaryConditions() const;
97 
98  private:
99  void initFinal();
100 
101  typename GridInterface::Vector
102  estimateCapPressureGradient(const FIt& f, const FIt& nbf) const;
103 
104  const GridInterface* pgrid_;
105  const ReservoirProperties* preservoir_properties_;
106  const BoundaryConditions* pboundary_;
107 
108  // Boundary id to face iterator mapping. May be mostly or completely empty.
109  // Obviously requires unique-face-per-bid grids.
110  std::vector<FIt> bid_to_face_;
111 
112  // Storing some cell iterators, so that we may use tbb for parallelizing.
113  std::vector<CIt> cell_iters_;
114 
115  // Precomputing the capillary pressures of cells saves a little time.
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_;
121  };
122 
123 } // namespace Opm
124 
125 #include "EulerUpstreamResidual_impl.hpp"
126 
127 
128 #endif // OPENRS_EULERUPSTREAMRESIDUAL_HEADER
Class for doing simple transport by explicit Euler upstream method for general grid.
Definition: EulerUpstreamResidual.hpp:58
EulerUpstreamResidual()
Definition: EulerUpstreamResidual_impl.hpp:356
EulerUpstreamResidual(const GridInterface &grid, const ReservoirProperties &resprop, const BoundaryConditions &boundary)
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43
Definition: EulerUpstreamResidual_impl.hpp:70