OilfieldToolsNet/opm-core
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

merge.

Browse Source
This commit is contained in:
Kjetil Olsen Lye 2012-04-12 15:48:59 +02:00
commit 34460f1711
2 changed files with 36 additions and 29 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

28
examples/spu_2p.cpp
View File

@ -94,6 +94,7 @@
#include <fstream> #include <fstream>
#include <iterator> #include <iterator>
#include <vector> #include <vector>
#include <numeric>
@ -272,6 +273,7 @@ main(int argc, char** argv)
create_directories(fpath); create_directories(fpath);
output_interval = param.getDefault("output_interval", output_interval); output_interval = param.getDefault("output_interval", output_interval);
} }
const int num_transport_substeps = param.getDefault("num_transport_substeps", 1);
// If we have a "deck_filename", grid and props will be read from that. // If we have a "deck_filename", grid and props will be read from that.
bool use_deck = param.has("deck_filename"); bool use_deck = param.has("deck_filename");
@ -358,11 +360,6 @@ main(int argc, char** argv)
if (!use_reorder) { if (!use_reorder) {
THROW("Cannot run implicit (non-reordering) transport solver with rock compressibility yet."); THROW("Cannot run implicit (non-reordering) transport solver with rock compressibility yet.");
} }
if (use_segregation_split) {
if (!use_gauss_seidel_gravity) {
THROW("For gravity segregation splitting, only use_gauss_seidel_gravity=true supports rock compressibility.");
}
}
nl_pressure_maxiter = param.getDefault("nl_pressure_maxiter", 10); nl_pressure_maxiter = param.getDefault("nl_pressure_maxiter", 10);
nl_pressure_tolerance = param.getDefault("nl_pressure_tolerance", 1.0); // in Pascal nl_pressure_tolerance = param.getDefault("nl_pressure_tolerance", 1.0); // in Pascal
} }
@ -560,32 +557,39 @@ main(int argc, char** argv)
// Solve transport. // Solve transport.
transport_timer.start(); transport_timer.start();
double stepsize = simtimer.currentStepLength();
if (num_transport_substeps != 1) {
stepsize /= double(num_transport_substeps);
std::cout << "Making " << num_transport_substeps << " transport substeps." << std::endl;
}
for (int tr_substep = 0; tr_substep < num_transport_substeps; ++tr_substep) {
if (use_reorder) { if (use_reorder) {
Opm::toWaterSat(state.saturation(), reorder_sat); Opm::toWaterSat(state.saturation(), reorder_sat);
reorder_model.solve(&state.faceflux()[0], &porevol[0], &reorder_src[0], reorder_model.solve(&state.faceflux()[0], &porevol[0], &reorder_src[0],
simtimer.currentStepLength(), &reorder_sat[0]); stepsize, &reorder_sat[0]);
Opm::toBothSat(reorder_sat, state.saturation()); Opm::toBothSat(reorder_sat, state.saturation());
Opm::computeInjectedProduced(*props, state.saturation(), src, simtimer.currentStepLength(), injected, produced); Opm::computeInjectedProduced(*props, state.saturation(), src, stepsize, injected, produced);
if (use_segregation_split) { if (use_segregation_split) {
if (use_column_solver) { if (use_column_solver) {
if (use_gauss_seidel_gravity) { if (use_gauss_seidel_gravity) {
reorder_model.solveGravity(columns, &porevol[0], simtimer.currentStepLength(), reorder_sat); reorder_model.solveGravity(columns, &porevol[0], stepsize, reorder_sat);
Opm::toBothSat(reorder_sat, state.saturation()); Opm::toBothSat(reorder_sat, state.saturation());
} else { } else {
colsolver.solve(columns, simtimer.currentStepLength(), state.saturation()); colsolver.solve(columns, stepsize, state.saturation());
} }
} else { } else {
std::vector<double> fluxes = state.faceflux(); std::vector<double> fluxes = state.faceflux();
std::fill(state.faceflux().begin(), state.faceflux().end(), 0.0); std::fill(state.faceflux().begin(), state.faceflux().end(), 0.0);
tsolver.solve(*grid->c_grid(), tsrc, simtimer.currentStepLength(), ctrl, state, linsolve, rpt); tsolver.solve(*grid->c_grid(), tsrc, stepsize, ctrl, state, linsolve, rpt);
std::cout << rpt; std::cout << rpt;
state.faceflux() = fluxes; state.faceflux() = fluxes;
} }
} }
} else { } else {
tsolver.solve(*grid->c_grid(), tsrc, simtimer.currentStepLength(), ctrl, state, linsolve, rpt); tsolver.solve(*grid->c_grid(), tsrc, stepsize, ctrl, state, linsolve, rpt);
std::cout << rpt; std::cout << rpt;
Opm::computeInjectedProduced(*props, state.saturation(), src, simtimer.currentStepLength(), injected, produced); Opm::computeInjectedProduced(*props, state.saturation(), src, stepsize, injected, produced);
}
} }
transport_timer.stop(); transport_timer.stop();
double tt = transport_timer.secsSinceStart(); double tt = transport_timer.secsSinceStart();

7
opm/core/fluid/SaturationPropsFromDeck.cpp
View File

@ -43,6 +43,7 @@ namespace Opm
} }
const int samples = 200; const int samples = 200;
double swco = 0.0; double swco = 0.0;
double swmax = 1.0;
if (phase_usage_.phase_used[Aqua]) { if (phase_usage_.phase_used[Aqua]) {
const SWOF::table_t& swof_table = deck.getSWOF().swof_; const SWOF::table_t& swof_table = deck.getSWOF().swof_;
if (swof_table.size() != 1) { if (swof_table.size() != 1) {
@ -58,7 +59,8 @@ namespace Opm
krocw_ = krow[0]; // At connate water -> ecl. SWOF krocw_ = krow[0]; // At connate water -> ecl. SWOF
swco = sw[0]; swco = sw[0];
smin_[phase_usage_.phase_pos[Aqua]] = sw[0]; smin_[phase_usage_.phase_pos[Aqua]] = sw[0];
smax_[phase_usage_.phase_pos[Aqua]] = 1.0; swmax = sw.back();
smax_[phase_usage_.phase_pos[Aqua]] = sw.back();
} }
if (phase_usage_.phase_used[Vapour]) { if (phase_usage_.phase_used[Vapour]) {
const SGOF::table_t& sgof_table = deck.getSGOF().sgof_; const SGOF::table_t& sgof_table = deck.getSGOF().sgof_;
@ -79,7 +81,8 @@ namespace Opm
} }
smax_[phase_usage_.phase_pos[Vapour]] = sg.back(); smax_[phase_usage_.phase_pos[Vapour]] = sg.back();
} }
smin_[phase_usage_.phase_pos[Liquid]] = 0.0; // These only consider water min/max sats. Consider gas sats?
smin_[phase_usage_.phase_pos[Liquid]] = 1.0 - swmax;
smax_[phase_usage_.phase_pos[Liquid]] = 1.0 - swco; smax_[phase_usage_.phase_pos[Liquid]] = 1.0 - swco;
} }