36 #ifndef OPM_STEADYSTATEUPSCALERMANAGERIMPLICIT_HEADER
37 #define OPM_STEADYSTATEUPSCALERMANAGERIMPLICIT_HEADER
42 #include <opm/core/utility/Units.hpp>
43 #include <opm/core/utility/SparseTable.hpp>
54 template <
class Istream>
55 void readControl(Istream& is, std::vector<double>& saturations, Opm::SparseTable<double>& all_pdrops)
59 std::vector<double> sat(num_sats);
60 Opm::SparseTable<double> ap;
61 std::vector<double> tmp_pd;
62 for (
int i = 0; i < num_sats; ++i) {
66 tmp_pd.resize(num_pdrops);
67 for (
int j = 0;
j < num_pdrops; ++
j) {
70 ap.appendRow(tmp_pd.begin(), tmp_pd.end());
73 saturations.swap(sat);
80 template <
class Ostream>
81 void writeControl(Ostream& os,
const std::vector<double>& saturations,
const Opm::SparseTable<double>& all_pdrops)
83 int num_sats = saturations.size();
84 os << num_sats <<
'\n';
85 for (
int i = 0; i < num_sats; ++i) {
86 os << saturations[i] <<
'\n';
87 int num_pdrops = all_pdrops[i].size();
89 for (
int j = 0;
j < num_pdrops; ++
j) {
90 os <<
' ' << all_pdrops[i][
j];
99 template<
class Ostream,
class Tensor>
100 void writeRelPerm(Ostream& os,
const Tensor& K,
double sat,
double pdrop,
bool success)
119 os << K(0,0) <<
'\t';
120 os << K(1,1) <<
'\t';
121 os << K(2,2) <<
'\t';
122 os << K(1,2) <<
'\t';
123 os << K(0,2) <<
'\t';
124 os << K(0,1) <<
'\t';
125 os << K(2,1) <<
'\t';
126 os << K(2,0) <<
'\t';
127 os << K(1,2) <<
'\t';
136 template <
class Upscaler>
140 void upscale(
const Opm::parameter::ParameterGroup& param)
143 std::vector<double> saturations;
144 Opm::SparseTable<double> all_pdrops;
145 bool from_file = param.has(
"sat_pdrop_filename");
147 std::string filename = param.get<std::string>(
"sat_pdrop_filename");
148 std::ifstream file(filename.c_str());
150 OPM_THROW(std::runtime_error,
"Could not open file " << filename);
155 int num_sats = param.getDefault(
"num_sats", 4);
156 double min_sat = param.getDefault(
"min_sat", 0.2);
157 double max_sat = param.getDefault(
"max_sat", 0.8);
158 saturations.resize(num_sats);
159 for (
int i = 0; i < num_sats; ++i) {
160 double factor = num_sats == 1 ? 0 : double(i)/double(num_sats - 1);
161 saturations[i] = (1.0 - factor)*min_sat + factor*max_sat;
164 int num_pdrops = param.getDefault(
"num_pdrops", 5);
165 double log_min_pdrop = std::log(param.getDefault(
"min_pdrop", 1e2));
166 double log_max_pdrop = std::log(param.getDefault(
"max_pdrop", 1e6));
167 std::vector<double> pdrops;
168 pdrops.resize(num_pdrops);
169 for (
int i = 0; i < num_pdrops; ++i) {
170 double factor = num_pdrops == 1 ? 0 : double(i)/double(num_pdrops - 1);
171 pdrops[i] = std::exp((1.0 - factor)*log_min_pdrop + factor*log_max_pdrop);
174 for (
int i = 0; i < num_sats; ++i) {
175 all_pdrops.appendRow(pdrops.begin(), pdrops.end());
178 int flow_direction = param.getDefault(
"flow_direction", 0);
179 bool start_from_cl = param.getDefault(
"start_from_cl",
true);
185 std::string kr_filename = param.getDefault<std::string>(
"kr_filename",
"upscaled_relperm");
186 std::string krwtmpname = kr_filename +
"_water";
187 std::string krotmpname = kr_filename +
"_oil";
188 std::string krw_filename = param.getDefault<std::string>(
"outputWater", krwtmpname);
189 std::string kro_filename = param.getDefault<std::string>(
"outputOil", krotmpname);
190 std::ofstream krw_out(krw_filename.c_str());
191 std::ofstream kro_out(kro_filename.c_str());
194 krw_out <<
"# Result from steady state upscaling" << std::endl;
195 krw_out <<
"# Pressuredrop Sw Krxx Kryy Krzz" << std::endl;
196 kro_out <<
"# Result from steady state upscaling" << std::endl;
197 kro_out <<
"# Pressuredrop Sw Krxx Kryy Krzz" << std::endl;
198 krw_out.precision(4); krw_out.setf(std::ios::scientific | std::ios::showpoint);
199 kro_out.precision(4); kro_out.setf(std::ios::scientific | std::ios::showpoint);
203 upscaler.init(param);
206 typedef typename Upscaler::permtensor_t permtensor_t;
207 permtensor_t upscaled_K = upscaler.upscaleSinglePhase();
208 permtensor_t upscaled_K_copy = upscaled_K;
209 upscaled_K_copy *= (1.0/(Opm::prefix::milli*Opm::unit::darcy));
210 std::cout.precision(15);
211 std::cout <<
"Upscaled K in millidarcy:\n" << upscaled_K_copy << std::endl;
212 std::cout <<
"Upscaled porosity: " << upscaler.upscalePorosity() << std::endl;
215 int num_cells = upscaler.grid().size(0);
216 int num_sats = saturations.size();
218 for (
int i = 0; i < num_sats; ++i) {
220 std::vector<double> init_sat;
223 upscaler.setToCapillaryLimit(saturations[i], init_sat);
225 init_sat.resize(num_cells, saturations[i]);
226 std::cout <<
"Failed to initialize with capillary limit for s = " << saturations[i]
227 <<
". Init with uniform distribution." << std::endl;
230 init_sat.resize(num_cells, saturations[i]);
233 const Opm::SparseTable<double>::row_type pdrops = all_pdrops[i];
234 int num_pdrops = pdrops.size();
235 for (
int j = 0;
j < num_pdrops; ++
j) {
236 upscaler.initSatLimits(init_sat);
237 double pdrop = pdrops[
j];
238 bool success =
false;
239 std::pair<permtensor_t, permtensor_t> lambda
240 = upscaler.upscaleSteadyState(flow_direction, init_sat, saturations[i], pdrop, upscaled_K, success);
241 double usat = upscaler.lastSaturationUpscaled();
243 std::cout <<
"Upscaling failed for " << usat << std::endl;
245 init_sat = upscaler.lastSaturationState();
247 std::cout <<
"\n\nTensor of upscaled relperms for initial saturation " << saturations[i]
248 <<
", real steady-state saturation " << usat
249 <<
" and pressure drop " << pdrop
250 <<
":\n\n[water]\n" << lambda.first
251 <<
"\n[oil]\n" << lambda.second << std::endl;
253 writeRelPerm(krw_out, lambda.first , usat, pdrop, success);
254 writeRelPerm(kro_out, lambda.second, usat, pdrop, success);
284 #endif // OPM_STEADYSTATEUPSCALERMANAGER_HEADER
Definition: SteadyStateUpscalerManagerImplicit.hpp:137
Definition: applier.hpp:18
void writeControl(Ostream &os, const std::vector< double > &saturations, const Opm::SparseTable< double > &all_pdrops)
Writes saturation and pressure drop data to an output stream.
Definition: SteadyStateUpscalerManager.hpp:81
void upscale(const Opm::parameter::ParameterGroup ¶m)
Definition: SteadyStateUpscalerManagerImplicit.hpp:140
void readControl(Istream &is, std::vector< double > &saturations, Opm::SparseTable< double > &all_pdrops)
Reads saturation and pressure drop data from an input stream.
Definition: SteadyStateUpscalerManager.hpp:55
void writeRelPerm(Ostream &os, const Tensor &K, double sat, double pdrop)
Definition: SteadyStateUpscalerManager.hpp:100
int j
Definition: elasticity_upscale_impl.hpp:658