OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-01-17 19:33:00 -06:00
Code Issues Packages Projects Releases Wiki Activity

Revising the function stablizeNewton().

Browse Source 
For the SOR type relaxation, we need the dx from the previous iteration.
This commit is contained in:
Kai Bao 2014-05-20 11:25:04 +02:00
parent 5b6c325a32
commit 8dcae05a7b
2 changed files with 25 additions and 5 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
opm/autodiff/FullyImplicitBlackoilSolver.hpp
View File

@ -276,7 +276,7 @@ namespace Opm {
const int it, const double relaxRelTol,
bool &oscillate, bool &stagnate ) const;
void stablizeNewton(V &dx, const bool &oscillate, const bool &stagnate, const double omega,
void stablizeNewton(V &dx, V &dxOld, const bool &oscillate, const bool &stagnate, const double omega,
const RelaxType relax_type) const;
};

28
opm/autodiff/FullyImplicitBlackoilSolver_impl.hpp
View File

@ -278,8 +278,25 @@ namespace {
<< std::setw(9) << it << std::setprecision(9)
<< std::setw(18) << r0 << std::endl;
// compute the size of the linear system
int linearSize = 0;
std::vector<ADB>::const_iterator massBalanceIt = residual_.material_balance_eq.begin();
const std::vector<ADB>::const_iterator endMassBalanceIt = residual_.material_balance_eq.end();
for (; massBalanceIt != endMassBalanceIt; ++massBalanceIt) {
linearSize += (*massBalanceIt).size();
}
linearSize += residual_.well_flux_eq.size();
linearSize += residual_.well_eq.size();
std::cout << " the size of the linear system is " << linearSize << std::endl;
std::cin.ignore();
V dxOld = V::Zero(linearSize);
while ((!converged) && (it < maxit)) {
const V dx = solveJacobianSystem();
V dx = solveJacobianSystem();
updateState(dx, x, xw);
@ -1707,9 +1724,12 @@ namespace {
template<class T>
void
FullyImplicitBlackoilSolver<T>::stablizeNewton( V &dx, const bool &oscillate, const bool &stagnate,
FullyImplicitBlackoilSolver<T>::stablizeNewton( V &dx, V &dxOld, const bool &oscillate, const bool &stagnate,
const double omega, const RelaxType relax_type) const {
const V tempDxOld = dxOld;
dxOld = dx;
switch (relax_type) {
case DAMPEN:
if (omega == 1.) {
@ -1721,8 +1741,8 @@ namespace {
if (omega == 1.) {
return;
}
const V dxold = dx;
dx = dx*omega + (1.-omega)*dxold;
// const V dxold = dx;
dx = dx*omega + (1.-omega)*tempDxOld;
return;
default:
OPM_THROW(std::runtime_error, "Can only handle DAMPEN and SOR relaxation type.");