2014-02-16 3 views
15

Ein klassischer Algorithmus Frage in 2D-Version der Regel alsDas maximale Volumen von Trapped Regen Wasser in 3D

beschrieben

Da n nicht negativen ganzen Zahlen eine Höhenkarte darstellt, wobei die Breite der einzelnen Balken 1 ist, berechnen, wie viel Wasser Es ist in der Lage, nach dem Regen zu fangen.

Zum Beispiel Da der Eingang

[0,1,0,2,1,0,1,3,2,1,2,1] 

der Rückgabewert

6 

enter image description here

Der Algorithmus, den ich verwenden, würde das obige 2D Problem zu lösen, ist

int trapWaterVolume2D(vector<int> A) { 
    int n = A.size(); 
    vector<int> leftmost(n, 0), rightmost(n, 0); 

    //left exclusive scan, O(n), the highest bar to the left each point 
    int leftMaxSoFar = 0; 
    for (int i = 0; i < n; i++){ 
     leftmost[i] = leftMaxSoFar; 
     if (A[i] > leftMaxSoFar) leftMaxSoFar = A[i]; 
    } 


    //right exclusive scan, O(n), the highest bar to the right each point 
    int rightMaxSoFar = 0; 
    for (int i = n - 1; i >= 0; i--){ 
     rightmost[i] = rightMaxSoFar; 
     if (A[i] > rightMaxSoFar) rightMaxSoFar = A[i]; 
    } 

    // Summation, O(n) 
    int vol = 0; 
    for (int i = 0; i < n; i++){ 
     vol += max(0, min(leftmost[i], rightmost[i]) - A[i]); 
    } 
    return vol; 
} 

Meine Frage ist, wie man den obigen Algorithmus auf die 3D-Version des Problems erweiterbar macht, um das Maximum an Wasser zu berechnen, das im realen 3D-Gelände gefangen ist.

enter image description here

Wir kennen die Höhe an jedem (x, y) Punkt und das Ziel ist die maximale Wassermenge zu berechnen, die in gefangen werden kann: das heißt, um

int trapWaterVolume3D(vector<vector<int> > A); 

Beispiel Graph zu implementieren die Form. Alle Gedanken und Referenzen sind willkommen.

+0

Es ist derselbe Algorithmus, der zur Berechnung der Fraktur in Metallen verwendet wurde. –

+1

Diese Frage ist wohl auch im Zusammenhang mit: http://stackoverflow.com/questions/15033555/tips-on-finding-the-volume-of-water-in-a-3d-chess-board/ – user3838351

Antwort

12

Für jeden Punkt auf dem Gelände alle Pfade von diesem Punkt bis zum Rand des Geländes berücksichtigen. Das Wasserniveau wäre das Minimum der maximalen Höhen der Punkte dieser Wege. Um es zu finden, müssen wir einen leicht modifizierten Dijkstra-Algorithmus ausführen, der die Wasserpegelmatrix von der Grenze aus füllt.

3

Der "leicht modifizierte Dijkstra-Algorithmus" von user3290797 ist näher am Algorithmus von Prim als Dijkstra. In minimalen Spannbaumbegriffen erstellen wir einen Graphen mit einem Eckpunkt pro Kachel, einem Eckpunkt für die Außenseite und Kanten mit Gewichten, die gleich der maximalen Höhe ihrer zwei angrenzenden Kacheln sind (die Außenseite hat die Höhe "minus unendlich").

Wenn in diesem Diagramm ein Pfad zum äußeren Scheitelpunkt angegeben wird, ist das maximale Gewicht einer Kante im Pfad die Höhe, die das Wasser erreichen muss, um entlang diesem Pfad zu entkommen. Die relevante Eigenschaft eines minimalen aufspannenden Baums ist, dass für jedes Paar von Vertices das maximale Gewicht einer Kante in dem Pfad in dem aufspannenden Baum das Minimum ist, das unter allen Wegen zwischen diesen Vertices möglich ist. Der minimale Spannbaum beschreibt somit die wirtschaftlichsten Fluchtwege für Wasser, und die Wasserhöhen können in linearer Zeit mit einer Traversierung extrahiert werden.

