ImpesTPFAAD.hpp
1 /*
2  Copyright 2013 SINTEF ICT, Applied Mathematics.
3  Copyright 2013 Statoil ASA.
4 
5  This file is part of the Open Porous Media Project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 #ifndef OPM_IMPESTPFAAD_HEADER_INCLUDED
22 #define OPM_IMPESTPFAAD_HEADER_INCLUDED
23 
24 #include <opm/autodiff/AutoDiffBlock.hpp>
25 #include <opm/autodiff/AutoDiffHelpers.hpp>
26 #include <opm/autodiff/BlackoilModelEnums.hpp>
27 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
28 
29 struct UnstructuredGrid;
30 struct Wells;
31 
32 namespace Opm {
33 
34  class DerivedGeology;
35  class LinearSolverInterface;
36  class BlackoilState;
37  class WellState;
38 
45  {
46  public:
48  ImpesTPFAAD(const UnstructuredGrid& grid,
49  const BlackoilPropsAdFromDeck& fluid,
50  const DerivedGeology& geo,
51  const Wells& wells,
52  const LinearSolverInterface& linsolver);
53 
58  void solve(const double dt,
59  BlackoilState& state,
60  WellState& well_state);
61  private:
62  // Disallow copying and assignment
63  ImpesTPFAAD(const ImpesTPFAAD& rhs);
64  ImpesTPFAAD& operator=(const ImpesTPFAAD& rhs);
65 
66  // Types
67  typedef AutoDiffBlock<double> ADB;
68  typedef ADB::V V;
69  typedef ADB::M M;
70  typedef Eigen::Array<double,
71  Eigen::Dynamic,
72  Eigen::Dynamic,
73  Eigen::RowMajor> DataBlock;
74  enum { Water = Opm::Water,
75  Oil = Opm::Oil,
76  Gas = Opm::Gas };
77 
78  // Data
79  const UnstructuredGrid& grid_;
80  const BlackoilPropsAdFromDeck& fluid_;
81  const DerivedGeology& geo_ ;
82  const Wells& wells_;
83  const LinearSolverInterface& linsolver_;
84  HelperOps ops_;
85  const M grav_;
86  ADB cell_residual_;
87  std::vector<ADB> well_flow_residual_;
88  ADB well_residual_;
89  ADB total_residual_;
90  std::vector<V> kr_;
91  std::vector<V> well_kr_;
92  ADB qs_;
93  V well_perf_dp_;
94 
95  // Methods for assembling and solving.
96  void computeExplicitData(const double dt,
97  const BlackoilState& state,
98  const WellState& well_state);
99  void assemble(const double dt,
100  const BlackoilState& state,
101  const WellState& well_state);
102  void solveJacobianSystem(BlackoilState& state,
103  WellState& well_state) const;
104  double residualNorm() const;
105  void computeFluxes(BlackoilState& state, WellState& well_state) const;
106 
107  // Fluid interface forwarding calls to correct methods of fluid_.
108  V fluidMu(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
109  ADB fluidMu(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
110  V fluidFvf(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
111  ADB fluidFvf(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
112  V fluidRho(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
113  ADB fluidRho(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
114  std::vector<V> fluidRelperm(const V& sw, const V& so, const V& sg, const std::vector<int>& cells) const;
115  V fluidKr(const int phase) const;
116  V fluidKrWell(const int phase) const;
117  };
118 
119 
120 } // namespace Opm
121 
122 #endif /* OPM_IMPESTPFAAD_HEADER_INCLUDED */
Contains vectors and sparse matrices that represent subsets or operations on (AD or regular) vectors ...
Definition: AutoDiffHelpers.hpp:44
AutoDiffMatrix is a wrapper class that optimizes matrix operations.
Definition: AutoDiffMatrix.hpp:43
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
void solve(const double dt, BlackoilState &state, WellState &well_state)
Solve forward in time.
Definition: ImpesTPFAAD.cpp:172
ImpesTPFAAD(const UnstructuredGrid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const Wells &wells, const LinearSolverInterface &linsolver)
Construct impes solver.
Definition: ImpesTPFAAD.cpp:146
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
Eigen::Array< Scalar, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99
Class for solving black-oil impes problems.
Definition: ImpesTPFAAD.hpp:44