00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020 #ifndef OPM_TRANSPORTSOLVERTWOPHASECOMPRESSIBLEPOLYMER_HEADER_INCLUDED
00021 #define OPM_TRANSPORTSOLVERTWOPHASECOMPRESSIBLEPOLYMER_HEADER_INCLUDED
00022
00023 #include <opm/polymer/PolymerProperties.hpp>
00024 #include <opm/core/transport/reorder/ReorderSolverInterface.hpp>
00025 #include <opm/core/utility/linearInterpolation.hpp>
00026 #include <vector>
00027 #include <list>
00028
00029 struct UnstructuredGrid;
00030
00031 namespace {
00032 class ResSOnCurve;
00033 class ResCOnCurve;
00034 }
00035
00036 namespace Opm
00037 {
00038
00039 class BlackoilPropertiesInterface;
00040
00044 class TransportSolverTwophaseCompressiblePolymer : public ReorderSolverInterface
00045 {
00046 public:
00047
00048 enum SingleCellMethod { Bracketing, Newton, NewtonC, Gradient};
00049 enum GradientMethod { Analytic, FinDif };
00050
00062 TransportSolverTwophaseCompressiblePolymer(const UnstructuredGrid& grid,
00063 const BlackoilPropertiesInterface& props,
00064 const PolymerProperties& polyprops,
00065 const SingleCellMethod method,
00066 const double tol,
00067 const int maxit);
00068
00070 void setPreferredMethod(SingleCellMethod method);
00071
00092 void solve(const double* darcyflux,
00093 const std::vector<double>& initial_pressure,
00094 const std::vector<double>& pressure,
00095 const std::vector<double>& temperature,
00096 const double* porevolume0,
00097 const double* porevolume,
00098 const double* source,
00099 const double* polymer_inflow_c,
00100 const double dt,
00101 std::vector<double>& saturation,
00102 std::vector<double>& surfacevol,
00103 std::vector<double>& concentration,
00104 std::vector<double>& cmax);
00105
00108 void initGravity(const double* grav);
00109
00121 void solveGravity(const std::vector<std::vector<int> >& columns,
00122 const double dt,
00123 std::vector<double>& saturation,
00124 std::vector<double>& surfacevol,
00125 std::vector<double>& concentration,
00126 std::vector<double>& cmax);
00127
00128
00129
00130
00131
00132 private:
00133
00134 const UnstructuredGrid& grid_;
00135 const BlackoilPropertiesInterface& props_;
00136 const PolymerProperties& polyprops_;
00137 const double* darcyflux_;
00138 const double* porevolume0_;
00139 const double* porevolume_;
00140 const double* source_;
00141 const double* polymer_inflow_c_;
00142 double dt_;
00143 double tol_;
00144 double maxit_;
00145 SingleCellMethod method_;
00146 double adhoc_safety_;
00147
00148 std::vector<double> saturation_;
00149 std::vector<int> allcells_;
00150 double* concentration_;
00151 double* cmax_;
00152 std::vector<double> fractionalflow_;
00153 std::vector<double> mc_;
00154 std::vector<double> visc_;
00155 std::vector<double> A_;
00156 std::vector<double> A0_;
00157 std::vector<double> smin_;
00158 std::vector<double> smax_;
00159
00160
00161 const double* gravity_;
00162 std::vector<double> trans_;
00163 std::vector<double> density_;
00164 std::vector<double> gravflux_;
00165 std::vector<double> mob_;
00166 std::vector<double> cmax0_;
00167
00168
00169 std::vector<double> s0_;
00170 std::vector<double> c0_;
00171
00172
00173 std::vector<int> ia_upw_;
00174 std::vector<int> ja_upw_;
00175 std::vector<int> ia_downw_;
00176 std::vector<int> ja_downw_;
00177
00178 struct ResidualC;
00179 struct ResidualS;
00180
00181 class ResidualCGrav;
00182 class ResidualSGrav;
00183
00184 class ResidualEquation;
00185 class ResSOnCurve;
00186 class ResCOnCurve;
00187
00188 friend class TransportSolverTwophaseCompressiblePolymer::ResidualEquation;
00189 friend class TransportSolverTwophaseCompressiblePolymer::ResSOnCurve;
00190 friend class TransportSolverTwophaseCompressiblePolymer::ResCOnCurve;
00191
00192
00193 virtual void solveSingleCell(const int cell);
00194 virtual void solveMultiCell(const int num_cells, const int* cells);
00195 void solveSingleCellBracketing(int cell);
00196 void solveSingleCellNewton(int cell, bool use_sc, bool use_explicit_step = false);
00197 void solveSingleCellGradient(int cell);
00198 void solveSingleCellGravity(const std::vector<int>& cells,
00199 const int pos,
00200 const double* gravflux);
00201 int solveGravityColumn(const std::vector<int>& cells);
00202
00203 void initGravityDynamic();
00204
00205 void fracFlow(double s, double c, double cmax, int cell, double& ff) const;
00206 void fracFlowWithDer(double s, double c, double cmax, int cell, double& ff,
00207 double* dff_dsdc) const;
00208 void fracFlowBoth(double s, double c, double cmax, int cell, double& ff,
00209 double* dff_dsdc, bool if_with_der) const;
00210 void computeMc(double c, double& mc) const;
00211 void computeMcWithDer(double c, double& mc, double& dmc_dc) const;
00212 void mobility(double s, double c, int cell, double* mob) const;
00213 void scToc(const double* x, double* x_c) const;
00214 #ifdef PROFILING
00215 class Newton_Iter {
00216 public:
00217 bool res_s;
00218 int cell;
00219 double s;
00220 double c;
00221
00222 Newton_Iter(bool res_s_val, int cell_val, double s_val, double c_val) {
00223 res_s = res_s_val;
00224 cell = cell_val;
00225 s = s_val;
00226 c = c_val;
00227 }
00228 };
00229
00230 std::list<Newton_Iter> res_counts;
00231 #endif
00232 };
00233
00234 }
00235
00236 #endif // OPM_TRANSPORTSOLVERTWOPHASECOMPRESSIBLEPOLYMER_HEADER_INCLUDED