River landscapes and optimal channel networks

Contributed by Andrea Rinaldo, May 10, 2018 (sent for review March 14, 2018; reviewed by Armaghan Abed-Elmdoust, Angelika Steger, and Ramarathnam Venkatesan)
June 11, 2018
115 (26) 6548-6553


Optimal channel networks (OCNs) are a well-studied static model of river network structures. We present exact results showing that every OCN is a natural river tree where a landscape exist such that the flow directions are always directed along its steepest descent. We characterize the family of natural river trees in terms of certain forbidden structures, called “k-path obstacles.” We thus determine conditions for river landscapes to imply a rooted tree as its network of topographic gradients. Our results are significant in particular for applications where OCNs may be used to produce statistically identical replicas of realistic matrices for ecological interactions.


We study tree structures termed optimal channel networks (OCNs) that minimize the total gravitational energy loss in the system, an exact property of steady-state landscape configurations that prove dynamically accessible and strikingly similar to natural forms. Here, we show that every OCN is a so-called natural river tree, in the sense that there exists a height function such that the flow directions are always directed along steepest descent. We also study the natural river trees in an arbitrary graph in terms of forbidden substructures, which we call k-path obstacles, and OCNs on a d-dimensional lattice, improving earlier results by determining the minimum energy up to a constant factor for every d2. Results extend our capabilities in environmental statistical mechanics.
River networks can be viewed as rooted trees that, when extracted from fluvial landscapes by topographic steepest descent directions, show deep similarities of their parts and the whole, often across several orders of magnitude, despite great diversities in their geology, exposed lithology, vegetation, and climate (1). The large related body of observational data provides quintessential examples of real physical phenomena that can be effectively modeled by using graph theory (2, 3). River networks are in fact the loopless patterns formed by fluvial erosion over a drainage basin. Such patterns, referred to as “spanning trees,” are such that a directed graph route exists for flows from every site of the catchment to an outlet and are strictly related to the topographical surface whose gradients determine the edges (the drainage directions). A specific class of spanning trees, called optimal channel networks (OCNs), was obtained by minimizing a specific functional (4, 5), later shown to be an exact property of the stationary solutions of the general equation describing landscape evolution (6, 7). The static properties and the dynamic origins of the scale-invariant structures of OCNs proved remarkable (1). OCNs are suboptimal (that is, dynamically accessible given initial conditions and quenched randomness frustrating the optimum search) configurations of a spanning network mimicking landscape evolution and network selection (SI Appendix). Empirical and theoretical works have generally focused on the two-dimensional (2D) case, although recently (inspired by vascular systems), three-dimensional (3D) settings have also been analyzed for an OCN embedded in a lattice (8). Several exact results have been derived for OCNs (615).
A model for a river network is obtained by taking a reasonably dense set of points on a terrain and then joining each point to a nearby point downhill. Experimental observations suggest that OCNs should satisfy the following properties: (i) each portion of the drainage basin has a single output; and (ii) a geomorphological stability condition holds. OCNs are defined on a finite graph G, consisting of a set V(G) of vertices (or nodes), which correspond to sites in a drainage river basin, and a set E(G) of edges between nearby sites. Moreover, each node v is endowed with a height, hv. Water from each site/node v is entirely directed along the steepest descent—that is, toward the neighbor u of v which maximizes hvhu—and we write
where N(v)={u:uvE(G)} denotes the set of neighbors of node v. It is easy to see that the set of edges along which water flows is acyclic; in a single river basin, it will be an oriented spanning tree, with all water flowing to a unique root. Note that the notation N(v){v} also includes the possibility of finding nodes with no outlet. The second empirical property provides a connection between the gradient of the descent and the landscape-forming flux of water flowing along the path. More precisely, in the (geological) steady state, if the landscape-forming water flow rate mv at a site v increases, then the heights of the sites are modified by water erosion in such a way as to keep the product mv1/2Δhv constant. The exponent 1/2 is experimentally determined by local slope-area plots where one assumes that landscape-forming discharges are proportional to total contributing area and refers to the fluvial domain, thereby neglecting unchanneled parts of the landscape (1). Such an exponent also emerges as the leading term of a general nonlocal erosion term under reparametrization invariance and a small-gradient approximation (7). Thus, the first rule microscopically ensures that no cycles are created, while the second rule introduces a feedback between river structure and landscape erosion (1, 1618).
Here, we present analytical results about the properties of the structure that minimizes the total gravitational energy loss. In particular, we will give a theoretical characterization of river networks in terms of forbidden substructures, which provides a simple way of determining which of the spanning trees of a graph are natural river trees. We will also study OCNs embedded in a d-dimensional lattice, proving, for each dimension d2, upper and lower bounds on the energy that differ by only a constant factor. In fact, river trees are 2D, while their landscapes are 3D objects where a drop in elevation is assigned along edges as a function of connectivity. d-dimensional optimal constructions are of notable interest, like, e.g., spanning vascular networks (d=3) delivering metabolites to a given volume (13). In a vascular system of size L, metabolic rates scale like Ld, while mass M scales as MLd+1, leading to metabolism scaling as Md/(d+1) (13). Constructions in d>3 may consider time as an extra dimension [see, e.g., the work on Abelian sandpile (19) models of the Bethe lattice (20)]. Great interest thus exists for OCNs on general graphs for which important applications may arise in the future. We thus begin with a formal description of OCNs, and then we present the analytical results concerning the characterization of river trees and the energy of an OCN embedded in lattices of any dimension d.
For each constant γ(0,1), we can associate to the (rooted) spanning trees of a graph G an energy function of the following form:
where Av(T) is the number of vertices uV(G) such that the path in T from u to the root contains the vertex v (in other words, the number of vertices in the subtree of T “rooted” at v). A spanning tree of G which minimizes the energy in Eq. 1 is called an OCN of G.
In the case of river networks, if we imagine an open system where injection occurs at each site (vertex) at rate 1, and flowing along the edges of T until it reaches the root, then Av(T) represents the total flow of water out of v. The energy function defined in Eq. 1 has its origins in the exact result that steady-state configurations of river networks minimize the total loss of gravitational energy. This corresponds to:
where mv is the mass of water leaving vertex v, so mv is proportional to Av(T). Furthermore, experimental evidence suggests the following relationship between Av(T) and the difference of heights of two adjacent nodes: ΔhvAv(T)1/2 (1). We therefore obtain:
Replacing the exponent 1/2 with the parameter γ(0,1), and ignoring the (unimportant) constant factor, gives Eq. 1. When γ=1 the class of directed networks that minimize the mean distance to the outlet is obtained, whose energy function is the same (6, 16). This is true if the functional in Eq. 1 is replaced by a more general one where T is substituted with the graph G itself, including all edges. For such functional, every tree is a local minimum (7), and minimization is carried out under the constraint that locally (for each node) and globally a conservation law holds matching injection rates and inflows/outflows (6, 7, 16).
Given a rooted graph G (with one or more roots), we say that T is a “rooted spanning forest” of G if T is an acyclic subgraph of G, with one root of G in each component of T, and all edges oriented toward the corresponding root. The following natural definition will play an important role in our investigation of the OCN(s) of a graph.


