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

Added lots of (inactive) experimental code.

Browse Source 
The code attempts to improve #iterations in the Gauss-Seidel-like
multicell solver by improving ordering. In general, experiment failed
to improve #iterations, except for one: totally random order was the
best (for the 100x100 case tried)!
This commit is contained in:
Atgeirr Flø Rasmussen 2012-02-17 09:39:45 +01:00
parent 8a2a3c3ff8
commit ccdacb791c
2 changed files with 82 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

80
opm/core/transport/reorder/TransportModelTwophase.cpp
View File

@ -20,7 +20,7 @@
#include <opm/core/transport/reorder/TransportModelTwophase.hpp>
#include <opm/core/fluid/IncompPropertiesInterface.hpp>
#include <opm/core/grid.h>
#include <opm/core/transport/reorder/nlsolvers.h>
#include <opm/core/transport/reorder/reordersequence.h>
#include <opm/core/utility/RootFinders.hpp>
#include <fstream>
@ -41,6 +41,8 @@ namespace Opm
dt_(0.0),
saturation_(0),
fractionalflow_(grid.number_of_cells, -1.0)
// ia_(grid.number_of_cells + 1, -1),
// ja_(grid.number_of_faces, -1)
{
if (props.numPhases() != 2) {
THROW("Property object must have 2 phases");
@ -57,6 +59,15 @@ namespace Opm
source_ = source;
dt_ = dt;
saturation_ = saturation;
// std::vector<int> seq(grid_.number_of_cells);
// std::vector<int> comp(grid_.number_of_cells + 1);
// int ncomp;
// compute_sequence_graph(&grid_, darcyflux_,
// &seq[0], &comp[0], &ncomp,
// &ia_[0], &ja_[0]);
reorderAndTransport(grid_, darcyflux);
}
@ -135,6 +146,71 @@ namespace Opm
{
// std::ofstream os("dump");
// std::copy(cells, cells + num_cells, std::ostream_iterator<double>(os, "\n"));
// Experiment: try a breath-first search to build a more suitable ordering.
// Verdict: failed to improve #iterations.
// {
// std::vector<int> pos(grid_.number_of_cells, -1);
// for (int i = 0; i < num_cells; ++i) {
// const int cell = cells[i];
// pos[cell] = i;
// }
// std::vector<int> done_pos(num_cells, 0);
// std::vector<int> upstream_pos;
// std::vector<int> new_pos;
// upstream_pos.push_back(0);
// done_pos[0] = 1;
// int current = 0;
// while (new_pos.size() < num_cells) {
// const int i = upstream_pos[current++];
// new_pos.push_back(i);
// const int cell = cells[i];
// for (int j = ia_[cell]; j < ia_[cell+1]; ++j) {
// const int opos = pos[ja_[j]];
// if (!done_pos[opos]) {
// upstream_pos.push_back(opos);
// done_pos[opos] = 1;
// }
// }
// }
// std::reverse(new_pos.begin(), new_pos.end());
// std::copy(new_pos.begin(), new_pos.end(), const_cast<int*>(cells));
// }
// Experiment: try a random ordering.
// Verdict: amazingly, reduced #iterations by approx. 25%!
// int* c = const_cast<int*>(cells);
// std::random_shuffle(c, c + num_cells);
// Experiment: implement a metric measuring badness of ordering
// as average distance in (cyclic) ordering from
// upstream neighbours.
// Verdict: does not seem to predict #iterations very well, if at all.
// std::vector<int> pos(grid_.number_of_cells, -1);
// for (int i = 0; i < num_cells; ++i) {
// const int cell = cells[i];
// pos[cell] = i;
// }
// double diffsum = 0;
// for (int i = 0; i < num_cells; ++i) {
// const int cell = cells[i];
// int num_upstream = 0;
// int loc_diffsum = 0;
// for (int j = ia_[cell]; j < ia_[cell+1]; ++j) {
// const int opos = pos[ja_[j]];
// if (opos == -1) {
// std::cout << "Hmmm" << std::endl;
// continue;
// }
// ++num_upstream;
// const int diff = (i - opos + num_cells) % num_cells;
// loc_diffsum += diff;
// }
// diffsum += double(loc_diffsum)/double(num_upstream);
// }
// std::cout << "Average distance from upstream neighbours: " << diffsum/double(num_cells)
// << std::endl;
double max_s_change = 0.0;
const double tol = 1e-9;
const int max_iters = 300;
@ -155,7 +231,7 @@ namespace Opm
const double old_s = saturation_[cell];
saturation_[cell] = s0[i];
solveSingleCell(cell);
// std::cout << "delta s = " << saturation_[cell] - old_s << std::endl;
// std::cout << "cell = " << cell << " delta s = " << saturation_[cell] - old_s << std::endl;
if (max_s_change < std::fabs(saturation_[cell] - old_s)) {
max_change_cell = cell;
}

4
opm/core/transport/reorder/TransportModelTwophase.hpp
View File

@ -57,6 +57,10 @@ namespace Opm
std::vector<double> fractionalflow_; // one per cell
const double* visc_;
// Storing the upwind graph for experiments.
// std::vector<int> ia_;
// std::vector<int> ja_;
struct Residual;
double fracFlow(double s, int cell) const;
};