750 likes | 929 Views
External Memory Graph Algorithms and Applications to GIS. Laura Toma Duke University July 14 2003. Massive Data. Massive datasets are being collected everywhere Storage management software is billion-$ industry. Examples : Geography: NASA satellites generate 1.2TB per day
E N D
External Memory Graph Algorithms and Applications to GIS Laura Toma Duke University July 14 2003
Massive Data • Massive datasets are being collected everywhere • Storage management software is billion-$ industry Examples: • Geography: NASA satellites generate 1.2TB per day • WEB: Web crawl of 200M pages and 2000M links, Akamai stores 7 billion clicks per day • Phone: AT&T 20TB phone call database • Consumer: WalMart 70TB database, buying patterns (supermarket checkout)
I/O Model [AV’88] N = problem size B = disk block size M = memory size • I/O-operation: • movement of one block of data from/to disk • Complexity measure: number of I/Os • Fundamental bounds: Block I/O M Scanning: scan(N) = I/Os Sorting: sort(N) = I/Os • In practice B and M are big
Outline • I/O-efficient graph algorithms • Problems, techniques and results • Algorithms for planar graphs using graph separation • A GIS application: TerraFlow
Adj(v1) Adj(v2) Adj(v3) … I/O-Efficient Graph Algorithms • Input: G = (V,E) • Assume edge-list representation of stored on disk • Basic problems: • BFS, DFS, CC, SSSP, MST • Hard in external memory! • Lower bound: Ω(min{V, sort(V)}) (practically Ω(sort(V)) • Standard internal memory algorithms for these problems use O(E) I/Os G
DFS(u) Mark u For every v in Adj(u) If v not marked DFS(v) BFS and DFS • Internal memory: O(V+E) • External memory: • one I/O per vertex to load adjacency list Ω(V ) I/Os • one I/O per edge to check if v is marked Ω(E) I/Os • O(V+E)= O(E) I/Os
v SSSP and MST • Dijkstra’s algorithm • Maintain p-queue on vertices not yet included in SSSP • Repeatedly • DeleteMin(v) and relax each adjacent edge (v,u) if d(s,u) > d(s,u) + wvu then DecreaseKey(u, d(s,u) + wvu) • External memory: • one I/O per vertex to load adjacency list Ω(V) I/Os • External p-queue: O(E) Insert/Delete/DeleteMin in O(sort(E)) I/Os • DecreaseKey:O(1) I/Os to read key of u Ω(E) I/Os O(V+E+sort(E))= O(E) I/Os
I/O-Efficient Graph Algorithms • Problems: • Random (unstructured) accesses to the adjacency lists of vertices as they are visited Ω(V) I/Os • Need to check if v has been already visited and/or read its key Ω(E) I/Os • o(E) algorithm: solve (2) • o(V) algorithm: solve (1) and (2)
u v v o(E) Algorithms • Store edges to previously seen vertices Undirected/directed BFS, DFS, SSSP • buffered repository tree (BRT) [BGVW’00] Insert(v, e), ExtractAll(v) • Process/update all adjacent edges without checking if necessary Undirected SSSP: • I/O-efficient tournament tree [KS’96] DecreaseKey(v,k) Undirected MST: O(V + sort(E)) [ABT’01] • Maintain a priority queue on edges incident to current MST • How to decide if v is in MST without doing one I/O? • If next edge returned by DeleteMin is the same then v already in MST
u3 u1 u3 u1 u2 u4 u4 u2 o(V) Algorithms • CC and MST: [MR’99, ABT’01] • graph contraction • Goal: reduce the problem to the same problem on a smaller graph by selecting disjoint subgraphs and contracting them • A contraction phase reduces nb of vertices by a constant fraction • Typically use a sequence of contraction steps G = G0 G1 G2 … Gi … • CC and MST algorithms: general idea • Use contraction steps • Use an O(V+sort(E)) algorithm on G’
o(V) Algorithms • Undirected BFS, SSSP [MM’02, MZ’03] • Clustering • partition graph into V/k subgraphs (clusters) of kvertices • BFS Idea: Keep a pool of hot clusters • A cluster is loaded in the pool once • A cluster stays in the pool until all its vertices have been visited
Upper Bounds • General undirected graphs • CC, MST: [MR’99, ABT’01] • BFS: [MM’02] • SSSP: [MZ’03] • DFS: [KS’96] • General directed graphs • BFS, DFS, SSSP: [BVWB’00] Topological sort
Upper BoundsSparse Graphs • Sparse graphs E=O(V) • CC, MST: O(sort(V)) if graph stays sparse under edge contraction • Undirected BFS: O(sort(V)) ? open • Undirected SSSP: O(sort(V)) ? open • Undirected DFS: O(V)o(V) ? open Directed BFS, DFS, SSSP • O(sort(N)) BFS, SSSP, (DFS) on special classes of sparse graphs • Planar • Outerplanar, grid, bounded-treewidth
DFS O(sort(N)) I/Os [AMTZ’01] O(sort(N)) I/Os [ABT’00] O(sort(N)) I/Os [ABT’00] separators BFS SSSP Planar Undirected Graphs • BFS, DFS, SSSP: O(sort(N)) I/Os • O(sort(N)) I/O-efficient reductions [ABT’00, AMTZ’01] • Separators can be computed in O(sort(N)) I/Os [MZ’02]
I/O-Efficient Graph Algorithms Our Contributions • An MST on general undirected graphs. • O(sort(N)) algorithms on planar graphs • Reducibility on planar undirected graphs • Planar digraphs: SSSP, BFS, directed ear decomposition and topological sort • An O(sort(N) log N) DFS algorithm for planar undirected graphs • O(sort(N)) cycle separator • All-pair-shortest-paths and diameter • Planar digraphs • General undirected graphs • Data structure for shortest path queries on planar digraphs • Trade-off space-query • GIS application: TerraFlow • Flow modeling on grid terrains • r.terraflow: Port into GRASS, the open source GIS
Outline • I/O-efficient graph algorithms • Problems, techniques and results • Algorithms for planar graphs using graph separation • Shortest paths (SSSP, BFS, APSP) • DFS • Topological sort on planar DAGs • Data structure for SP queries • A GIS application: TerraFlow
R R R R R R R R R Planar graph separation:R-division • A partition of a planar graph using a set S of separator vertices into . subgraphs (clusters) Giof at most R vertices each such that: • There are separators vertices in total • There is no edge between a vertex in Giand a vertex in Gj • Each cluster is adjacent to separator vertices
R-division • Boundary vertices Bnd(Gi) of Gi • The separator vertices adjacent to Gi • Boundary set • Maximal subset of separator vertices that are adjacent to the same clusters • Lemma [Frederickson’87]: • R-division of a planar graph of bounded degree has boundary sets.
R-divisions and Planar Graph Algorithms • R-divisions [Frederickson’87] dynamic graph algorithms [GI’91,KS’93], faster SP algorithms [HKRS’97], SP data structures • In external memory choose R = B2 • O(N/B) separator vertices • O(N/B2) clusters of O(B2) vertices each and O(B) boundary vertices • O(N/B2) boundary sets • Can be computed in O(sort(N)) I/Os [MZ’02] • B2-division SSSP, BFS, DFS, topological sort, APSP, diameter, SP data structures,..
t s Planar SSSP B2 1. Compute a B2-division of G 2. Construct asubstitute graph GR on the separator vertices such that it preserves SP in G between any u,v in S • replace each subgraph Gi with a complete graph on Bnd(Gi) • for any u, v onBnd(Gi), the weight of edge (u,v) is δGi(u,v) GR has O(N/B2)· O(B2)=O(N) edges and O(N/B) vertices 3. Compute SSSP on GR 4. Compute SSSP to vertices inside clusters
v Planar SSSP SSSP on GRwith O(N/B) vertices and O(N) edges • Dijkstra’s algorithm with I/O-efficient p-queue • Access to adjacency list of each vertex takes O(N/B) I/Os • O(N) Insert/Delete/DeleteMin in O(sort(N)) I/Os [A95] • But..need dist(s,u) for all u in Adj(v) • Keep list LS={dist(s,u), for any u in S} • For each vertex v read from LS the current distances of adjacent vertices O(N) edges => O(N) accesses toLS O(N) I/Os
Planar SSSP SSSP on GR Idea: use boundary sets Store LS so that vertices in the same boundary set are consecutive • There are O(N/B2) boundary sets • Vertices in same boundary set have same O(B) neighbors in GR assuming G has bounded degree • Each boundary set is accessed once by each neighbor in GR • Each boundary set has size O(B) O(N/B2) x O(B) = O(N/B)I/Os
α u v β Gj Gi Planar APSP • Straightforward bound: O(N sort(N)) = O(sort(N2)) • Improved to optimalO(scan(N2)) • Idea: compute SP from all vertices in a cluster while cluster is in memory • For each cluster Gi • For any α in Bnd(Gi) compute SSSP(α) in GR • For each cluster Gj • load in memory Gj, Bnd(Gj) and δ(Bnd(Gi), Bnd(Gj)) • compute the shortest paths between all vertices in Gi and Gj d(u,v)=min{δGj (u,α) + δGR(α, β) + δGi(β,v) | α in Bnd(Gi), β in Bnd(Gj)} • write the output • O(N/B2) clusters O(sort(N2)/B) [compute] + O(scan(N2)) [output] • Diameter: O(sort(N2)/B)
General AP-BFS The APSP idea (compute SP from all vertices of a cluster while the cluster is in main memory) can be generalized to other algorithms which use clustering, like the BFS algorithm [MM’02] on general undirected graphs. Theorem: • AP-BFS of a general undirected graph and its unweighted diameter can be computed in O(V sort(E)) I/Os. Note: • general undirected BFS is O(sort(E)) amortized over V vertices
Planar DFS Idea: Partition the faces of G into levels around a source face containing s and grow DFS level-by-level • Levels can be obtained from BFS in dual graph • Structure of levels is simple (bicomps are cycles) • Rooting/Attaching: use that a spanning tree is a DFS-tree if and only if it has no cross edges A DFS-tree of a planar graph can be computed in O(sort(N)) I/Os
Planar Graphs • Shortest paths • generalize to digraphs: compute B2-division on the underlying graph • BFS, SSSP in O(sort(N)) • APSP (transitive closure) in O(scan(N2)) • diameter in O(sort(N2)/B) • DFS • Undirected • O(sort(N)) using BFS in the dual [O(sort(N) log N) direct algorithm using cycle separators] • Directed • The planar undirected DFS algorithms do not extend to digraphs • O(sort(N)) DFS? open
Outline • I/O-efficient graph algorithms • Problems, techniques and results • Algorithms for planar graphs using graph separation • Shortest paths (SSSP, BFS, APSP) • DFS • Topological sort on planar DAGs • O(sort(N)) using directed ear decomposition (DED) of its dual • Simplified algorithm using B2-division • Data structure for SP queries • A GIS application: TerraFlow
Directed Ear Decomposition (DED) • A directed ear decomposition of a graph G is a partition of G into simple directed paths P0, P1, …, Pk such that: • P0 is a simple cycle • endpoints of each Pii>0 are in lower-indexed paths Pj, Pl, j,l<i • internal vertices of each Pii>0 are not in any Pjj<i • G has a directed ear decomposition if and only if it is strongly connected (exist directed cycle containing each pair of vertices u,v). • Planar DED: O(sort(N)) I/Os
Planar Topological Sort using DED • Theorem [KK’79]: The directed dual of a planar DAG is strongly connected and therefore has a directed ear decomposition. • Idea: • Place vertices to the left of P0 before vertices to the right • Sort two sets recursively • Used in PRAM topological sort algorithm [KK93,K93] • PRAM simulation O(sort(N)log N) I/Os • Improved to O(sort(N)) by defining and utilizing ordered ear decomposition tree [ATZ’03]
B2 O(sort(N)) Topological Sort using B2-division Same idea as in planar SSSP algorithm • Construct a substitute graph GR using B2-division • edge from v to uon boundary of Gi if exists path from v to u in Gi • Topologically sort GR(separator vertices in G): • Store in-degree of each vertex in list L • Maintain list of in-degree zero vertices • Repeatedly: • Number an in-degree zero vertex v • Consider all edges (v,u) and decrement in-degree of u in L • analysis exactly as in SSSP algorithm O(scan(N)) if B2-division is given v
B2 t s O(sort(N)) Topological Sort using B2-division 5 • Problem: • Not clear how to incorporate removed vertices from G in topological order of separator vertices (GR) • Solution (assuming only one in-degree zero vertex s for simplicity): • Longest-path-from-s order is a topological order • Longest paths to removed vertices locally computable from longest-paths to boundary vertices F E C D B A 4 3 1 2
O(sort(N)) Topological Sort using B2-division • Compute a B2-division of G • Construct substitute graph GR using • Weight of edge between v and u on boundary of Gi equal to length of longest path from v to u in Gi 2. Compute longest path to each vertex in GR (same as in G): • Maintain list L of longest paths seen to each vertex • Repeatedly: • Obtain longest path for next vertex v in topological order • Consider all edges (v,u) and update longest path to u 3. Find longest path to vertices inside clusters analysis exactly as for planar SSSP algorithm O(scan(N)) if B2-division is given v
Outline • I/O-efficient graph algorithms • Problems, techniques and results • Algorithms for planar graphs using graph separation • Shortest paths (SSSP, BFS, APSP) • DFS • Topological sort on planar DAGs • Data structure for SP queries • A GIS application: TerraFlow
Data Structure for SP Queries on Planar Digraphs • Problem: pre-process a planar digraph into a data structure in order to answer efficiently distance (shortest path) queries between arbitrary vertices • Trade-off space-query: O(S) space, query = ? • The two extreme straightforward solutions: • O(N) space, O(sort(N)) I/O query • O(N2) space, O(1) I/O query • Related work: • Planar graphs: [Arikati et al, Djidjev, 1996] [Chen & Xu, 2000] • Space-query trade-off: forany S in [N, N2], S x Q = O(N2) • General graphs: • approx shortest paths [Cohen, Halperin, Zwick, …] • I/O-model : space, query [HMZ’99]
Data Structure for SP Queries on Planar Digraphs • Basic data structure [Arikati et al, Djidjev]: • Recursively, compute a separator and store for each vertex u in G the shortest path from u to all separator vertices. Space , query time, I/Os [HMZ’99] • Generalized to any S in [N, N2]: O(S) space,Q=O(N2/S) • Use R-division • S in [N, N3/2]: Store shortest paths between the separator vertices and compute shortest path in each cluster on the fly. • S in [N3/2, N2]: Pre-process each cluster as a basic data structure and for any vertex u in G store shortest paths from u to all separator vertices. • I/O-model • S in [N, N3/2]: ? • S in [N3/2, N2]: O(S) space, query using [HMZ’99]
α u v β Gi Gj Data Structure for SP Queries on Planar Digraphs • General framework: Compute an R-division. Store APSP between separator vertices. This uses space O(N2/R). Query: δ(u,v)=min{δGj (u,α) + δGR(α, β) + δGi(β,v) | α in Bnd(Gi), β in Bnd(Gj)} • Problems • Store APSP between separator vertices so that the O(R) distances δ(Bnd(Gi), Bnd(Gj)) can be retrieved efficiently in O(scan(R)) I/Os • Compute δGj (u,v) in O(scan(R)) I/Os Pre-process each cluster recursively • Compute δGj (u, Bnd(Gi)) in O(scan(R)) I/Os Pre-process each cluster into a data structure for answering all-boundary-SP queries
Data Structure for SP Queries on Planar Digraphs • Let G be a planar graph of size N and Bnd(G) its boundary of size O(N1/2). There exists a data structure that uses space O(N lg N) and answers all-boundary-shortest-path queries in O(N/B) I/Os. • Theorem: For any S in [N, N2/B] there exists a data structure which answers distance queries in I/Os and can be built in I/Os. The size is . if and if . • For S =Θ(N): O(N log2N) space and O(N/B) query • For any S/N = Ω(Nε) or S = Ω(N1 +ε) for some ε in (0,1] There exists a data structure of size O(S) which answers distance queries in I/Os and can be built in I/Os.
Outline • I/O-efficient graph algorithms • Problems, techniques and results • Algorithms for planar graphs using graph separation • Shortest paths • DFS • Topological sort on planar DAGs • Data structure for SP queries • A GIS application: TerraFlow
TerraFlow 3 3 3 3 2 2 2 2 4 4 4 4 7 7 7 7 5 5 5 5 8 8 8 8 7 7 7 7 1 1 1 1 9 9 9 9 DEM Representations TIN Grid Contour lines Sample points • Grids DEMs grid graphs • On grid graphs: BFS, SSSP, CC in O(sort(N)) I/Os
TerraFlow Example: LIDAR Terrain Data • Massive (irregular) point sets (1-10m resolution) • Relatively cheap and easy to collect Example: Jockey’s ridge (NC coast)
TerraFlow Modeling Flow on Terrains • What happens when it rains? • Predict areas susceptible to floods. • Predict location of streams. • Flow is modeled by computing two basic attributes from the DEM of the terrain: • Flow Direction (FD) • The direction water flows at a point • Flow Accumulation (FA) • Total amount of water that flows through a point if water is distributed according to the flow directions
TerraFlow Flow Accumulationof Panama
TerraFlow Panama Flow Accumulation: zoom
TerraFlow GIS Performance on Massive Data • GRASS (open source GIS) • Killed after running for 17 days on a 6700 x 4300 grid (approx 50 MB dataset) • TARDEM (research, U. Utah) • Killed after running for 20 days on a 12000 x 10000 grid (appox 240MB dataset) • CPU utilization5%, 3GB swap file • ArcInfo (commercial GIS) • Can handle the 240MB dataset • Doesn’t work for datasets bigger than 2GB
TerraFlow Flow Direction (FD) on Grids • On grids: Approximated using 3x3 neighborhood Problem: flat areas - Plateas and sinks • Goal: compute FD grid • Every cell has flow direction • Flow directions do not induce cycles • Every cell has a flow path outside the terrain
TerraFlow FD on Flat Areas • Plateaus • A cell flows towards the nearest spill point on the boundary of the plateau • Compute FD on plateaus using CC and BFS • Sinks • Route the water uphill out of the sink by modeling flooding:uniformly pouring water on terrain until steady-state is reached • Flooding removes (fills) sinks Assign uphill flow directions on the original terrain by assigning downhill flow directions on the flooded terrain
TerraFlow Flooding • Watershed: part of the terrain that flows into a sink • Sinks partition of terrain into watersheds watershed graph GT • Vertices are watersheds; add vertex for the “outside” watershed • Edge (u,v) if watersheds u,v are adjacent • Edge (u,v) labeled with lowest height on boundary between u and v • Flooding: Compute for each watershed u to the height hu of the lowest-height path in GT from u to the “outside” watershed. • the height of a path is the height of the highest edge on path
TerraFlow Flooding • Plane-sweep algorithm with a Union-Find structure • Initially only the outside watershed is done • Sweep watershed graph bottom-up with a horizontal plane • When hit edge (u,v) • If both watersheds u and v are done, ignore • If none is done, union them • If precisely one is not done, raise it at h(u,v) and mark it done • Theorem: Flooding and the FD grid can be computed in O(sort(N)) I/Os on a grid DEM of size N.
TerraFlow Flow Accumulation (FA) on Grids FA models water amount of flow through each cell with “uniform rain” • Initially one unit of water in each cell • Water distributed from each cell to neighbors pointed to by its FD • Flow conservation: If several FD, distribute proportionally to height difference • Flow accumulation of cell is total flow through it Goal: compute FA for every cell in the grid (FA grid) Theorem: The FA grid can be computed in O(sort(N)) I/Os.
TerraFlow http://www.cs.duke.edu/geo*/terraflow TerraFlow • TerraFlow: implementation of I/O-efficient FD and FA algorithms • Significantly faster on very large grids than existing GIS software • Scalable: 1 billion elements!! (>2GB data) • Allows multiple methods flow modeling • Implementation • C++, uses TPIE (Transparent Parallel I/O Environment) • Library of I/O-efficient modules developed at Duke • Experimental platform • TerraFlow, ArcInfo: 500MHz Alpha, FreeBSD 4.0, 1GB RAM • GRASS/TARDEM: 500MHz Intel PIII, FreeBSD/Windows, 1GB RAM