00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020 #ifndef OPM_TRANSPORTSOLVERTWOPHASEPOLYMER_HEADER_INCLUDED
00021 #define OPM_TRANSPORTSOLVERTWOPHASEPOLYMER_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 Opm
00032 {
00033
00034 class IncompPropertiesInterface;
00035
00039 class TransportSolverTwophasePolymer : public ReorderSolverInterface
00040 {
00041 public:
00042
00043 enum SingleCellMethod { Bracketing, Newton, Gradient, NewtonSimpleSC, NewtonSimpleC};
00044 enum GradientMethod { Analytic, FinDif };
00045
00056 TransportSolverTwophasePolymer(const UnstructuredGrid& grid,
00057 const IncompPropertiesInterface& props,
00058 const PolymerProperties& polyprops,
00059 const SingleCellMethod method,
00060 const double tol,
00061 const int maxit);
00062
00064 void setPreferredMethod(SingleCellMethod method);
00065
00080 void solve(const double* darcyflux,
00081 const double* porevolume,
00082 const double* source,
00083 const double* polymer_inflow_c,
00084 const double dt,
00085 std::vector<double>& saturation,
00086 std::vector<double>& concentration,
00087 std::vector<double>& cmax);
00088
00100 void solveGravity(const std::vector<std::vector<int> >& columns,
00101 const double* porevolume,
00102 const double dt,
00103 std::vector<double>& saturation,
00104 std::vector<double>& concentration,
00105 std::vector<double>& cmax);
00106
00107 public:
00108 virtual void solveSingleCell(const int cell);
00109 virtual void solveMultiCell(const int num_cells, const int* cells);
00110 void solveSingleCellBracketing(int cell);
00111 void solveSingleCellNewton(int cell);
00112 void solveSingleCellGradient(int cell);
00113 void solveSingleCellNewtonSimple(int cell,bool use_sc);
00114 class ResidualEquation;
00115
00116 void initGravity(const double* grav);
00117 void solveSingleCellGravity(const std::vector<int>& cells,
00118 const int pos,
00119 const double* gravflux);
00120 int solveGravityColumn(const std::vector<int>& cells);
00121 void scToc(const double* x, double* x_c) const;
00122
00123 #ifdef PROFILING
00124 class Newton_Iter {
00125 public:
00126 bool res_s;
00127 int cell;
00128 double s;
00129 double c;
00130
00131 Newton_Iter(bool res_s_val, int cell_val, double s_val, double c_val) {
00132 res_s = res_s_val;
00133 cell = cell_val;
00134 s = s_val;
00135 c = c_val;
00136 }
00137 };
00138
00139 std::list<Newton_Iter> res_counts;
00140 #endif
00141
00142
00143 private:
00144 const UnstructuredGrid& grid_;
00145 const double* porosity_;
00146 const double* porevolume_;
00147 const IncompPropertiesInterface& props_;
00148 const PolymerProperties& polyprops_;
00149 std::vector<double> smin_;
00150 std::vector<double> smax_;
00151 double tol_;
00152 double maxit_;
00153
00154 const double* darcyflux_;
00155 const double* source_;
00156 const double* polymer_inflow_c_;
00157 double dt_;
00158 std::vector<double> saturation_;
00159 double* concentration_;
00160 double* cmax_;
00161 std::vector<double> fractionalflow_;
00162 std::vector<double> mc_;
00163 const double* visc_;
00164 SingleCellMethod method_;
00165 double adhoc_safety_;
00166
00167
00168 std::vector<double> gravflux_;
00169 std::vector<double> mob_;
00170 std::vector<double> cmax0_;
00171
00172 std::vector<double> s0_;
00173 std::vector<double> c0_;
00174
00175 struct ResidualC;
00176 struct ResidualS;
00177
00178 class ResidualCGrav;
00179 class ResidualSGrav;
00180
00181
00182 void fracFlow(double s, double c, double cmax, int cell, double& ff) const;
00183 void fracFlowWithDer(double s, double c, double cmax, int cell, double& ff,
00184 double* dff_dsdc) const;
00185 void fracFlowBoth(double s, double c, double cmax, int cell, double& ff,
00186 double* dff_dsdc, bool if_with_der) const;
00187 void computeMc(double c, double& mc) const;
00188 void computeMcWithDer(double c, double& mc, double& dmc_dc) const;
00189 void mobility(double s, double c, int cell, double* mob) const;
00190 };
00191
00192 }
00193
00194 #endif // OPM_TRANSPORTSOLVERTWOPHASEPOLYMER_HEADER_INCLUDED