Let G be a rooted graph, and T be a rooted spanning forest of G. We say that T is a “natural river spanning forest” of G if there is a height function f:V(G)R such that if vV(G) is a nonroot, and uN(v), then vu is an edge of T if and only if
We will assume throughout that the values of f are distinct, and hence that the minimum in Eq. 4 occurs at a unique vertex x. We remark that we will usually be interested in the case when G has exactly one root, in which case we say that T is a “natural river spanning tree” of G.
Note that (i) every nonroot vertex of T has out-degree exactly 1, (ii) f(v) strictly decreases along every directed path in T, and (iii) the function f attains its minimum at one of the roots. In Theorem 1, we will provide a simple characterization of the natural river spanning trees of a graph in terms of forbidden substructures. We will also show (Theorem 3) that every OCN in a graph is also a natural river spanning tree of that graph. As a consequence, to study (either analytically or numerically) the tree structures that minimize the energy function Eγ(T) defined in Eq. 1, it suffices to consider natural river trees.

Which Networks Are River Trees?

The definition of a natural river spanning tree creates various constraints on its possible shape and structure; in this section we shall consider these constraints.


Let G be a rooted graph, and let T a spanning tree of G. A k-path obstacle for T consists of a sequence of directed paths (η1,,ηk) in T, with the following property: For each iZk, the last vertex of ηi is connected to the first vertex of ηi+1 by an edge of G that is not in T.
In other words, if ηi is the path aici(1)ci(i)bi in T, and
then we say that T contains a k-path obstacle (Fig. 1). We emphasize that the paths ηi do not need to be vertex-disjoint (and it will be convenient in the proofs below to allow them to intersect); however, if k is minimal such that T contains a k-path obstacle, then the paths in any k-path obstacle will be disjoint. Observe that we will always consider the indices modulo k, so that (for example) b0=bk and b1=bk+1, and that if k=1 then necessarily 1>0.
Fig. 1.
k-path obstacles. (A) A three-path obstacle. For i{1,2,3}, ηi, is the directed path aici(1)ci(i)bi in T, where 1,2,3 are arbitrary nonnegative integers. Dashed red lines represent the edges b1a2, b2a3 and b3a1 belonging to G but not to T. (B) An example of a two-path obstacle (a1c1b1 and a2c2b2), where G (dashed red lines) is a regular 2D lattice and T (black arrows) is a spanning tree.
It is not difficult to show that a natural river tree cannot contain a k-path obstacle for any kN; our first theorem states that this simple necessary condition is also sufficient, and hence characterizes the natural river trees.

