00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #ifndef OPM_IMPESTPFAAD_HEADER_INCLUDED
00022 #define OPM_IMPESTPFAAD_HEADER_INCLUDED
00023
00024 #include <opm/autodiff/AutoDiffBlock.hpp>
00025 #include <opm/autodiff/AutoDiffHelpers.hpp>
00026 #include <opm/autodiff/BlackoilModelEnums.hpp>
00027 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
00028
00029 struct UnstructuredGrid;
00030 struct Wells;
00031
00032 namespace Opm {
00033
00034 class DerivedGeology;
00035 class LinearSolverInterface;
00036 class BlackoilState;
00037 class WellState;
00038
00044 class ImpesTPFAAD
00045 {
00046 public:
00048 ImpesTPFAAD(const UnstructuredGrid& grid,
00049 const BlackoilPropsAdFromDeck& fluid,
00050 const DerivedGeology& geo,
00051 const Wells& wells,
00052 const LinearSolverInterface& linsolver);
00053
00058 void solve(const double dt,
00059 BlackoilState& state,
00060 WellState& well_state);
00061 private:
00062
00063 ImpesTPFAAD(const ImpesTPFAAD& rhs);
00064 ImpesTPFAAD& operator=(const ImpesTPFAAD& rhs);
00065
00066
00067 typedef AutoDiffBlock<double> ADB;
00068 typedef ADB::V V;
00069 typedef ADB::M M;
00070 typedef Eigen::Array<double,
00071 Eigen::Dynamic,
00072 Eigen::Dynamic,
00073 Eigen::RowMajor> DataBlock;
00074 enum { Water = Opm::Water,
00075 Oil = Opm::Oil,
00076 Gas = Opm::Gas };
00077
00078
00079 const UnstructuredGrid& grid_;
00080 const BlackoilPropsAdFromDeck& fluid_;
00081 const DerivedGeology& geo_ ;
00082 const Wells& wells_;
00083 const LinearSolverInterface& linsolver_;
00084 HelperOps ops_;
00085 const M grav_;
00086 ADB cell_residual_;
00087 std::vector<ADB> well_flow_residual_;
00088 ADB well_residual_;
00089 ADB total_residual_;
00090 std::vector<V> kr_;
00091 std::vector<V> well_kr_;
00092 ADB qs_;
00093 V well_perf_dp_;
00094
00095
00096 void computeExplicitData(const double dt,
00097 const BlackoilState& state,
00098 const WellState& well_state);
00099 void assemble(const double dt,
00100 const BlackoilState& state,
00101 const WellState& well_state);
00102 void solveJacobianSystem(BlackoilState& state,
00103 WellState& well_state) const;
00104 double residualNorm() const;
00105 void computeFluxes(BlackoilState& state, WellState& well_state) const;
00106
00107
00108 V fluidMu(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
00109 ADB fluidMu(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
00110 V fluidFvf(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
00111 ADB fluidFvf(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
00112 V fluidRho(const int phase, const V& p, const V& T, const std::vector<int>& cells) const;
00113 ADB fluidRho(const int phase, const ADB& p, const ADB& T, const std::vector<int>& cells) const;
00114 std::vector<V> fluidRelperm(const V& sw, const V& so, const V& sg, const std::vector<int>& cells) const;
00115 V fluidKr(const int phase) const;
00116 V fluidKrWell(const int phase) const;
00117 };
00118
00119
00120 }
00121
00122 #endif