Use volume weighted arithmetic average of NTG for cells merged by minpv.

This is a big hack to make the TRAN[X,Y] like flow_legacy and a bit more
like Eclipse.
This commit is contained in:
Tor Harald Sandve
2017-05-12 12:39:03 +02:00
parent bcdc4e5e38
commit 5fb005d87b

View File

@@ -108,6 +108,10 @@ public:
const std::vector<double>& ntg =
eclState.get3DProperties().getDoubleGridProperty("NTG").getData();
// modify the ntg values for the cells merged with minpv
std::vector<double> ntg_mod (ntg);
minPvFillNtg_(ntg_mod);
unsigned numElements = elemMapper.size();
extractPermeability_();
@@ -209,8 +213,8 @@ public:
axisCentroids),
permeability_[outsideElemIdx]);
applyNtg_(halfTrans1, insideFaceIdx, insideCartElemIdx, ntg);
applyNtg_(halfTrans2, outsideFaceIdx, outsideCartElemIdx, ntg);
applyNtg_(halfTrans1, insideFaceIdx, insideCartElemIdx, ntg_mod);
applyNtg_(halfTrans2, outsideFaceIdx, outsideCartElemIdx, ntg_mod);
// convert half transmissibilities to full face
// transmissibilities using the harmonic mean
@@ -398,6 +402,44 @@ private:
}
}
void minPvFillNtg_(std::vector<double>& ntg) const {
// compute volume weighted arithmetic average of NTG for
// cells merged as an result of minpv.
const auto& eclState = gridManager_.eclState();
const auto& eclGrid = eclState.getInputGrid();
const auto& porv = eclState.get3DProperties().getDoubleGridProperty("PORV").getData();
const auto& actnum = eclState.get3DProperties().getIntGridProperty("ACTNUM").getData();
const auto& cartMapper = gridManager_.cartesianIndexMapper();
const auto& cartDims = cartMapper.cartesianDimensions();
assert(dimWorld > 1);
const size_t nxny = cartDims[0] * cartDims[1];
for (size_t cartesianCellIdx = 0; cartesianCellIdx < ntg.size(); ++cartesianCellIdx) {
const double cellVolume = eclGrid.getCellVolume(cartesianCellIdx);
ntg[cartesianCellIdx] *= cellVolume;
double totalCellVolume = cellVolume;
// Average properties as long as there exist cells above
// that has pore volume less than the MINPV threshold
int cartesianCellIdxAbove = cartesianCellIdx - nxny;
while ( cartesianCellIdxAbove >= 0 &&
actnum[cartesianCellIdxAbove] > 0 &&
porv[cartesianCellIdxAbove] < eclGrid.getMinpvValue() ) {
// Volume weighted arithmetic average of NTG
const double cellAboveVolume = eclGrid.getCellVolume(cartesianCellIdxAbove);
totalCellVolume += cellAboveVolume;
ntg[cartesianCellIdx] += ntg[cartesianCellIdxAbove]*cellAboveVolume;
cartesianCellIdxAbove -= nxny;
}
ntg[cartesianCellIdx] /= totalCellVolume;
}
}
const GridManager& gridManager_;
std::vector<DimMatrix> permeability_;
std::unordered_map<std::uint64_t, Scalar> trans_;