Theorem 1.

Let G be a finite graph with a single root. A spanning tree T of G is a natural river spanning tree of G if, and only if, it has no path obstacle.
Looking at Fig. 1, it may be deduced that a natural river tree by definition cannot contain any three-path obstacle (i.e., cannot contain the red dashed edges). The same reasoning can be extended to any natural river tree and value of k. We will give here only a formal proof of the necessary and sufficient conditions for T to be a natural river network.


We will show first that, for each k1, a natural river tree T cannot contain a k-path obstacle. Indeed, let f:V(G)R be the height function associated to T, and suppose that there exist paths η1,,ηk in T, where ηi is the path aici(1)ci(i)bi in T, and b1a2,b2a3,,bka1E(G)\E(T). Since aici(1) (or aibi if i=0) is in T, and aibi1E(G), it follows from the definition that f(bi)f(ci(1))<f(bi1) for each iZk, where b0bk. But then,
and this contradiction proves the claim.
We now turn to our main task: showing that if T has no path obstacle, then it is a natural river tree. Given such a rooted spanning tree T in G, we introduce a directed graph D with vertex set V(G) and
E(D)={xy:eitherxyis an edge ofT, orzsuch thatxzE(G)\E(T)and(zy)E(T)}.


The directed graph D is acyclic.


Suppose that D contains a directed cycle with exactly k edges that are not edges of T, and note that k1, since T is a tree. Let us denote the vertices of the cycle as follows:
where b1,,bk are the vertices whose out-neighbor in the cycle is not an edge of T, and 1,,k0. Thus, there exist vertices a1,,ak such that, for each iZk, we have bi1aiE(G)\E(T) (where b0bk), and aici(1) (or aibi if i=0) is an edge of T. This gives a k-path obstacle, contradicting our assumption.
It follows from the claim that the directed graph D can be extended to a linear ordering < on V(G) (i.e., a linear ordering with x>y for every edge xy of D). Choose a height function f such that f(x)<f(y) if and only if x<y; we claim that T is a natural river spanning tree with this height function.
Indeed, let vV(G) be a nonroot, and note that there is a unique uN(v) such that vu is an edge of T, since all edges of T are directed toward the root. Note that vu is an edge of D, by definition, so f(u)<f(v). Now let uwN(v), and observe that wu is an edge of D (since wvE(G)\E(T) and vu is an edge of T), so f(u)<f(w). Hence,
as required. It follows that a spanning tree of G with no path obstacle is a natural river spanning tree of G, completing the proof of Theorem 1.
Recall that an OCN (with parameter γ) in a graph G is a spanning tree of G that minimizes the energy Eγ(T). Our next theorem states that such a tree does not contain a path obstacle.

