Merge pull request #1035 from michal-toth/partitioning/lower-bigM-constant

Make well edge weight equal to sum of grid edges instead of limits::max

Author: blattms

opm/grid/common/MetisPartition.cpp

@@ -275,10 +275,13 @@ metisSerialGraphPartitionGridOnRoot(const CpGrid& cpgrid,
 
         if( wells )
         {
+            // well edge weight for partitioning, big enough that wells should not get split
+            idx_t weWeight = sumOfGridEdges<idx_t>(cpgrid, *gridAndWells);
+
             int neighborCounter = 0;
             for( int cell = 0; cell < n;  cell++ )
             {
-                fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(*gridAndWells, lids[cell], gids, neighborCounter, adjncy, adjwgt);
+                fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(*gridAndWells, lids[cell], gids, neighborCounter, adjncy, adjwgt, weWeight);
             }
         }
         else

opm/grid/common/ZoltanGraphFunctions.cpp

@@ -179,6 +179,23 @@ void getNullEdgeList(void *cpGridPointer, int sizeGID, int sizeLID,
     *err = ZOLTAN_OK;
 }
 
+template <typename EdgeWeightType>
+EdgeWeightType sumOfGridEdges(const Dune::CpGrid& grid,
+                              const CombinedGridWellGraph& graph)
+{
+    EdgeWeightType weWeight = 0;
+    for (int edge=0; edge<grid.numFaces(); ++edge) {
+        // assumes that transmissibility is positive (too obvious to leave an assert)
+        const auto& addend = graph.transmissibility(edge);
+        if (weWeight >= std::numeric_limits<EdgeWeightType>::max()-addend) {
+            weWeight = std::numeric_limits<EdgeWeightType>::max();
+            break;
+        }
+        weWeight += graph.transmissibility(edge);
+    }
+    return weWeight;
+}
+
 template<typename ID>
 void fillNBORGIDForSpecificCellAndIncrementNeighborCounter(const Dune::CpGrid& grid, int localCellId, ID globalID, int& neighborCounter, ID& nborGID) {
     for ( int local_face = 0 ; local_face < grid.numCellFaces(localCellId); ++local_face )
@@ -245,14 +262,14 @@ void getCpGridEdgeList(void *cpGridPointer, int sizeGID, int sizeLID,
 #endif
 }
 template<typename ID, typename weightType>
-void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(const CombinedGridWellGraph& graph, const int localCellId, ID globalID, int& neighborCounter, ID& nborGID, weightType *ewgts) {
+void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(const CombinedGridWellGraph& graph, const int localCellId, ID globalID, int& neighborCounter, ID& nborGID, weightType *ewgts, const weightType& weWeight) {
     const Dune::CpGrid&  grid = graph.getGrid();
     // First the strong edges of the well completions.
     auto wellEdges = graph.getWellsGraph()[localCellId];
     for( auto edge : wellEdges)
     {
         nborGID[neighborCounter] = edge;
-        ewgts[neighborCounter++] = std::numeric_limits<weightType>::max();
+        ewgts[neighborCounter++] = weWeight;
     }
 
     // Now the ones of the grid that are not handled by the well completions
@@ -308,9 +325,12 @@ void getCpGridWellsEdgeList(void *graphPointer, int sizeGID, int sizeLID,
 #endif
     int neighborCounter = 0;
 
+    // well edge weight for partitioning, big enough that wells should not get split
+    float weWeight = sumOfGridEdges<float>(grid, graph);
+
     for( int cell = 0; cell < numCells;  cell++ )
     {
-        fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(graph, localID[cell], globalID, neighborCounter, nborGID, ewgts);
+        fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(graph, localID[cell], globalID, neighborCounter, nborGID, ewgts, weWeight);
 #ifndef NDEBUG
         assert(neighborCounter-oldNeighborCounter==numEdges[cell]);
         oldNeighborCounter = neighborCounter;
@@ -416,12 +436,16 @@ void setCpGridZoltanGraphFunctions(Zoltan_Struct *zz,
 template
 void fillNBORGIDForSpecificCellAndIncrementNeighborCounter(const Dune::CpGrid&, int, int*, int&, int*& nborGID);
 template
-void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(const CombinedGridWellGraph&, const int, int*, int&, int*&, int*);
+void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(const CombinedGridWellGraph&, const int, int*, int&, int*&, int*, const int&);
+template
+int sumOfGridEdges(const Dune::CpGrid& grid, const CombinedGridWellGraph& graph);
 
 template
 void fillNBORGIDForSpecificCellAndIncrementNeighborCounter(Dune::CpGrid const&, int, long*, int&, long*&);
 template
-void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(Dune::cpgrid::CombinedGridWellGraph const&, int, long*, int&, long*&, long*);
+void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(Dune::cpgrid::CombinedGridWellGraph const&, int, long*, int&, long*&, long*, const long&);
+template
+long sumOfGridEdges(const Dune::CpGrid& grid, const CombinedGridWellGraph& graph);
 
 #endif
 

opm/grid/common/ZoltanGraphFunctions.hpp

@@ -244,9 +244,16 @@ class CombinedGridWellGraph
 /// \brief Get the number of edges of the graph of the grid and the wells for one cell
 int getNumberOfEdgesForSpecificCellForGridWithWells(const CombinedGridWellGraph& graph, int localCellId);
 
+/// \brief Iterate over the grid and get the sum of all edge weights
+///
+/// Used as a weight for edges between cells of a well
+template <typename EdgeWeightType>
+EdgeWeightType sumOfGridEdges(const Dune::CpGrid& grid,
+                              const CombinedGridWellGraph& graph);
+
 /// \brief Get the list of edges and weights for one cell of a grid with wells
 template<typename ID, typename weightType>
-void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(const CombinedGridWellGraph& graph, const int localCellId, ID globalID, int& neighborCounter, ID& nborGID, weightType *ewgts);
+void fillNBORGIDAndWeightsForSpecificCellAndIncrementNeighborCounterForGridWithWells(const CombinedGridWellGraph& graph, const int localCellId, ID globalID, int& neighborCounter, ID& nborGID, weightType *ewgts, const weightType& weWeight);
 
 #ifdef HAVE_ZOLTAN
 /// \brief Sets up the call-back functions for ZOLTAN's graph partitioning.