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

Added Newton iterations to column solver.

Browse Source
This commit is contained in:
Atgeirr Flø Rasmussen 2012-03-02 13:45:33 +01:00
parent 560cb08313
commit dc80e2a38f
1 changed files with 21 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

25
opm/core/transport/GravityColumnSolver_impl.hpp
View File

@ -78,11 +78,28 @@ namespace Opm
StateWithZeroFlux state(s); // This holds s by reference.
JacSys sys(grid_.number_of_cells);
model_.initStep(state, grid_, sys);
model_.initIteration(state, grid_, sys);
std::map<int, std::vector<int> >::const_iterator it;
for (it = columns.begin(); it != columns.end(); ++it) {
solveSingleColumn(it->second, dt, s, sys.vector().writableSolution());
const int max_iter = 40;
const double tol = 1e-4;
int iter = 0;
double max_delta = 1e100;
while (iter < max_iter) {
model_.initIteration(state, grid_, sys);
std::map<int, std::vector<int> >::const_iterator it;
for (it = columns.begin(); it != columns.end(); ++it) {
solveSingleColumn(it->second, dt, s, sys.vector().writableSolution());
}
const double maxelem = *std::max_element(sys.vector().solution().begin(), sys.vector().solution().end());
const double minelem = *std::min_element(sys.vector().solution().begin(), sys.vector().solution().end());
max_delta = std::max(maxelem, -minelem);
std::cout << "Iteration " << iter << " max_delta = " << max_delta << std::endl;
if (max_delta < tol) {
break;
}
++iter;
}
if (max_delta >= tol) {
THROW("Failed to converge!");
}
// Finalize.
// model_.finishIteration(); // Doesn't do anything in th 2p model.