Theorem 2.

Let G be a finite graph with a single root, let T be a spanning tree of G, and let kN and γ(0,1). If T is an OCN in G with parameter γ, then T does not contain a k-path obstacle.
The proof of this theorem is presented in Materials and Methods. Theorems 1 and 2 immediately imply the following fundamental fact.

Theorem 3.

Every OCN is a natural river spanning tree.
In other words, if G is a finite graph with a single root, and T is an OCN in G with parameter γ(0,1), then T contains no path obstacle, by Theorem 2, and hence, by Theorem 1, is a natural river spanning tree of G.

Energy of OCNs in d-Dimensional Grids.

Here, we will study in more detail the energy of an OCN in a d-dimensional lattice. In particular, we will focus on the case of the d-dimensional grid [n]d, with root (1,,1), and edges between vertices at -distance 1, although it will be clear from the proofs below that our results can be easily extended to other lattices. For each n,dN and γ(0,1), define
T(n,d,γ)min{γ(T):Tis a spanning tree of[n]d},
the energy of an OCN in the d-dimensional grid. The following theorem determines T(n,d,γ) up to a constant factor (depending on d and γ) for all values of the parameters. It extends results of Colaiori et al. (12), who proved the lower bounds for d=2 and γ1/2 (in fact, on an orthogonal grid), and improves a result of Briggs and Krishnamoorthy (8), who obtained the lower bound 3n2/2O(n) in the case d=2 and γ=1/2.

Theorem 4.

For every integer d2, we have
as n.


Let us first prove the upper bounds, which follow from a straightforward recursive construction. We may assume for simplicity that n=2k, and consider 2d copies of the construction Tn/2 for [n/2]d, oriented so that 2d1 of them have their root near the center of [n]d, and the last one has its root at (1,,1), Fig. 2. Join the vertex (n/2,,n/2) in this last subcube to each of the roots of the other subcubes, which are each of the form (n/2,,n/2)+x for some x{0,1}d, so as to form a spanning tree Tn in [n]d.
Fig. 2.
The upper bound construction in the proof of Theorem 4. The tree is formed (in two dimensions) by four copies of the construction for [n/2]2, placed so that squares 1, 2, and 3 have their root adjacent to the center of [n]2, and the last one has its root at (1,1) (the black dot in the left bottom corner).
We claim that the subtree Tn of [n]d constructed in this way contains the diagonal edge (x+1,,x+1)(x,,x) for each 1x<n. Indeed, this holds (trivially) when n=1, and by induction it follows for all n=2k, since we add the diagonal edge (n/2+1,,n/2+1)(n/2,,n/2) in the induction step, and all other diagonal edges are in Tn by the induction hypothesis. It follows that all of the mass from 2d1 copies of [n/2]d flows down the diagonal of the final copy, and hence we have the following recursive formula:
Since (y+z)γyγ+zγ, and noting that (n/2)nd(n/2)dγn1+dγ, it follows that
and, hence, noting that Eγ(T1)=1, and recalling that n=2k, we obtain
This inequality immediately implies the claimed upper bounds. Indeed, if 1+dγ<d then the right-hand side of Eq. 8 is an increasing geometric series, so
If 1+dγ>d then it is a decreasing geometric series, so
Finally, if 1+dγ=d then all terms are equal, so
as required.
We next turn to the lower bounds. Observe first that the lower bound
holds trivially for every spanning tree T of [n]d, and every γ(0,1), since we have Av(T)1 for every vertex v. Hence, T(n,d,γ)nd holds always.
For the other lower bounds, we will find it useful to work in a more general setting, assuming for simplicity that n=3k. Let [n]d denote a copy of the d-dimensional grid in which every vertex on the boundary is a root and define
F(n,d,γ)minEγ(F):Fis a rootedspanning forest of[n]d},
where each component of F is assumed to have a single root (which lies on the boundary of [n]d), and the energy Eγ(F) is defined to be the sum of the energies of the components. Note that F(n,d,γ)T(n,d,γ), since a spanning tree in [n]d contains a rooted spanning forest of [n]d. We claim that
To prove Eq. 9, we partition [3n]d into 3d subcubes of size [n]d. In each subcube we must send the mass to the boundary of this subcube, and this gives a contribution of at least 3dF(n,d,γ) to the energy. Additionally, the middle subcube has to send its mass from its boundary to the boundary of the large cube [3n]d, and therefore each drop of water must pass through at least n additional vertices on its way to a root. Note that, since Av(T)(3n)d and γ<1, we have
for any x>0. Therefore, moving the mass of the center subcube to the boundary increases the energy by at least
and, hence, we obtain Eq. 9 (Fig. 3).
Fig. 3.
The proof of the lower bounds in Theorem 4. Black dots represent the root on the boundary of [3n]d for each component of the spanning forest. The red square highlights the central subcube; the water from this subcube must travel to a root on the boundary, and to do so it must first cross its own boundary and then pass through at least n additional vertices.
Setting C=γ3(d+1), and noting that F(1,d,γ)=1, we obtain
where in the penultimate term =k2. If 1+dγ<d then this gives only a minor improvement over the trivial lower bound nd, but for 1+dγ>d the improvement is more substantial. Indeed, it implies that
as n. Finally, if 1+dγ=d, then we obtain
as required.
An interesting (and likely difficult) challenge would be to decrease the implicit constant factors between the upper and lower bounds in Theorem 4 to factors of 1+o(1), and hence to determine T(n,d,γ) asymptotically.


