36 #ifndef OPM_STEADYSTATEUPSCALERMANAGER_HEADER
37 #define OPM_STEADYSTATEUPSCALERMANAGER_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)
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';
136 template <
class Traits>
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);
185 typedef typename Upscaler::permtensor_t permtensor_t;
187 upscaler.
init(param);
190 permtensor_t upscaled_K = upscaler.upscaleSinglePhase();
191 permtensor_t upscaled_K_copy = upscaled_K;
192 upscaled_K_copy *= (1.0/(Opm::prefix::milli*Opm::unit::darcy));
193 std::cout.precision(15);
194 std::cout <<
"Upscaled K in millidarcy:\n" << upscaled_K_copy << std::endl;
195 std::cout <<
"Upscaled porosity: " << upscaler.upscalePorosity() << std::endl;
198 std::string kr_filename = param.getDefault<std::string>(
"kr_filename",
"upscaled_relperm");
199 std::string krw_filename = kr_filename +
"_water";
200 std::string kro_filename = kr_filename +
"_oil";
201 std::ofstream krw_out(krw_filename.c_str());
202 std::ofstream kro_out(kro_filename.c_str());
203 krw_out <<
"# Result from steady state upscaling" << std::endl;
204 krw_out <<
"# Pressuredrop Sw Krxx Kryy Krzz" << std::endl;
205 kro_out <<
"# Result from steady state upscaling" << std::endl;
206 kro_out <<
"# Pressuredrop Sw Krxx Kryy Krzz" << std::endl;
209 krw_out.precision(15); krw_out.setf(std::ios::scientific | std::ios::showpoint);
210 kro_out.precision(15); kro_out.setf(std::ios::scientific | std::ios::showpoint);
214 int num_cells = upscaler.grid().size(0);
215 int num_sats = saturations.size();
216 for (
int i = 0; i < num_sats; ++i) {
218 std::vector<double> init_sat(num_cells, saturations[i]);
219 const Opm::SparseTable<double>::row_type pdrops = all_pdrops[i];
220 int num_pdrops = pdrops.size();
221 for (
int j = 0;
j < num_pdrops; ++
j) {
222 double pdrop = pdrops[
j];
223 std::pair<permtensor_t, permtensor_t> lambda
224 = upscaler.upscaleSteadyState(flow_direction, init_sat, saturations[i], pdrop, upscaled_K);
225 double usat = upscaler.lastSaturationUpscaled();
226 std::cout <<
"\n\nTensor of upscaled relperms for initial saturation " << saturations[i]
227 <<
", real steady-state saturation " << usat
228 <<
" and pressure drop " << pdrop
229 <<
":\n\n[water]\n" << lambda.first
230 <<
"\n[oil]\n" << lambda.second << std::endl;
232 init_sat = upscaler.lastSaturationState();
250 #endif // OPM_STEADYSTATEUPSCALERMANAGER_HEADER
void upscale(const Opm::parameter::ParameterGroup ¶m)
Definition: SteadyStateUpscalerManager.hpp:140
Definition: applier.hpp:18
A class for doing steady state upscaling.
Definition: SteadyStateUpscaler.hpp:51
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
Definition: SteadyStateUpscalerManager.hpp:137
void init(const Opm::parameter::ParameterGroup ¶m)
Initializes the upscaler from parameters.
Definition: UpscalerBase_impl.hpp:64
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