From 7ad5a4d29b709cfe131b01993f829649330a6c13 Mon Sep 17 00:00:00 2001 From: Tor Harald Sandve Date: Mon, 15 Oct 2018 09:44:42 +0200 Subject: [PATCH] Support modifying trans from deck --- ebos/ecltransmissibility.hh | 78 +++++++++++++++++++++++++++++++++++++ 1 file changed, 78 insertions(+) diff --git a/ebos/ecltransmissibility.hh b/ebos/ecltransmissibility.hh index 9b11f61ca..862f22efa 100644 --- a/ebos/ecltransmissibility.hh +++ b/ebos/ecltransmissibility.hh @@ -348,6 +348,9 @@ public: trans_[isId_(elemIdx, outsideElemIdx)] = trans; } } + + // potentially overwrite and/or modify transmissibilities based on input from deck + updateFromEclState_(); } /*! @@ -389,6 +392,81 @@ public: { return thermalHalfTransBoundary_.at(std::make_pair(insideElemIdx, boundaryFaceIdx)); } private: + + void updateFromEclState_(){ + const auto& gridView = vanguard_.gridView(); + const auto& cartMapper = vanguard_.cartesianIndexMapper(); + const auto& cartDims = cartMapper.cartesianDimensions(); + const auto& properties = vanguard_.eclState().get3DProperties(); + +#if DUNE_VERSION_NEWER(DUNE_GRID, 2,6) + ElementMapper elemMapper(gridView, Dune::mcmgElementLayout()); +#else + ElementMapper elemMapper(gridView); +#endif + + const auto& inputTranx = properties.getDoubleGridProperty("TRANX"); + const auto& inputTrany = properties.getDoubleGridProperty("TRANY"); + const auto& inputTranz = properties.getDoubleGridProperty("TRANZ"); + + // compute the transmissibilities for all intersections + auto elemIt = gridView.template begin(); + const auto& elemEndIt = gridView.template end(); + + for (; elemIt != elemEndIt; ++elemIt) { + const auto& elem = *elemIt; + auto isIt = gridView.ibegin(elem); + const auto& isEndIt = gridView.iend(elem); + for (; isIt != isEndIt; ++ isIt) { + // store intersection, this might be costly + const auto& intersection = *isIt; + if (!intersection.neighbor()) + continue; // intersection is on the domain boundary + + unsigned c1 = elemMapper.index(intersection.inside()); + unsigned c2 = elemMapper.index(intersection.outside()); + + if (c1 > c2) + continue; // we only need to handle each connection once, thank you. + + auto isId = isId_(c1, c2); + + int gc1 = std::min(cartMapper.cartesianIndex(c1), cartMapper.cartesianIndex(c2)); + int gc2 = std::max(cartMapper.cartesianIndex(c1), cartMapper.cartesianIndex(c2)); + + if (gc2 - gc1 == 1) { + if (inputTranx.deckAssigned()) { + // set simulator internal transmissibilities to values from inputTranx + trans_[isId] = inputTranx.iget(gc1); + } else { + // Scale transmissibilities with scale factor from inputTranx + trans_[isId] *= inputTranx.iget(gc1); + } + } + else if (gc2 - gc1 == cartDims[0]) { + if (inputTrany.deckAssigned()) { + // set simulator internal transmissibilities to values from inputTrany + trans_[isId] = inputTrany.iget(gc1); + } else { + // Scale transmissibilities with scale factor from inputTrany + trans_[isId] *= inputTrany.iget(gc1); + } + } + else if (gc2 - gc1 == cartDims[0]*cartDims[1]) { + if (inputTranz.deckAssigned()) { + // set simulator internal transmissibilities to values from inputTranz + trans_[isId] = inputTranz.iget(gc1); + } else { + // Scale transmissibilities with scale factor from inputTranz + trans_[isId] *= inputTranz.iget(gc1); + } + } + //else.. We don't support modification of NNC at the moment. + } + } + } + + template void computeFaceProperties( const Intersection& intersection, const int insideElemIdx,