When one considers only the structure of the network (without imposing a height function on the vertices), the requirement that a river network at stationarity minimizes total energy dissipation takes a simple and elegant form: We need to minimize the function
We have studied the energy and the structure of the OCN (the tree that minimizes Eq. 10), first on a general graph and then (in more detail) for the particular case of a discrete lattice (Fig. 4 and SI Appendix). Every OCN is a natural river tree, meaning that there exists a height function which determines the tree, and we characterized the natural river trees of a graph in terms of forbidden substructures. In the case of a d-dimensional lattice, we determined the energy of the OCN up to a constant factor for every d2 and γ(0,1), extending and improving previous work (8, 12).
Fig. 4.
OCN and landscapes. (A) An OCN derived within a 64 × 64 regular 2D lattice by means of a simulated annealing approach (SI Appendix) assuming γ=1/2. The algorithm ensures that the state is a local, dynamically accessible minimum, not necessarily the absolute minimum of the function Eγ(T). (B) The corresponding landscape (i.e., the elevation hv) computed assuming ΔhvAv(T)1/2. We have verified numerically that the obtained OCN is a natural river spanning tree.
We have also verified that, while it is natural to expect that any OCN on a discrete lattice should take the highly symmetric form of a Peano structure (1, 21), this is in fact not the case: For example, on the square lattice, several suboptimal minimum energy structures are obtained by breaking the Peano structure symmetry (1) (SI Appendix). Much effort has been devoted to reconciling the features of the ground state (known exactly since ref. 10) with those of feasible configurations. The resulting exponents associated with the global minimum did not match either the observational data or the numerical simulations (1). Because every local minimum of the OCN functional is a stationary solution of the general landscape evolution equation, any self-organization of the fluvial landscape corresponds to the dynamical settling of optimal structures into suboptimal niches of their fitness landscape. Thus, feasible optimality (i.e., the search for optima that are accessible to the dynamics given the initial conditions and path obstacles) is the rule in the fluvial landscape, and this might apply to a broad spectrum of interfaces in nature (16). Different network shapes give different values of Eγ(T) on a discrete lattice. Such a case imposes sharp conditions on the minimizing structures; in particular, a specific subset of all of the possible trees was selected as the natural rivers restricting the quest of minimum in this set.
As stated previously, great interest is placed in OCN on general graphs. In a first application of the theorems presented here, we have addressed the constructability of elevation fields and topographies compatible with the planar imprinting of OCNs (Fig. 4 and SI Appendix). Such results imply the possibility of building replicas of statistically identical matrices for ecological interactions. This is due to the flexibility of random search algorithms in a “greedy” mode (SI Appendix), that cannot access the absolute minimum but rather suboptimal configurations that are known to reproduce accurately, differently from the ground state, those of natural fluvial landforms forms (16). Interestingly, an application of limit scaling properties as a test of optimality has been proposed for foodwebs (22) and a direct use of OCN landscapes has been made to explain elevational gradients of biodiversity (17).

