From f01c4dac10bd5659e7e3c45ee14c7049930d9ba2 Mon Sep 17 00:00:00 2001 From: Liu Ming Date: Mon, 27 Jan 2014 16:53:45 +0800 Subject: [PATCH] fix eigen bug when compute the capPress. --- .../fullyimplicit/IncompPropsAdFromDeck.cpp | 50 +++++++++---------- 1 file changed, 25 insertions(+), 25 deletions(-) diff --git a/opm/polymer/fullyimplicit/IncompPropsAdFromDeck.cpp b/opm/polymer/fullyimplicit/IncompPropsAdFromDeck.cpp index 4ca86247b..45a946c8b 100644 --- a/opm/polymer/fullyimplicit/IncompPropsAdFromDeck.cpp +++ b/opm/polymer/fullyimplicit/IncompPropsAdFromDeck.cpp @@ -137,42 +137,42 @@ namespace Opm const ADB& so, const Cells& cells) const { - - const int n = cells.size(); - const int np = numPhases(); - const int num_blocks = so.numBlocks(); - Block s_all(n, np); - assert(sw.size() == n && so.size() == n); + const int numCells = cells.size(); + const int numActivePhases = numPhases(); + const int numBlocks = so.numBlocks(); + assert(sw.value().size() == numCells); + assert(so.value().size() == numCells); + Block s_all(numCells, numActivePhases); s_all.col(0) = sw.value(); s_all.col(1) = so.value(); - Block pc(n, np); - Block dpc(n, np * np); + Block pc(numCells, numActivePhases); + Block dpc(numCells, numActivePhases*numActivePhases); + satprops_.capPress(numCells, s_all.data(), cells.data(), pc.data(), dpc.data()); - satprops_.capPress(n, s_all.data(), cells.data(), pc.data(), dpc.data()); - - std::vector capPressures; - capPressures.reserve(2); + std::vector adbCapPressures; + adbCapPressures.reserve(2); const ADB* s[2] = { &sw, &so}; - for (int phase1 = 0; phase1 < 3; ++phase1) { - const int phase1_pos = phase1; - std::vector jacs(num_blocks); - for (int block = 0; block < num_blocks; ++block) { - jacs[block] = ADB::M(n, s[phase1]->derivative()[block].cols()); + for (int phase1 = 0; phase1 < 2; ++phase1) { + const int phase1_pos = phase1; + std::vector jacs(numBlocks); + for (int block = 0; block < numBlocks; ++block) { + jacs[block] = ADB::M(numCells, s[phase1]->derivative()[block].cols()); } - for (int phase2 = 0; phase2 < 3; ++phase2) { + for (int phase2 = 0; phase2 < 2; ++phase2) { const int phase2_pos = phase2; - // Assemble dpc1/ds2. - const int column = phase1_pos + phase2_pos; // Recall: Fortran ordering from props_.relperm() + // Assemble dpc1/ds2. + const int column = phase1_pos + numActivePhases*phase2_pos; // Recall: Fortran ordering from props_.relperm() ADB::M dpc1_ds2_diag = spdiag(dpc.col(column)); - for (int block = 0; block < num_blocks; ++block) { + for (int block = 0; block < numBlocks; ++block) { jacs[block] += dpc1_ds2_diag * s[phase2]->derivative()[block]; } - } - capPressures.emplace_back(ADB::function(pc.col(phase1_pos), jacs)); - } - return capPressures; + } + adbCapPressures.emplace_back(ADB::function(pc.col(phase1_pos), jacs)); + } + return adbCapPressures; } + } //namespace Opm