Da der Graph eben ist, gibt es einen linearen Algorithmus zur Berechnung des minimalen Spannbaums, bestehend aus alternierenden Boruvka-Pässen und Vereinfachungen. Dies verbessert die O (n log n) Laufzeit von Prim.

+0

Nun, im Dijkstra-Algorithmus wird der nächste Eckpunkt basierend auf den Kosten des Pfades zu diesem Eckpunkt ausgewählt (genau wie in meinem), und im Prim-Algorithmus wird die nächste Kante basierend auf ihren eigenen Kosten ausgewählt; Deshalb glaube ich, dass mein Algorithmus Dijkstra eher ähnelt als Prim. – Anton

2

Dieses Problem ist sehr nahe an der Konstruktion der morphologischen Wasserscheide eines Graustufenbildes (http://en.wikipedia.org/wiki/Watershed_(image_processing)).

Ein Ansatz ist wie folgt (Flutverfahren):

  • sortieren alle Pixel durch Erhöhung erhöht.

  • Arbeit schrittweise, durch Erhebungen zunimmt, Etiketten auf die Pixel pro Einzugsgebiet zuordnen.

  • für eine neue Erhebung Ebene, müssen Sie einen neuen Satz von Pixeln beschriften:

    • Einige Nachbarn nicht markiert haben, sie eine lokale Minimalkonfiguration bilden und ein neues Einzugsgebiet beginnen. Einige
    • haben nur Nachbarn mit dem gleichen Label, können sie in ähnlicher Weise gekennzeichnet werden (sie erstrecken sich ein Einzugsgebiet). Einige
    • haben Nachbarn mit unterschiedlichen Etiketten. Sie gehören nicht zu einem bestimmten Einzugsgebiet und definieren die Wasserscheidegrenzen.

Sie müssen den Standard-Wendepunkt-Algorithmus zu verbessern, um der Lage sein, das Wasservolumen zu berechnen. Sie können dies tun, indem Sie den maximalen Wasserstand in jedem Becken bestimmen und die Bodenhöhe auf jedem Pixel ableiten. Der Wasserspiegel in einem Becken wird durch die Höhe des niedrigsten Wasserstrichpixels um ihn herum bestimmt.

Sie können jedes Mal, handeln Sie einen Wendepunkt Pixel entdecken: wenn ein benachbartes Becken noch eine Ebene nicht zugeordnet wurde, kann das Becken das aktuelle Niveau stehen, ohne zu lecken.

1

Um in 3D Klopfen Wasserproblem zu erreichen heißt das maximale Volumen der eingeschlossenen regen Wasser zu berechnen Sie etwas tun können:

#include<bits/stdc++.h> 
using namespace std; 

#define MAX 10 

int new2d[MAX][MAX]; 
int dp[MAX][MAX],visited[MAX][MAX]; 


int dx[] = {1,0,-1,0}; 
int dy[] = {0,-1,0,1}; 


int boundedBy(int i,int j,int k,int in11,int in22) 
{ 
    if(i<0 || j<0 || i>=in11 || j>=in22) 
     return 0; 

    if(new2d[i][j]>k) 
     return new2d[i][j]; 

    if(visited[i][j]) return INT_MAX; 

    visited[i][j] = 1; 

    int r = INT_MAX; 
    for(int dir = 0 ; dir<4 ; dir++) 
    { 
     int nx = i + dx[dir]; 
     int ny = j + dy[dir]; 
     r = min(r,boundedBy(nx,ny,k,in11,in22)); 
    } 
    return r; 
} 

void mark(int i,int j,int k,int in1,int in2) 
{ 
    if(i<0 || j<0 || i>=in1 || j>=in2) 
     return; 

    if(new2d[i][j]>=k) 
     return; 

    if(visited[i][j]) return ; 

    visited[i][j] = 1; 

    for(int dir = 0;dir<4;dir++) 
    { 
     int nx = i + dx[dir]; 
     int ny = j + dy[dir]; 
     mark(nx,ny,k,in1,in2); 
    } 
    dp[i][j] = max(dp[i][j],k); 
} 

struct node 
{ 
    int i,j,key; 
    node(int x,int y,int k) 
    { 
     i = x; 
     j = y; 
     key = k; 
    } 
}; 

bool compare(node a,node b) 
{ 
    return a.key>b.key; 
} 
vector<node> store; 

int getData(int input1, int input2, int input3[]) 
{ 

    int row=input1; 
    int col=input2; 
    int temp=0; 
    int count=0; 

     for(int i=0;i<row;i++) 
     { 
      for(int j=0;j<col;j++) 
      { 
       if(count==(col*row)) 
       break; 
      new2d[i][j]=input3[count]; 
      count++; 
      } 
     } 

    store.clear(); 
    for(int i = 0;i<input1;i++) 
     { 
      for(int j = 0;j<input2;j++) 
      { 
       store.push_back(node(i,j,new2d[i][j])); 
      } 
     } 
    memset(dp,0,sizeof(dp)); 

     sort(store.begin(),store.end(),compare); 

     for(int i = 0;i<store.size();i++) 
     { 
      memset(visited,0,sizeof(visited)); 

      int aux = boundedBy(store[i].i,store[i].j,store[i].key,input1,input2); 
      if(aux>store[i].key) 
      { 

       memset(visited,0,sizeof(visited)); 
       mark(store[i].i,store[i].j,aux,input1,input2); 
      } 

     } 

     long long result =0 ; 

     for(int i = 0;i<input1;i++) 
     { 

      for(int j = 0;j<input2;j++) 
      { 
       result = result + max(0,dp[i][j]-new2d[i][j]); 
      } 
     } 

     return result; 

} 


int main() 
{ 
    cin.sync_with_stdio(false); 
    cout.sync_with_stdio(false); 
     int n,m; 
     cin>>n>>m; 
     int inp3[n*m]; 
     store.clear(); 

      for(int j = 0;j<n*m;j++) 
      { 
       cin>>inp3[j]; 

      } 

    int k = getData(n,m,inp3); 
     cout<<k; 

    return 0; 
} 
2

Dieses Problem gelöst werden kann, die Priority-Flood-Algorithmus. Es ist schon einige Male in den letzten Jahrzehnten (und wieder von anderen Menschen diese Frage zu beantworten) entdeckt und veröffentlicht, obwohl die spezifische Variante Sie suchen ist nicht, mein Wissen in der Literatur.

Sie können eine Übersichtsarbeit des Algorithmus und seine Varianten here finden. Seitdem dieses Papier veröffentlicht wurde, wurde eine noch schnellere Variante entdeckt (link), sowie Verfahren, um diese Berechnung an Datensätzen von Billionen von Zellen durchzuführen (link). Ein Verfahren zum selektiven Durchbrechen niedriger/enger Teilungen wird diskutiert here. Kontaktieren Sie mich, wenn Sie Kopien von diesen Papieren wünschen.

Ich habe ein Repository here mit vielen der oben genannten Varianten; Weitere Implementierungen finden Sie unter here.

Ein einfaches Skript Volumen berechnen die RichDEM Bibliothek ist wie folgt:

#include "richdem/common/version.hpp" 
#include "richdem/common/router.hpp" 
#include "richdem/depressions/Lindsay2016.hpp" 
#include "richdem/common/Array2D.hpp" 

/** 
    @brief Calculates the volume of depressions in a DEM 
    @author Richard Barnes ([email protected]) 

    Priority-Flood starts on the edges of the DEM and then works its way inwards 
    using a priority queue to determine the lowest cell which has a path to the 
    edge. The neighbours of this cell are added to the priority queue if they 
    are higher. If they are lower, then they are members of a depression and the 
    elevation of the flooding minus the elevation of the DEM times the cell area 
    is the flooded volume of the cell. The cell is flooded, total volume 
    tracked, and the neighbors are then added to a "depressions" queue which is 
    used to flood depressions. Cells which are higher than a depression being 
    filled are added to the priority queue. In this way, depressions are filled 
    without incurring the expense of the priority queue. 

    @param[in,out] &elevations A grid of cell elevations 

    @pre 
    1. **elevations** contains the elevations of every cell or a value _NoData_ 
     for cells not part of the DEM. Note that the _NoData_ value is assumed to 
     be a negative number less than any actual data value. 

    @return 
    Returns the total volume of the flooded depressions. 

    @correctness 
    The correctness of this command is determined by inspection. (TODO) 
*/ 
template <class elev_t> 
double improved_priority_flood_volume(const Array2D<elev_t> &elevations){ 
    GridCellZ_pq<elev_t> open; 
    std::queue<GridCellZ<elev_t> > pit; 
    uint64_t processed_cells = 0; 
    uint64_t pitc   = 0; 
    ProgressBar progress; 

    std::cerr<<"\nPriority-Flood (Improved) Volume"<<std::endl; 
    std::cerr<<"\nC Barnes, R., Lehman, C., Mulla, D., 2014. Priority-flood: An optimal depression-filling and watershed-labeling algorithm for digital elevation models. Computers & Geosciences 62, 117–127. doi:10.1016/j.cageo.2013.04.024"<<std::endl; 

    std::cerr<<"p Setting up boolean flood array matrix..."<<std::endl; 
    //Used to keep track of which cells have already been considered 
    Array2D<int8_t> closed(elevations.width(),elevations.height(),false); 

    std::cerr<<"The priority queue will require approximately " 
      <<(elevations.width()*2+elevations.height()*2)*((long)sizeof(GridCellZ<elev_t>))/1024/1024 
      <<"MB of RAM." 
      <<std::endl; 

    std::cerr<<"p Adding cells to the priority queue..."<<std::endl; 

    //Add all cells on the edge of the DEM to the priority queue 
    for(int x=0;x<elevations.width();x++){ 
    open.emplace(x,0,elevations(x,0)); 
    open.emplace(x,elevations.height()-1,elevations(x,elevations.height()-1)); 
    closed(x,0)=true; 
    closed(x,elevations.height()-1)=true; 
    } 
    for(int y=1;y<elevations.height()-1;y++){ 
    open.emplace(0,y,elevations(0,y) ); 
    open.emplace(elevations.width()-1,y,elevations(elevations.width()-1,y)); 
    closed(0,y)=true; 
    closed(elevations.width()-1,y)=true; 
    } 

    double volume = 0; 

    std::cerr<<"p Performing the improved Priority-Flood..."<<std::endl; 
    progress.start(elevations.size()); 
    while(open.size()>0 || pit.size()>0){ 
    GridCellZ<elev_t> c; 
    if(pit.size()>0){ 
     c=pit.front(); 
     pit.pop(); 
    } else { 
     c=open.top(); 
     open.pop(); 
    } 
    processed_cells++; 

    for(int n=1;n<=8;n++){ 
     int nx=c.x+dx[n]; 
     int ny=c.y+dy[n]; 
     if(!elevations.inGrid(nx,ny)) continue; 
     if(closed(nx,ny)) 
     continue; 

     closed(nx,ny)=true; 
     if(elevations(nx,ny)<=c.z){ 
     if(elevations(nx,ny)<c.z){ 
      ++pitc; 
      volume += (c.z-elevations(nx,ny))*std::abs(elevations.getCellArea()); 
     } 
     pit.emplace(nx,ny,c.z); 
     } else 
     open.emplace(nx,ny,elevations(nx,ny)); 
    } 
    progress.update(processed_cells); 
    } 
    std::cerr<<"t Succeeded in "<<std::fixed<<std::setprecision(1)<<progress.stop()<<" s"<<std::endl; 
    std::cerr<<"m Cells processed = "<<processed_cells<<std::endl; 
    std::cerr<<"m Cells in pits = " <<pitc   <<std::endl; 

    return volume; 
} 

template<class T> 
int PerformAlgorithm(std::string analysis, Array2D<T> elevations){ 
    elevations.loadData(); 

    std::cout<<"Volume: "<<improved_priority_flood_volume(elevations)<<std::endl; 

    return 0; 
} 

int main(int argc, char **argv){ 
    std::string analysis = PrintRichdemHeader(argc,argv); 

    if(argc!=2){ 
    std::cerr<<argv[0]<<" <Input>"<<std::endl; 
    return -1; 
    } 

    return PerformAlgorithm(argv[1],analysis); 
} 

Es sollte geradlinig sein, dies zu Format, was auch immer 2d Array anzupassen sind Sie

In Pseudo-Code verwenden, die findet auf den vorhergehenden äquivalent:

Let PQ be a priority-queue which always pops the cell of lowest elevation 
Let Closed be a boolean array initially set to False 
Let Volume = 0 
Add all the border cells to PQ. 
For each border cell, set the cell's entry in Closed to True. 
While PQ is not empty: 
    Select the top cell from PQ, call it C. 
    Pop the top cell from PQ. 
    For each neighbor N of C: 
    If Closed(N): 
     Continue 
    If Elevation(N)<Elevation(C): 
     Volume += (Elevation(C)-Elevation(N))*Area 
     Add N to PQ, but with Elevation(C) 
    Else: 
     Add N to PQ with Elevation(N) 
    Set Closed(N)=True 
+0

Ich habe [Python-Implementierung] (http://www.jianshu.com/p/540ee7ec725b) des im Pseudo-Code in Ihrer Antwort beschriebenen Algorithmus (Text in Chinesisch) gesehen. Hier ist eine [leicht modifizierte Version] (https://ru.stackoverflow.com/a/729857/23044) (Text auf Russisch) – jfs

0

Hier ist der einfache Code für die gleichgeschlechtlichen

#include<iostream> 
using namespace std; 
int main() 
{ 
    int n,count=0,a[100]; 
    cin>>n; 
    for(int i=0;i<n;i++) 
    { 
     cin>>a[i]; 
    } 
    for(int i=1;i<n-1;i++) 
    { 
     ///computing left most largest and Right most largest element of array; 
     int leftmax=0; 
     int rightmax=0; 
     ///left most largest 
     for(int j=i-1;j>=1;j--) 
     { 
      if(a[j]>leftmax) 
      { 
       leftmax=a[j]; 
      } 
     } 
     ///rightmost largest 

     for(int k=i+1;k<=n-1;k++) 
     { 
      if(a[k]>rightmax) 
      { 
       rightmax=a[k]; 
      } 
     } 
     ///computing hight of the water contained- 

     int x=(min(rightmax,leftmax)-a[i]); 
     if(x>0) 
     { 
      count=count+x; 
     } 

    } 
    cout<<count; 
    return 0; 
} 
+0

Diese Frage betrifft die Erweiterung des Algorithmus auf 3D. – Blastfurnace

+0

Das Absenden von Blastofen in der Annahme "mach ... erweiterbar" ist eine unglückliche Formulierung, in der * erweitert * gewünscht ist: das beantwortet die Frage nicht. Schön von Ihnen, kommentierten Code zu präsentieren (schauen Sie sich Konventionen und Tools an - [dogygen] (http://www.stack.nl/~dimitri/doxygen/manual/docblocks.html), zum einen), aber mit einer Frage ohne eine Programmiersprache zu erwähnen, sollten Sie dies tun. Code-only-Antworten sind verpönt: Bitte fassen Sie Ihre Idee/Vorgehensweise zusammen (siehe [Wie stelle ich eine gute Frage?] (Https://stackoverflow.com/help/how-to-ask)). – greybeard