Materials and Methods

Proof of Theorem 2.

Let T be an OCN in G with parameter γ, and suppose (for a contradiction) that T contains a path obstacle. Let k be minimal such that there is a k-path obstacle in T, so there exist paths η1,,ηk in T, where ηi is aici(1)ci(i)bi, and
Now, for each iZk, define a new tree Ti by removing the edge aici(1) (or aibi if i=0) from T, and adding the edge aibi1 (where b0bk). If k=1, then T1 is a spanning tree of G with lower energy than T (since 1>0), contradicting our assumption that T is an OCN. Let us therefore assume from now on that k2.
To show that Ti is a tree, we need to use the minimality of k. Indeed, if there is a cycle in Ti then it must contain the edge aibi1, and therefore there exists a path
in T. But now we can construct a (k1)-path obstacle in T by replacing the paths ηi1 and ηi by the path
contradicting the minimality of k. Hence, Ti is a tree, and therefore, since T is an OCN, we have Eγ(T)Eγ(Ti) for each iZk.
We next bound Eγ(Ti) from above by introducing a new tree, Ti, by removing the edge aibi1 from Ti, and adding the edge aici1(1). Note that this may not be a subgraph of G; however, it is a tree [since, as above, there is no path in T from ci1(1) to ai, by the minimality of k], and therefore its energy Eγ(Ti) may still be defined via Eq. 1. Moreover, we have Eγ(Ti)Eγ(Ti), since Av(Ti)Av(Ti) for every vertex vV(G), and, hence, by the observations above,
For each iZk, let us write ξi for the (unique) path in T from ci(1) to the root of T. Observe that Av(Ti)=Av(T)wi for each vV(ξi)\V(ξi1), where wiAai, and that Av(Ti)=Av(T)+wi for each vV(ξi1)\V(ξi), so
Now, let us write Sv{iZk:vV(ξi)} for each vertex vV(G), and define VR{vV(G):Sv=R} for each set RZk. Note that the vertices in VR form a subpath of each ξi such that iR, and set
Since Eγ(T)Eγ(Ti), it follows from Eq. 11 that
Now, observe that fR is a concave function of x for any γ(0,1) and RZk, and note that fR(0)=0. It follows that there exists a constant αRR such that fR(x)αRx for every xR, and moreover fR(x)<αRx whenever x0 and VR. Hence, it follows from Eq. 12, and the fact that wi>0, that
and moreover, the inequality is strict if VR for any set R included in either sum. Hence, adding the αR with {i1,i}R to both sides of Eq. 14, we obtain
since ci(1)V{i}, so V{i}. Setting αiR:iRαR, it follows that α1>α2>>αk>α1, which is the desired contradiction.


P.B. and B.B. are partly supported by NSF Grant DMS 1600742. J.B. is partly supported by NSF Grant DMS-1500121 and an Arnold O. Beckman Research Award (University of Illinois Urbana–Champaign Campus Research Board 15006). B.B. and G.C. are partially supported by MULTIPLEX Grant 317532. R. Morris is partially supported by Brazilian National Council for Scientific and Technological Development Grant 303275/2013-8, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro Grant 201.598/2014, and European Research Council (ERC) Starting Grant 680275 MALIG. A.R. and E.B. have been supported by ERC Advanced Grant RINEC 22761. Part of the research in this paper was carried out while J.B. was a Visiting Fellow Commoner at Trinity College, Cambridge.

Supporting Information

Appendix (PDF)


I Rodríguez-Iturbe, A Rinaldo Fractal River Basins: Chance and Self-Organization (Cambridge Univ Press, Cambridge, UK, 2001).
B Bollobás Graph Theory, An Introductory Course (Springer, 1st Ed, New York, 1979).
B Bollobás Random Graphs (Academic, London, 1985).
I Rodríguez-Iturbe, A Rinaldo, R Rigon, R Bras, E Ijjasz-Vasquez, Fractal structures as least energy patterns: The case of river networks. Geophys Res Lett 19, 889–892 (1992).
A Rinaldo, et al., Minimum energy and fractal structure of drainage networks. Water Resour Res 28, 2183–2190 (1992).
J Banavar, F Colaiori, A Flammini, A Maritan, A Rinaldo, Topology of the fittest transportation network. Phys Rev Lett 84, 4745–4748 (2000).
J Banavar, F Colaiori, A Flammini, A Maritan, A Rinaldo, Scaling, optimality and landscape evolution. J Stat Phys 104, 1–33 (2001).
LA Briggs, M Krishnamoorthy, Exploring network scaling through variations on optimal channel networks. Proc Natl Acad Sci USA 110, 19295–300 (2013).
A Rinaldo, I Rodríguez-Iturbe, R Rigon, E Ijjasz-Vasquez, R Bras, Self-organized fractal river networks. Phys Rev Lett 70, 822–825 (1993).
A Maritan, F Colaiori, A Flammini, M Cieplak, J Banavar, Disorder, river patterns and universality. Science 272, 984–988 (1996).
A Rinaldo, A Maritan, F Colaiori, A Flammini, R Rigon, Thermodynamics of fractal networks. Phys Rev Lett 76, 3364–3367 (1996).
F Colaiori, A Flammini, A Maritan, JR Banavar, Analytical and numerical study of optimal channel networks. Phys Rev E 55, 1298–1310 (1997).
JR Banavar, A Maritan, A Rinaldo, Size and form in efficient transportation networks. Nature 399, 130–132 (1999).
F Santambrogio, Optimal channel networks, landscape function and branched transport. Inter Free Boundaries 9, 149–169 (2007).
A Abed-Elmdoust, A Singh, Z Yang, Emergent spectral properties of river network topology: An optimal channel network approach. Sci Rep 7, 11486 (2017).
A Rinaldo, R Rigon, JR Banavar, A Maritan, I Rodriguez-Iturbe, Evolution and selection of river networks: Statics, dynamics, and complexity. Proc Natl Acad Sci USA 111, 2417–2424 (2014).
E Bertuzzo, et al., Geomorphic controls on elevational gradients of species richness. Proc Natl Acad Sci USA 113, 1737–1742 (2016).
A Abed-Elmdoust, M Miri, A Singh, Reorganization of river networks under changing spatiotemporal precipitation patterns: An optimal channel network approach. Water Resour Res 52, 8845–3049 (2016).
P Bak, C Tang, K Wieselfeld, Self-organized criticality. An explanation of 1/f noise. Phys Rev Lett 59, 381–384 (1987).
D Dhar, S Majumdar, Abelian sandpile model on Bethe lattice. J Phys A: Math Gen 23, 4333–4350 (1990).
A Flammini, F Colaiori, Exact analysis of the Peano basin. J Phys A Math Gen 29, 6701–6708 (1996).
D Garlaschelli, G Caldarelli, L Pietronero, Universal scaling relations in food webs. Nature 423, 165–168 (2003).

Information & Authors


Published in

Go to Proceedings of the National Academy of Sciences
Go to Proceedings of the National Academy of Sciences
Proceedings of the National Academy of Sciences
Vol. 115 | No. 26
June 26, 2018
PubMed: 29891709


Submission history

Published online: June 11, 2018
Published in issue: June 26, 2018


  1. spanning trees
  2. graph theory
  3. slope-area law
  4. landscape evolution


P.B. and B.B. are partly supported by NSF Grant DMS 1600742. J.B. is partly supported by NSF Grant DMS-1500121 and an Arnold O. Beckman Research Award (University of Illinois Urbana–Champaign Campus Research Board 15006). B.B. and G.C. are partially supported by MULTIPLEX Grant 317532. R. Morris is partially supported by Brazilian National Council for Scientific and Technological Development Grant 303275/2013-8, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro Grant 201.598/2014, and European Research Council (ERC) Starting Grant 680275 MALIG. A.R. and E.B. have been supported by ERC Advanced Grant RINEC 22761. Part of the research in this paper was carried out while J.B. was a Visiting Fellow Commoner at Trinity College, Cambridge.



Paul Balister
Department of Mathematical Sciences, University of Memphis, Memphis, TN 38152;
József Balogh
Department of Mathematics, University of Illinois, Urbana, IL 61801;
Enrico Bertuzzo
Dipartimento di Scienze Ambientali e Informatica, Università Cà Foscari Venezia, 30123 Venice, Italy;
Béla Bollobás
Department of Mathematical Sciences, University of Memphis, Memphis, TN 38152;
Trinity College, Cambridge CB2 1TQ, United Kingdom;
London Institute for Mathematical Sciences, London W1K 2XF, United Kingdom;
Guido Caldarelli
Networks Unit, IMT (Istituzioni, Mercati, Tecnologie) School for Advanced Studies Lucca, 55100 Lucca, Italy;
Institute of Complex Systems, Consiglio Nazionale delle Ricerche, 00185 Rome, Italy;
European Centre for Living Technology, Università Cà Foscari Venezia, 30124 Venice, Italy;
Amos Maritan
Dipartimento di Fisica e Astronomia, University of Padua, 35122 Padua, Italy;
Istituto Nazionale di Fisica Nucleare, University of Padua, 35122 Padua, Italy;
Rossana Mastrandrea
Networks Unit, IMT (Istituzioni, Mercati, Tecnologie) School for Advanced Studies Lucca, 55100 Lucca, Italy;
Robert Morris
Instituto de Matemática Pura e Aplicada, Jardim Botânico, Rio de Janeiro, RJ 22460-320, Brazil;
Andrea Rinaldo1 [email protected]
Laboratory of Ecohydrology, École Polytechnique Fédérale de Lausanne, Lausanne CH-1015, Switzerland;
Dipartimento di Ingegneria Civile, Edile e Ambientale, University of Padua, 35122 Padua, Italy


To whom correspondence should be addressed. Email: [email protected].
Author contributions: P.B., J.B., B.B., and R. Morris obtained the mathematical proofs; and P.B., J.B., E.B., B.B., G.C., A.M., R. Mastrandrea, R. Morris, and A.R. determined the problem, discussed the results, and wrote the paper.
Reviewers: A.A.-E., University of Texas at Austin; A.S., ETH Zurich; and R.V., Microsoft Research Redmond.

Competing Interests

The authors declare no conflict of interest.

Metrics & Citations


Note: The article usage is presented with a three- to four-day delay and will update daily once available. Due to ths delay, usage data will not appear immediately following publication. Citation information is sourced from Crossref Cited-by service.

Citation statements



If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

Cited by


    View Options

    View options

    PDF format

    Download this article as a PDF file


    Get Access

    Login options

    Check if you have access through your login credentials or your institution to get full access on this article.

    Personal login Institutional Login

    Recommend to a librarian

    Recommend PNAS to a Librarian

    Purchase options

    Purchase this article to get full access to it.

    Single Article Purchase

    River landscapes and optimal channel networks
    Proceedings of the National Academy of Sciences
    • Vol. 115
    • No. 26
    • pp. 6513-E6098







    Share article link

    Share on social media