Articles | Volume 26, issue 2
Hydrol. Earth Syst. Sci., 26, 355–374, 2022
Hydrol. Earth Syst. Sci., 26, 355–374, 2022

Research article 24 Jan 2022

Research article | 24 Jan 2022

Hydrology without dimensions

Hydrology without dimensions
Amilcare Porporato1, Amilcare Porporato
  • 1Department of Civil and Environmental Engineering and High Meadows Environmental Institute, Princeton University, Princeton, USA
  • Invited contribution by Amilcare Porporato, recipient of the EGU John Dalton Medal 2020.

Correspondence: Amilcare Porporato (


By rigorously accounting for dimensional homogeneity in physical laws, the Π theorem and the related self-similarity hypotheses allow us to achieve a dimensionless reformulation of scientific hypotheses in a lower-dimensional context. This paper presents applications of these concepts to the partitioning of water and soil on terrestrial landscapes. For such processes, their complexity and lack of first principle formulation make dimensional analysis an excellent tool to formulate theories that are amenable to empirical testing and analytical developments. The resulting scaling laws help reveal the dominant environmental controls for these partitionings. In particular, we discuss how the dryness index and the storage index affect the long-term rainfall partitioning, the key nonlinear control of the dryness index in global datasets of weathering rates, and the existence of new macroscopic relations among average variables in landscape evolution statistics. The scaling laws for the partitioning of sediments, the elevation profile, and the spectral scaling of self-similar topographies also unveil tantalizing analogies with turbulent flows.

1 Introduction

Galileo is credited as the first scientist to have used dimensional analysis and scaling. In his 1638 “Dialogues Concerning Two New Sciences” (Galilei1914), he deduced that geometrically similar objects are not equally strong under their own weight: “A small dog could probably carry on his back two or three dogs of his own size; but I believe that a horse could not carry even one of his own size”. Since this discovery of scaling laws for complex biological materials, dimensional analysis has continued to fascinate many scientists, from Fourier and Maxwell to Reynolds, Rayleigh, Kolmogorov, and Taylor, and contributed to numerous new results in several fields (e.g., Barenblatt1996; Szirtes2007; Bolster et al.2011; Katul et al.2019). These methods have been extremely useful in complex problems lacking closed-form solutions due to nonlinear interactions and multi-physics as well as in design of experiments and numerical simulations and the rational interpretation of their results.

Looking at some of the most striking applications of the Π theorem and self-similarity (e.g., the famous example of the atomic bomb of Taylor; Taylor1950a, b; Barenblatt1996, and the Kolmogorov spectrum of turbulence; Kolmogorov1941; Katul et al.2019), it is easy to be lured by the promise of a mathematical structure offering a powerful dimension reduction in the space of variables and parameters, even for vaguely formulated problems. In spite of the straightforward steps for its application, a brute force approach to dimensional analysis rarely leads to useful results (and too often leaves one with too many dimensionless groups). These failures are probably at the origin of some of the backlash in the literature, hailing these methods as nothing more than fancy tricks, capable only of recasting already-known solutions. What is true is that applications of dimensional analysis cannot be done automatically but require careful consideration of the problem at hand. The results, which follow from the initial hypotheses, have to be scrutinized using the available data or simulations. In other words, when performing dimensional analysis, one cannot avoid the necessary iterative process of hypothesis formulation and subsequent verification, which forms the basis of any mathematical modeling (Logan2013). Its results, like a good dish, depend on its ingredients as much as on its recipe.

The presence of emerging scaling laws in geophysics has been widely recognized for a long time, and the related power laws have been observed in rainfall and streamflow statistics, in landscape and river-network geometry, as well as in the aggregation properties of soils and aquifers (among many others, see, e.g., Rodríguez-Iturbe and Rinaldo2001, Gagnon et al.2006, and Sposito2008, and references therein). Notwithstanding these numerous examples, the presence of only a few applications of the Π theorem in geophysics appears to be at odds with the early conclusion of Strahler (1958) that “dimensional analysis will become increasingly useful in empirical, quantitative studies in geomorphology by offering a systematic means of describing and comparing the form elements of the landscape”.

In this paper, we revisit a series of fundamental problems in hydrology and geomorphology and scrutinize them under the common lens of dimensional analysis with the goal of sharpening our intuition of the underlying physical processes. After a brief review of the concepts of dimensional analysis and scaling in Sect. 2, in Sect. 3 we consider three different instances of the partitioning of water and soils in natural landscapes. Because of their complexity, availability of global data, and the lack of governing equations from first principles, these phenomena provide a particularly fertile test bed for dimensional techniques.

2 From group theory to street-fighting hydrology

Dimensional analysis is a chapter of the more general group theory (Gilmore2012), based on the elegant formalism of generalized homogeneous functions (Barenblatt1996) and their underlying linear algebra (Logan2013). It formalizes the principle of dimensional homogeneity, which expresses the fact that physical equations should be valid independently of the system of units chosen. Taking this principle to its rigorous consequences allows us to exploit the dimensional symmetries of the problem in a way that may lead to useful results.

On the theoretical front, the Π theorem provides a systematic recipe to find self-similar solutions of partial differential equations (PDEs); these underlie the local Lie groups under which the PDEs are invariant (Barenblatt1996; Bluman and Cole2012). On the more practical side, in data analysis and design of experiments (Barenblatt1996; Szirtes2007; Shen et al.2014), dimensional analysis helps make the so-called “Fermi reasoning” (i.e., privileging good reasoning to accuracy to achieve fair estimates) more methodical (Bhaskar and Nigam1990; Persico2004; Efthimiou and Llewellyn2007). By offering a rigorous way to organize physical hypotheses about a problem, perform dimension reduction in terms of fundamental governing groups, and verify and refine the hypotheses with available information, dimensional analysis is an asset for “street-fighting” mathematical guessing (Mahajan2010).

2.1 Scaling and power laws

The word scaling, in the sense used here, refers to the existence of a property that allows one to go from one scale to another (upscaling or downscaling) using the same mathematical law and, therefore, it is related to the absence of a preferred scale or unit of measure. From a mathematical viewpoint, this property is linked to the fact that the solutions of dimensional physical problems are generalized homogeneous functions (Hankey and Stanley1972; Widom2009) and that dimension functions are power-law monomials (Barenblatt1996).

The defining property of a homogeneous function of degree n,

(1) f ( λ x ) = λ n f ( x ) ,

makes it evident that, apart from a scale coefficient, λn, the same function is used to go from one scale (x) to another scale (λx). The linkage to power laws becomes evident when setting λ=1/x, which transforms the previous equation in

(2) f ( x ) = x n f ( 1 ) = a x n ,

where a=f(1) is a constant. Power-law functions, f(x)=axn, satisfy Eq. (1)1.

As we will see, scaling laws appear naturally in the application of the Π theorem and are not limited to space and time variables but can involve any dimensional quantity. Scaling is often connected to fractal properties of the underlying processes and also appears in the related fields of critical phenomena and anomalous scaling (Hankey and Stanley1972; Sornette2006), renormalization group theory (Goldenfeld2018), and complete and incomplete self-similarity (Barenblatt1996). Sornette (2006) has several examples of power laws ensuing from different mechanisms.

2.2 The powerful dimension reduction of the Π theorem

The starting point of a dimensional analysis is the formulation of a “physical law”,

(3) a = f a 1 , , a n ,

which relates a quantity of interest, a, the so-called governed quantity, to n other governing quantities, a1, …, an. Equation (3) embeds formally our scientific hypothesis about the physical problem and serves as a mathematical placeholder to collect the initial ingredients, on which the recipe of the Π theorem then operates.

With Eq. (3) established, the next step is to take stock of the dimension (i.e., the factor by which a physical quantity changes upon passage from the original system of units to another) of the quantities present in the physical law (Eq. 3). According to Maxwell's convention, the dimension of a variable q is indicated as [q]. For example, the dimension of the velocity v is [v]=LT-1 (as is well known, in mechanics the dimensions of length mass and time – LMT – are used as a class of system of units). As a result of dimensional homogeneity, dimension functions are always power-law monomials (Barenblatt1996); thus, in mechanics a generic variable q has dimension function [q]=LαMβTγ. The latter is a specific instance of Bridgman's equation (Bridgman1922; Panton2006), in which LMT are the chosen principal dimensions. It follows that the dimension of the argument of transcendental functions is unity (i.e., they are dimensionless quantities), so that their numerical value is identical in all systems of units.

The number k of fundamental quantities in the physical law (Eq. 3) plays a key role in the Π theorem. Such fundamental quantities, also called repeating variables, must be dimensionally independent (Barenblatt1996): that is, none of their dimensions can be represented in terms of a product of powers of the dimensions of the remaining k−1 quantities (in turn, this is true only if the determinant formed with the exponents of the dimension functions is different from zero). Often k is equal to the fewest independent dimensions required to specify the dimensions of all quantities involved in Eq. (3), but in a few cases the number of primary dimensions differs when the variables are expressed in terms of different systems of dimensions (e.g., LMT, FLT – force, length, and time – or any other combination). In any case, the value of k is given by the rank of the dimensional matrix, formed by listing all the exponents of the primary dimensions of each variable in Eq. (3); see, e.g., Panton (2006) and Logan (2013). Once it has been ascertained that the problem admits k dimensionally independent quantities among the governing quantities, the k repeating variables may be chosen and the physical law (Eq. 3) may be formally re-arranged as

(4) a = f a 1 , , a k ; a k + 1 , , a n ,

where conventionally the semicolon separates the repeating variables from the other ones.

Enforcing dimensional homogeneity in Eq. (4) leads to the main outcome of the Π theorem (Barenblatt1996), namely, that the physical law, if true, can be written in the form

(5) Π = φ Π k + 1 , , Π n

in terms of the dimensionless Π groups

(6) Π = a a 1 α 1 , a 2 α 2 , , a k α k , Π k + 1 = a k + 1 a 1 α 1 , k + 1 , a 2 α 2 , k + 1 , , a k α k , k + 1 , Π n = a n a 1 α 1 , n , a 2 α 2 , n , , a k α k , n .

A comparison of Eqs. (3) and (5) clearly shows the great achievement of the Π theorem in re-expressing a general mathematical relationship between a quantity of interest and n dimensional quantities as a new relationship between nk dimensionless quantities in a more manageable, lower-dimensional space. The most striking applications in the literature (e.g., Barenblatt1996) are in fact linked to a drastic dimension reduction (in the problem of the pendulum n=2 and k=2, while for the atomic bomb n=3 and k=3, so that n-k=0 and one is left with a dimensionless function which is a constant). In nonlinear PDEs this may allow them to be transformed into ordinary differential equations (ODEs), with much greater chances of finding a solution, either analytically or numerically. While this happens more frequently in one spatial dimension (e.g., Barenblatt1996; Daly and Porporato2004b; Eggers and Fontelos2008), sometimes it also works in higher dimensions (Hills and Moffatt2000; Xue and Stone2020).

Depending on the governing variables involved in the physical law, there may be freedom in choosing the k governing repeating variables. Each admissible choice leads to dimensionless groups that are related to the ones obtained from a different admissible set. As we will see, while in general the different ensuing representations are equivalent, some combination may be more revealing of the underlying dynamics and afford a more parsimonious representation. In PDEs, this freedom may lead to different phase-space representations, some of which could be more amenable to analysis (Gratton and Minotti1990; Daly and Porporato2004b).

In summary, starting from a physically meaningful law (1), the Π theorem not only provides a mathematically more specific and elegant expression for it, but also helps reveal the actual physical controls of the problem, which emerge through it in the form of the Π numbers obtained: these, and not the single variables listed in the original law, are the actual quantities governing the physical phenomenon (for example, the Reynolds, Mach, and Froude numbers and the many other ones, including those that we will see later in this article).

2.2.1 Self-similarity

When one or more of the Π groups attain very large or very small values, the function φ in Eq. (5) may reach an asymptotic form related to a self-similar regime. In the simplest form of self-similarity, called complete or of the first kind (see Barenblatt1996), the function φ reaches a constant plateau for either very small or very large values of the governing groups. As a result, the physical problem does not change even if the values of these groups change, and the “self-similar group” can then be eliminated from Eq. (5), allowing for further dimension reduction. In more complicated cases the self-similarity is of the second type, or incomplete, because the function φ still depends on the self-similar Π group according to a power law with an exponent, which is not directly related to any of the dimension functions of the governing quantities.

Thus, assuming for example self-similarity with respect to the group Πk+1, Eq. (5) becomes

(7) Π = lim Π k + 1 0 or φ Π k + 1 , Π k + 2 , , Π n = Π 1 β ψ Π k + 2 , , Π n ,

where β=0 in the case of complete self-similarity. Thus complete self-similarity allows us to further reduce the dimensionality of the problem by as many dimensions as there are self-similar groups. This type of similarity is frequently encountered in near-wall turbulence, where the global Reynolds number and the dimensionless distance from the wall appear as self-similar groups. As we will see later on, complete self-similarity is also present in landscape channelization, when fluvial erosion dominates over soil creep.

Incomplete self-similarity seems to be more rare and is often difficult to distinguish from the case of complete self-similarity, especially for experimental or numerical problems where there are transitions to different regimes (Spagnoli2005) or in which the available data are not sufficiently precise (Barenblatt et al.1997; Smits et al.2011; Yin et al.2019). In self-similar solutions of PDEs, incomplete self-similarity typically results from a nonlinear eigenvalue problem and leads to power laws that control the temporal or spatial behavior of the problem with exponents that are irrational numbers that are not directly related to the dimension functions of any of the variables or initial and boundary conditions (e.g., Gratton and Minotti1990; Aronson and Graveleau1993; Barenblatt1996; Burton and Taborek2007; Zheng et al.2014). It is also worth noting that, starting from different physical laws with or without a certain variable, one sometimes arrives at different final forms of self-similar behavior. As we will see in the example of weathering in Sect. 3.2, it appears that such cases are actually related to problems of complete self-similarity, because the alternative formulation as an incomplete self-similarity problem is characterized by integer exponents, which allow simplifications among variables, which in turn lead to the complete self-similarity obtained with the other choice of variables.

2.3 Augmented and directional dimensional analysis

It is not infrequent in the literature to come across a point of view, made explicit by Bridgman (1922), that “there is nothing sacrosanct” about the choice of primary dimensions and that “dimensional analysis is merely a man-made tool that may be manipulated at will”. This is indeed in line with modern physics, which links length, time, mass, and energy to the different descriptions of reality depending on the scale of observation. Once a free choice of the primary units is accepted, then the question remains of what choice will be of maximum utility (Moon and Spencer1949).

Compared to the aforementioned (Sect. 2.2) freedom of choosing repeating variables or classes of systems of units (e.g., the length, force, and time, LFT, instead of the length, mass, and time, LMT), the augmented dimensional analysis refers to a more drastic freedom of “inventing” the type and number of primary dimensions. It is related to the original observation that “dimensions of quantities do not always afford a test of their identity” (Lodge1888). Since this ultimately affects the number of dimensionless groups and the extent of dimension reduction, the practical implications may be significant.

As the reader may imagine, this line of reasoning has attracted both stern skepticism and enthusiastic support, leading to interesting debates and controversies, with defenders of rigor and objectivity on one side and advocates of a flexible approach on the other. In the writer's experience, the subject is always a source of interesting discussion, if not else because, when one writes in favor of it, the editors somehow always manage to select a reviewer who belongs to the skeptical camp.

The most emblematic case of controversy is perhaps the well-known Rayleigh–Riabouchinsky controversy (e.g., Gibbings1982; Butterfield2001). Dealing with a problem of heat transfer between a body and a fluid stream, Rayleigh (1915a) solved it within the domain of thermodynamics, i.e., including temperature as a primary dimension and obtaining one governing dimensionless group. Riabouchinsky (1915) objected that, adopting the more advanced point of view of statistical mechanics, according to which temperature is the mean molecular kinetic energy, one could more parsimoniously limit the primary dimensions to length, time, and energy without involving temperature. As a result, he obtained two dimensionless groups, reaching the paradoxical conclusion that apparently more detailed knowledge of the problem yields a less informative result.

The resolution of this controversy (see Appendix A for more details) shows that it is not a matter of the arbitrariness of which dimensions are considered but of the level of description intended for a problem. Whether to include the specifics of the molecular motion depends on the size of the considered object. If one, following Rayleigh, accepts a thermodynamic approach, then the details of the disorderly molecular energy are irrelevant and temperature can be treated as an independent quantity compared to the kinetic energy of the mean motion. Formally, choosing a greater number of fundamental units (i.e., the number of primary dimensions) is made possible by the addition of corresponding dimensional unifiers (Panton2006). The confusion often arises in those cases where the structure of the equation is such that the dimensional unifier can be tacitly eliminated or taken for granted.

Of a similar nature, and perhaps even more subtle and controversial, is directional dimensional analysis, which is based on the fact that in some cases distinguishing between vertical and horizontal dimensions provides more informative results (i.e., fewer dimensionless groups). Williams (1890) argued that “owing to the dimensions of space, the unit of length is involved in different ways, according to the different relative directions in which it may be taken.” While for some authors it remains controversial (see, e.g., Barr1984; Gibbings1980; Kader and Yaglom1990, and references therein), others have worked to link these extensions of dimensional analysis more solidly to group theory (Moran and Marshek1972). Directional dimensional analysis works, provided that the ratio of the length dimensions does not play a role in the physics of the problem (e.g., in the original equations); this happens, for example, in cases that assume incompressibility where one of the lengths may be simply related to mass flow because the mass dimension has been canceled out of the equations. As we will see, the variable z in the landscape evolution model, analyzed in Sect. 3.3, falls within this category.

Directional dimensional analysis and its generalizations have been used successfully in a variety of problems, including applications in mechanics and atmospheric sciences (e.g., Huntley1967; Moran and Marshek1972; Kader and Yaglom1990; Siano1985; Daly and Porporato2004a, b; Dimitrakopoulos and DeJong2012; Bonetti et al.2020; Hooshyar et al.2020; Sun2020; Hooshyar et al.2021). Notwithstanding we do not have rigorous criteria to assess whether and when the “tricks” of augmented dimensional analysis can be applied, there are certainly several cases in which extending the set of primary dimensions is useful. Even in case of failure, the negative results can always be used to sharpen our physical laws and improve our starting point for dimensional analysis.

3 Water and soil mineral partitioning in the critical zone

Some of the most important questions of terrestrial geophysics are related to the partitioning of water and soil minerals at the land surface. Figure 1 shows the three cases that we will analyze in this section: the partitioning of rainfall into evapotranspiration and percolation plus runoff, the soil partitioning of minerals either lost by chemical dissolution (weathering) or transported away, and the related geomorphologic partitioning of soil over the landscape due to soil creep and fluvial erosion. Since their complexity prevents us from writing the governing equations and boundary conditions in detail, these problems are excellent candidates for dimensional analysis, applied in combination with available data and numerical simulations of simplified models. For simplicity, we focus on the dominant components of such partitionings and consider only long-term averages, assuming stationary conditions.

Figure 1The three partitionings analyzed in this paper using dimensional analysis: the rainfall partitioning taking place at the land surface, the soil mineral partitioning by chemical dissolution (weathering) and transport processes, and the related partitioning of soil sediments responsible for landscape evolution and the formation of drainage networks.


The scaling laws obtained from dimensional analysis shed light on the dominant soil, vegetation and climate controls for these hydrologic partitionings. Specifically, in the first application it will be seen how the dimensionless dryness and storage indices determine the long-term rainfall partitioning. The second application will reveal the key nonlinear control of the dryness index on weathering rates, while the third application will analyze new macroscopic relations among average variables' landscape evolution, uncovering an intriguing analogy between self-similar topographies and turbulent flow fields.

3.1 Rainfall partitioning

The rainfall reaching the soil surface is either lost by runoff or infiltrates into the soil, where in turn it is lost either by evapotranspiration or percolation. The fate of this partitioning is essentially controlled by the properties of the soil–plant system, which – in a sense – acts as a geophysical valve, not only for the entire hydrologic cycle, but also for the energy and carbon cycles. This fundamental hydrological problem presents robust behaviors for its macroscopic (i.e., averaged) patterns as well as rich and complex controls at the detail level, where the role of temporal fluctuations and spatial heterogeneities becomes important.

We will indicate the mean rainfall rate as R and consider together the main losses: the evapotranspiration rate (ET) and the mean rate of percolation plus runoff (LQ); see Fig. 1. Assuming stationarity over long timescales and referring to a given area, the input balances the outputs:

(8) R = ET + LQ .

It should be clear that, depending on the area chosen and whether the latter is homogeneous in terms of land cover, soil properties, and topography, the actual meaning of the term LQ in Eq. (8) may be very different. In particular, if the control volume for our water balance is the rooting zone of a small homogeneous plot of vegetated soil, then LQ will be the average of a very intermittent term given by the sum of surface runoff and percolation to soil layers below the rooting zone, while if we consider an entire river basin, it will be mostly the average of the streamflow draining the area. How these processes scale with the land area is an interesting open question, which is outside of our scope (see Fig. 4 by Yin et al.2019, for an analysis of this effect). We specifically focus on ET, as many have done before, but here we adopt the special lens offered by the Π theorem to explore the implications of different hypotheses in the physical laws used as starting points.

3.1.1 Turc and Budyko spaces: dryness and humidity as the main controls of rainfall partitioning

While they did not use dimensional analysis, Turc and Budyko started their work by making what is perhaps the simplest hypothesis of a physical law for the rainfall partitioning (see Dooge1992; Daly et al.2019, and references therein),

(9) ET = f TB R , ET max ,

where ETmax is a reference (e.g., potential or maximum) evapotranspiration. The rank of the dimension matrix

(10) ET R ET max L 1 1 1 T - 1 - 1 - 1

is 1, which leads to a total of two dimensionless groups, one governed and one governing. Note that the same result is obtained adopting a system of units with the only dimension of flux, say Φ, for which the dimension matrix

(11) ET R ET max Φ 1 1 1

also has rank 1.

If one chooses R as the repeating variable, applying the Π theorem then gives the so-called Budyko hypothesis

(12) Π B = φ B D I ,


(13) Π B = ET R and D I = ET max R

is Budyko's dryness (or aridity) index. If instead one chooses ETmax, then Turc's hypothesis follows (Daly et al.2019):

(14) Π T = φ T H I ,


(15) Π T = ET ET max and H I = R ET max = 1 D I

is the humidity index. Which representation to use is mostly a matter of convenience, since they are related by ΠT=ΠB/DI and φT=φB/DI. However, Budyko's hypothesis may be more suitable for emphasizing the dry end of the hydrologic spectrum, while Turc's privileges the humid end.

Budyko's hypothesis is tested in Fig. 2 using the MOPEX data for several river basins in the continental USA and plotted along with the original Budyko interpolating curve. Clearly, as shown many times before, the main controls on the rainfall partitioning are captured by the dryness index.

Figure 2Rainfall partitioning. (a) Test of Budyko's hypothesis using the MOPEX data (Porporato and Yin2022). The solid line shows the original semi-empirical curve used by Budyko, ΠB=(DItanhDI-1(1-e-DI))1/2. (b) Different solutions of the mean partitioning of the minimalist soil moisture balance model of Porporato et al. (2004), showing the role of the storage index, γ=w0/α, and seasonality (Feng et al.2012). The red dotted line is the original curve of Budyko.


3.1.2 Storage index: the role of the hydrologic active depth and the variability timescale

While the dryness index captures the main variability of the data, the scatter in Fig. 2 also suggests additional controls, which prompt us to revisit the physical law (Eq. 9) by considering additional governing quantities. There are several quantities that capture the different properties of the hydrologic system and its climatic forcing, such as storage capacity, areal extension, vegetation type, and seasonality. Among these candidates, the storage depth (the storage volume divided by the area) is probably the first one to consider, as suggested by the dramatic effects of paving a vegetated field on the rainfall partitioning.

To account for this effect, we may add a quantity w0 with dimensions of length (L), which is a volume divided by an area. This leads to a physical law of the type

(16) ET = f ET R , ET max , w 0 .

However, when trying to apply the Π theorem, it becomes immediately apparent that, if w0 is chosen as a repeating variable, it actually drops from the dimensionless formulation (similarly to the way the mass disappears from the list of variables when applying the Π theorem to the pendulum – Barenblatt1996). Physically, this tells us that to see a difference in the partitioning among the hydrologic balance of a parking lot, a crusted soil, and a deep vegetated soil, the role of the soil depth must be associated with other variables, which include a timescale to account for the variability of the hydrologic balance, for otherwise all the other influences are already contained in the dryness index.

Thus, to achieve a further refinement in the Pi-theorem formulation, in combination with the storage depth w0, Porporato et al. (2004) introduced a timescale related to the mean time between rainfall events or – better – its inverse, the frequency of rainfall events, say λ. Since the frequency of rainfall times the mean rainfall depth per event α is equal to the mean rainfall rate R=αλ, the presence of the mean rainfall rate in the list becomes redundant if one already includes the mean rainfall depth α. With these additions, the physical law becomes

(17) ET = f ET λ , α , ET max , w 0 ,

with the dimension matrix of rank 2 given by

(18) ET λ α ET max w 0 L 1 0 1 1 1 T - 1 - 1 0 - 1 0 .

It is instructive to explore the full range of choices related to the three possible dimensionless groups given by the five variables minus the rank 2 of the matrix. In turn, this implies five possibilities (of the six potential pairs of repeating variables, (αw0) obviously must be excluded because of dimensional dependence), which can all be brought back to the dryness index DI and the storage index γ=w0/α, as shown in the following table:

(19) Repeating var . Pi theorem Group relation ( α , λ ) ET α λ = φ 1 E max α λ , w 0 α Π B = φ 1 D I , γ λ , E max ET E max = φ 2 α λ E max , w 0 λ E max Π B D I - 1 = φ 2 D I - 1 , γ D I - 1 α , E max ET E max = φ 3 E max α λ , w 0 α Π B D I - 1 = φ 3 D I - 1 , γ λ , w 0 ET λ w 0 = φ 4 E max λ w 0 , α w 0 Π B γ - 1 = φ 4 D I γ - 1 , γ - 1 E max , w 0 ET E max = φ 5 α w 0 , λ w 0 E max Π B D I - 1 = φ 5 γ - 1 , γ D I - 1 .

These five hydrologic spaces provide different, if related, scaling laws that allow us to focus on the role of specific combinations of parameters and emphasize different hydrologic conditions. It is useful to note that the dimensionless group ϑ=γDI-1=λw0EETmax appears often as an independent variable. The latter can also be written as ϑ=λη, which uncovers an interpretation of it as a ratio of timescales, one related to the mean time between rainfall occurrence, 1/λ, and one related to the time of depletion of the soil store reservoir of depth w0 at the maximum evapotranspiration rate, i.e., η=Emax/w0. Note also that ϑ=λη appears naturally in the normalized form of the evolution equation of the probability distributions of soil moisture dynamics in minimalist stochastic models (Porporato et al.2004; Rodríguez-Iturbe and Porporato2005). Figure 2b reports different partitioning curves for different values of the storage index γ, obtained by solving a minimalist stochastic soil moisture model with only four dimensional parameters as in Eq. (17).

Figure 3Transient trajectories of the time-varying ratio ETt/Rt as a function of the time-varying dryness index Dt and time of year for Mediterranean and tropical dry climates. Brackets indicate ensemble averages for a given time of year, taken over the rainfall variability due to a time-dependent (i.e., seasonal) marked Poisson process of rainfall and seasonal potential evapotranspiration. Each loop is derived from the same climate inputs, with the exception of the phase difference between rainfall and potential evapotranspiration, with the thick grey lines showing results for out-of-phase and thick black lines for in-phase. The +1 on the time-of-year axis corresponds to when maximum rainfall occurs (wet season), and −1 corresponds to timing of minimum rainfall (dry season). The annual average value for each loop falls on a single point on the annual Budyko curves, shown as dashed lines (black and grey) for those that account for climate seasonality and as thin black lines for the classical curve which considers only annually averaged climate values. After Feng et al. (2015).

With the goal of analyzing the suitability of different hydrologic spaces for capturing the information available in global datasets on hydrologic partitioning, Daly et al. (2019) analyzed a physical law focused on hydrologic fluxes (i.e., rates) only. Accordingly, they considered an extension of Budyko and Turc physical laws (Eq. 9) by hypothesizing an additional control by a governing variable ϕ with the dimensions of a flux, [ϕ]=LT-1=Φ,

(20) ET = f ET R , ET max , ϕ .

Since each of the fluxes (L/T=Φ) can be used as a repeating variable, three hydrologic spaces are obtained. Assuming further that ϕ is a flux related to the storage capacity and the frequency of rainfall, ϕ=λw0, the three spaces (Daly et al.2019) are found to coincide with specific cases in Eq. (19): case 1 if R is chosen as a repeating variable, case 2 if ETmax is chosen, and case 4 if ϕ is chosen.

The analysis by Daly et al. (2019) shows that by accounting for the ability of a catchment to store water to supply evapotranspiration, the flux ϕ (referred to as the maximum storage rate) serves as a modulator for the relationships of ET with R and ETmax for very dry and very wet catchments. In this way, accounting for it, it allows us to expand the Budyko and Turc frameworks, suggesting that they are not equivalent, as often assumed in parametric models, unless ϕ→∞. Thus the variable ϕ helps group catchments with similar evapotranspiration rates with respect to different combinations of the dimensionless groups ϕ/ETmax=(λw0)/ETmax=ϑ=γDI-1 and ϕ/R=w0/λ=γ, which account for key catchment and hydrologic characteristics. It also facilitates the analysis of the partitioning in catchments with intermediate values of the dryness index, while Budyko's and Turc's hypotheses help in the analysis of very dry and wet catchments, respectively (Daly et al.2019).

3.1.3 Seasonality, variable coevolution, and higher-order effects

It is logical to wonder about the effects of adding other potentially important variables to the physical law of rainfall partitioning. These could include variables describing seasonality, soil and vegetation properties, and other details of the climatic forcing. With the goal of investigating the role of seasonality in rainfall and evapotranspiration, Feng et al. (2012) considered the duration of the wet season and the intensity of rainfall seasonality. The effect of the new dimensionless groups obtained on the long-term partitioning is shown in Fig. 2b. Depending on the degree and type of seasonality, a reduction of evapotranspiration is typically observed due to an increase in percolation during the wet season compared to the homogeneous case with no seasonality. At the level of approximation afforded by the long-term averaging, the effects of seasonality are hardly distinguishable from those of a decrease in the storage index γ. To disentangle the various effects and more clearly see the effects of seasonality, the temporal variability in the rainfall partitioning must be considered. An example of this is presented in Fig. 3, where the Budyko curve is plotted parametrically as a function of time and the non-uniqueness related to seasonal storage becomes evident (Feng et al.2015).

It is clear that more detailed analyses to disentangle the role of different ecohydrological variables in the rainfall partitioning should focus on specific aspects of the space–time variability of evapotranspiration, which instead are lost in the spatially lumped, long-term evapotranspiration rates of Budyko's type analyses. The combination with simple physically based models may be a valuable way to sharpen the hypotheses of the physical laws, especially when trying to unravel the effects of the covariation of some of the variables. These regard questions of how the rooting depth, and thus the storage capacity of the active soil layer, may depend on the dryness index and whether such covariations may imply some adaptation of vegetation and soil properties to the hydroclimatic characteristics. Along these lines, the minimalist stochastic model developed by Porporato et al. (2004) pointed to a tendency, for the points lying on the Budyko curve, to converge towards typical values of the storage index of 5–6, in turn suggesting a root-depth adaptation to the intensity of rainfall. Similarly, Or and Lehmann (2019) attributed the convergence of the rainfall partitioning to typical evaporative depths, while Hunt (2021) connected a dependence of the storage index with hydroclimatic characteristics, γ(DI), to an adaptation brought about by constraints linked to percolation statistics.

Shedding light on the role of higher-order controls in rainfall partitioning requires charting of more detailed model–data investigations, branching off the beaten path of Budyko-type analyses. Long-term averages may only contain a weak signal of those interactions, which may be easily overwhelmed by noise and other data limitations. This suggests a need for adapting dimensional analysis to these other aspects of the rainfall partitioning and formulating more sophisticated physical laws that capture the more subtle controls by the soil–plant–atmosphere system (see Feng et al.2018, for an example).

3.2 Soil formation, weathering, and loss

A partitioning of the soil mineral component also takes place, on much larger timescales, on the land surface (Riebe et al.2004). The loss of minerals at the land surface is controlled by a series of tightly connected chemico-physical processes with interesting geomorphology and hydrologic interactions (Maher and Chamberlain2014), which in turn have a crucial impact on the surface energy and water balances as well as several biogeochemical, ecological, and climate processes (Garrels1983; Richter and Billings2015).

The input of solid material to the soil from breaking up parent rock is called rock denudation. While the dissolved minerals and the very fine particles resulting from denudation can be transported away by water, the remaining particles stay in the ground, where they continue to be weathered until they too can be transported away. As the soil ages and transforms chemically, the soil may also be deformed by the action of internal stresses and move as creeping flow. As a result, if one assumes, for simplicity, that the input and output balance (Riebe et al.2004), denudation D equals soil formation rate. These in turn balance the losses by weathering W and erosion plus creep, which are lumped here into a single term E (Fig. 1), yielding the following soil partitioning:

(21) D = W + E .

To analyze this partitioning and help resolve the complex interaction between chemical weathering, climate, and the hydrologic cycle, Calabrese and Porporato (2020) used dimensional analysis. This allowed them to obtain a theoretical framework to organize the existing data on weathering rates. Because of the crucial role of leaching of dissolved inorganic carbon on the sites of weathering, the problem strongly depends on the rainfall partitioning at the surface discussed in the previous section. Accordingly, focusing on the weathering rate as the governed quantity, Calabrese and Porporato (2020) wrote a physical law assuming that weathering rates are a function of the input given by the denudation rate, a maximum weathering rate (which includes the effects of type of parent material, temperature, etc.), the concentration of dissolved inorganic carbon [DIC], and the percolation flux LQ (see Eq. 3). Further assuming that both [DIC] and LQ only depend on the surface rainfall partitioning through the dryness index, as in Budyko's hypothesis (Eq. 12), they wrote

(22) W = f W D , W max , D I .

Formally, apart from the presence of the dimensionless dryness index, the situation is similar to the one of Eq. (3), since W, D, and Wmax all have the same dimensions of a flux of mineral per time. Choosing the input D as the repeating variable in Eq. (21), as Budyko did in the rainfall partitioning, the Π theorem gives

(23) W D = φ W W max D , D I .

To verify this functional relationship, Calabrese and Porporato (2020) used literature data including granitic, basaltic, and shale terrains and encompassing a broad range of environmental settings. After normalizing the data by their maximum rate and assuming for the range of data available that WD, the results (Fig. 4a) show a remarkable linear trend of the type

(24) W D = W max D m ψ W D I .

Since the exponent m of the power law is very close to 1, the expression can be further simplified to the final form

(25) W W max = ψ W D I .

Figure 4(a) Analysis of the self-similarity hypothesis in Eq. (24); (b) specific weathering rates as a function of the dryness index; dashed line shows Eq. (26). After Calabrese and Porporato (2020).

Both these expressions are suggestive of self-similarity. At a first glance, looking at Eq. (24), one would be led to conclude that we are in the presence of self-similarity of the second type with respect to Wmax/D. However, an exponent practically equal to one, unlike the typical irrational numbers expected in incomplete self-similarity (Barenblatt1996), leads to a simplification which eliminates the governing variable D. As a result, the final expression (Eq. 25) can be interpreted as a case of complete self-similarity with respect to D/Wmax, when Wmax is chosen as a repeating variable in Eq. (22) instead of D. Alternatively, the same expression can be obtained directly from the Π theorem starting from an abbreviated physical law W=fWa(Wmax,DI). The precise meaning of this self-similarity remains an open problem to be investigated. The physical message however is that the denudation rate is not a relevant quantity compared to the maximum weathering.

Equation (25) provides a mathematical structure to analyze how the water cycle affects the chemical depletion fraction through the dryness index. The empirical data fit provides an equation with only one parameter,

(26) ψ W D I = 1 - ln D I α + 1 1 + ln D I α + 1 ,

which emphasizes the strongly nonlinear relation between water availability and weathering rate (see Fig. 4b). Starting from high DI values, the normalized weathering rate increases only slightly up until DI∼2, after which it steeply increases before plateauing at DI∼0.5. The increase region corresponds to the establishment of grassland, savanna, and shrubland ecosystems, emphasizing the important role of vegetation in acidifying the soil and in turn enhancing weathering rates. A similar nonlinear dependence on wetness has been observed in the global distribution of soil organic carbon (Kramer and Chadwick2018).

The scatter shown by the data in Fig. 4b is not unexpected, because of the large datasets used, which are obtained from different locations and different methods. While it is likely that part of it may be attributable to noise, it would also be interesting to consider whether adding other quantities in the physical law might help explain more of this variability. A logical starting point could be the analysis of the effects of the soil depth and water-storage capacity, which could be explored by correlating the data with estimates of the storage index γ. It is likely that this will reveal further linkages between surface hydrology, soil formation, and weathering.

3.3 Landscape self-similarity

Besides the macroscopic effects of wetness on weathering discussed in the previous section, the loss of minerals driven by the coupled water and sediment fluxes on the topographic surface is also responsible for the formation of complex topographic patterns, which in turn impact ecohydrologic processes (Dietrich and Perron2006). In this section, we focus on the dendritic morphology of interlocked ridges and valleys, which often displays a self-similar scaling linked to the fractal behavior of landscapes (Rodríguez-Iturbe and Rinaldo2001). While on the one hand the complexity of these patterns hinders analytical results, on the other hand it also helps with the investigation of macroscopic (i.e., averaged) behavior, giving rise, as we will see, to emergent scaling behaviors for the large-scale statistics of the sediment budget, the mean elevation profile, and the landscape spectral properties.

Figure 5Steady-state simulation of landscape elevation given by Eqs. (27) and (28) for m=0.5 and different channelization indices, Eq. (33), 𝒞I=103 in (a) and 𝒞I=104 in (b). The elevation along A–A and B–B transects is shown in (c) and (d). After Hooshyar et al. (2021).


Towards the goal of performing dimensional and self-similarity analyses, here we will consider only (spatially) averaged quantities in idealized geometries for which the solutions may be expected to depend only on a limited number of dimensionless quantities. The effort of formulating meaningful physical laws for average landscape quantities is facilitated by the existence of simplified, semi-empirical landscape evolution models (LEMs) (Chen et al.2014), which can be inspected to infer the main governing quantities of the problems. We specifically refer to a minimalist landscape evolution model in detachment-limited conditions (Howard1994), although many considerations also apply to transport-limited and other intermediate formulations (Davy and Lague2009; Pelletier2012). Accordingly, the equation for the evolution of the landscape elevation z(x,y,t) is (Chen et al.2014)

(27) z t = δ 2 z - K a m | z | + U ,

where t is time, δ is the diffusion coefficient used to represent the intensity of soil creep, K is the fluvial erosion coefficient, a(x,y,t) is the specific drainage area, m is a dimensionless coefficient (for simplicity we assume the exponent of the gradient in the erosion term to be equal to 1), and U is the surface-growth term due to tectonic uplift. Equation (27) is coupled to the “conservation” equation for the specific drainage area:

(28) a z | z | = 1 .

The latter was derived by Bonetti et al. (2018) as an idealized representation of the water flow, following the steepest descent lines of a topographic surface with a given characteristic speed (see also Bonetti et al.2020). Mathematically, a=limw0A/w, where A is the contributing area and w is a finite portion of a contour (iso-elevation); hence, the specific drainage area has dimensions of length, [a]=L, and is defined at a point on the landscape, whereas the contributing area A, [A]=L2, is actually zero when considered at a point, unless the topographic surface is discontinuous (see Bonetti et al.2018, for details). The use of the nonlocal quantity A instead of a in Eq. (27) is therefore incorrect and leads to grid-dependent results in mathematical codes. Since the equation of a is not very well known, and given the misuse of A in the literature, we present a short derivation of Eq. (28) in Appendix B.

The coupled Eqs. (27) and (28) form a closed system once the initial and boundary conditions are specified. Figure 5 shows simulation results producing complex patterns (Anand et al.2020) that resemble real landscapes, with characteristic branching and channelization instability. A dimensional analysis of terms in a LEM in detachment-limited conditions was presented by Theodoratos et al. (2018), although they use A instead of a in the evolution equation for h. Throughout this section, we will also see how the statistics of landscape elevation bear a striking resemblance to those of streamwise velocity in near-wall turbulence (Hooshyar et al.2020). For comparison, Appendix C briefly reports the Navier–Stokes equations for the case of channel flow between two parallel plates, which corresponds to the case of the landscape evolution between two parallel boundaries at fixed elevation (e.g., Fig. 5).

3.3.1 Sediment budget and soil partitioning

In steady state, spatial averaging leads to the partitioning of soil material into soil creep (SC) and the soil loss by stream erosion (SL),

(29) U = CR + SL .

This equation is obviously related to Eq. (21), since at steady state the rate of uplift needs to equal the denudation rate, while the losses may be grouped by focusing on different features of the underlying processes. Equation (21) emphasizes chemical weathering while lumping all the other losses (mainly erosion and creep) into the term E, while here Eq. (29) draws attention to the water-erosion term SL (possibly including also weathering), with the rest lumped into a diffusive term, CR, representing soil creep and other smoothing processes.

Averaging spatially Eqs. (27) and (28), the mean sediment balance equation can be written as (see Hooshyar and Porporato2021, for details)

(30) δ 2 l y - 1 z Ω - CR - K ( 1 - m ) κ c a m z + m a m - 1 z SL + U = 0 ,

where ly is the domain width, z is the component of the elevation gradient normal to the domain contour Ω, κc is the plane or contour curvature, the brackets 〈⋅〉 indicate the spatial average, and 〈⋅〉Ω is the average over the domain contour. While this is an exact result, the terms in brackets are unknown because of the complexity of the landscape surface. This is where dimensional analysis comes in handy.

Bonetti et al. (2020) tackled the problem of these spatial averages by adopting a directional dimensional analysis (see Sect. 2.3), considering independent dimensions for horizontal lengths, L, and the vertical length, Lz. Physically, a justification for the use of directional dimensional analysis may be found in the fact that Eq. (27) is the sediment budget equation at a point written in terms of volume of sediments per unit ground area. For sediments of constant density, the same equation can be written in terms of mass per unit ground area, say ρz, which then can be considered a new variable with its independent dimensions. Indicating the unknown average slope at the boundary as S*=zΩ, the physical law can be written as

(31) S * = f l y , K , U ; δ , m ,

where the various quantities are suggested by inspection of the governing equations and the type of boundary conditions. For this problem, [S*]=LzL-1, [U]=LzT-1, [δ]=L2T-1, [K]=L1-mT-1, and [m]=1. Choosing ly, K, and U as repeating variables, the Π theorem then yields

(32) δ S * U l y = 2 CR U = φ CR C I - 1 , m = ψ CR C I , m ,


(33) C I = K l m + 1 δ

is the global channelization index. Figure 6 shows the plot of the function (Eq. 32) obtained from a set of simulations with different CI and m values.

These results suggest a tantalizing analogy with the analysis of wall-bounded turbulent shear flows. In fact, the regular behavior of φ(CI,m) reminds one of the behavior of the Darcy friction factor plotted in the Moody diagram for pipe flow (Munson et al.2013). For a detailed description of dimensional analysis in turbulence, see Barenblatt (1996) and Katul et al. (2019). In this analogy, the slope of the elevation profile and the slope of the velocity profile play a similar role. As is well known, in turbulence, the mean velocity profile at the wall is proportional to the wall shear stress and relates to the partitioning of viscous and turbulent stresses; this is analogous here to the partitioning of soil coming from uplift into creep and stream erosion – see Eq. (30). We will highlight further links with the dimensional analysis of wall-bounded turbulent flows in the subsequent developments.

It is also useful to note that in Bonetti et al. (2020) the Π theorem was applied to Eq. (31) considering ly, U, and δ to be repeating variables, which instead leads to Πδ=φδ(CI), with ΠδCRCI and φδ(CI)=ψCR(CI)CI-1. In the turbulence analogy, the present analysis corresponds to choosing viscosity instead of density (which is the typical choice; see Katul et al.2019) as one of the repeating variables in the physical law for the wall-shear stress. Thus, for the wall-shear stress, τ=μΣ*, where μ is the dynamic viscosity and Σ* is the slope of the streamwise velocity profile at the wall, the physical law is τ=fτ(V,L,ρ,μ,ϵ), where ρ is the density, V the mean velocity, L the characteristic lateral dimension, and ϵ the roughness height. Using V, L, and μ as repeating variables, one has τLμV=ReτρV2=RefDarcyRe,ϵL, where Re=ρVLμ is the bulk Reynolds number.

Figure 6Average sediment partitioning according to Eq. (32) as a function of the channelization index for different values of m. Points refer to numerical simulations and solid line refers to the analytical solution for the unchannelized regime (courtesy of Sara Bonetti).

A second interesting fact is the plateauing of the curves at large values of CI for both the analytical solution in the unchannelized case and the numerical results in the channelized regime. This suggests the existence of complete self-similarity for CI→∞, again analogous to the self-similarity observed in the Moody diagram for the friction factor in the fully rough regime (Munson et al.2013). Thus, while erosion dominates over creep, the effect of diffusion does not fully disappear but remains present, probably concentrated in an area of zero measure corresponding to the network of ridges and valleys. Again, this singular limit bears similarities to the much investigated limit of the viscous Navier–Stokes equations for Reynolds numbers tending to infinity as well as the hypothesis of dissipative anomalies in the inviscid Euler equations (e.g., Eyink2008).

3.3.2 Mean elevation profile

Hooshyar et al. (2020) analyzed the profile of the mean elevation, z(y)=limlx1lx0lxz(x,y)dx, of channelized landscapes. Based on the inspection of the governing equations and again using directional dimensional analysis, the physical law for the slope of the profile was assumed to be

(34) d z d y = f d z d y y , δ , z * ; l y , K , U , m ,

where y is the distance from the boundary and z* is an elevation scale. Choosing y, δ, and z* as repeating variables and with simple manipulation of the dimensionless groups, the Π theorem gives

(35) ( m + 1 ) η d φ d η = f 3 η , C I , ζ , m ,

where φ=zz*. In addition to the global channelization index, CI, reflecting the relative impact of fluvial erosion to diffusive transport, η=Kym+1δ is a local variable with a similar form to 𝒞I but capturing the local relative contribution of those two processes, while ζ=Uly2δz* describes the relative impact of tectonic uplift on diffusive transport.

In a system with relatively small diffusive transport and dominated by erosion and uplift, 𝒞I and ζ take high values. The same argument also applies to η except for locations at an intermediate distance from the boundary. Thus, when the variables η, 𝒞I, and ζ reach such an asymptotic condition, one may assume complete self-similarity (Barenblatt1996) according to which the function f3 is independent of these quantities

(36) η d φ d η = κ ( m ) ,

where κ is only a function of m. Integrating Eq. (36) yields

(37) φ = κ ( m ) ln η + C ,

where C is independent of η but may still depend on m, 𝒞I, and ζ. Equation (37) describes the logarithmic scaling of the mean-elevation profile with respect to η. The emergence of such a logarithmic profile was confirmed in numerical simulations, laboratory experiments, and real landscapes as well as in other models of complex branching such as the optimal channel networks and directed percolation (Hooshyar et al.2020). Figure 7 shows the flattening of mean elevation profiles with increasing CI (recall that a similar effect is also observed in turbulent velocity profiles) along with the values of the coefficient κ of the logarithmic profile.

Figure 7Mean elevation profiles for increasing values of CI (from blue to red lines) and, on the right, coefficient of the logarithmic profile as a function of the exponent m for simulations, the XLE laboratory experiments (Hooshyar et al.2019), and the landscape at the Calhoun Critical Zone Observatory (CCZO). Modified after Hooshyar et al. (2020).

As observed by Hooshyar et al. (2020), the presence of a mean logarithmic profile of elevation at an intermediate distance from the domain boundaries is similar to the classic results for the turbulent velocity profile (see also Barenblatt1996). In wall-bounded turbulence (Katul et al.2019), the logarithmic profile for the mean velocity profile is obtained from the Π theorem applied to the physical law for the mean gradient, dudy=f(y,u*,ρ,L,μ,ϵ), where u*=τρ is the friction velocity (the other quantities have the same meaning as in the previous subsection) and then assuming complete self-similarity with respect to both the bulk and local Reynolds numbers. From a phenomenological point of view, the analogy between the two phenomena is found in the resemblance between progressive penetration to smaller scales of ridges and valleys into the landscape and the intensification of vorticity producing smaller and smaller vortices (i.e., smaller Kolmogorov scales). The increased turbulent mixing and the progressive land-surface dissection with sharper sequences of ridges and valley surface dissection with increasing channelization index and Reynolds number, respectively, produce in both cases a flattening of the mean profile observed in the turbulence mean velocity and the mean elevation (see Fig. 5a) as well as a logarithmic scaling. The fact that a logarithmic region is found in different contexts, including directed percolation (Hooshyar et al.2020), hints at the generality of such self-similar scaling as a robust outcome in dynamically different complex systems, appearing as a dimensional consequence of length-scale independence in spatially bounded complex systems.

3.3.3 Spectral analysis of elevation fluctuations

Our final example follows Hooshyar et al. (2021), who analyzed the power spectral density (PSD) of elevation transects at a given y (see Fig. 5),

(38) P y ( ω ) = 1 l x | z ^ y ( ω ) | 2 ,

for lx→∞, where

(39) z ^ y ( ω ) = 0 l x z y ( x ) e - 2 π i ω d x

is the Fourier transform of zy(x) and ω is the longitudinal wavenumber or inverse scale. These PSDs peak at a wavenumber (i.e., the most energetic mode) which provides a characterization of the typical valley spacing. Beyond such wavenumbers, a power-law scaling is typically visible in simulations, producing an asymptotic behavior which can be collapsed onto a single curve at high 𝒞I.

Hooshyar et al. (2021) analyzed the PSD spectral scaling with the aid of directional dimensional analysis. They argued that for a longitudinal elevation series at a given y, the amount of “energy” (i.e., variance of the elevation fluctuations) at a wavenumber ω must depend on the following variables:

(40) P y ( ω ) = g 1 ω , y , l y , δ , K , U , m ,

where g1 is the physical law. The energy 𝒫y(ω) is defined over the fluctuations along the x axis (see Eq. 38) and has dimension Lz2L. The wavenumber along the x direction has dimension [ω]=L-1, while δ, K, and U have dimensions L2T−1, L1-mT-1, and LzT−1, respectively; y and ly have dimension L.

Given three dimensions Lz, L, and T and seven dimensional governing variables and choosing K, U, and ω as repeating variables, the Π theorem yields

(41) P y ( ω ) K 2 ω 3 - 2 m U 2 = g 2 K ω - ( m + 1 ) δ , ω l y , ω y , m .

A manipulation of Eq. (41) leads to

(42) P y ( ω ) K 2 ω 3 - 2 m U 2 = g 3 C I , η , η ω , m .

The quantity η=Kym+1/δ has the same form as that of 𝒞I but defined locally at y distance from the boundary, and ηω=Kω-(m+1)/δ is equivalent to η but defined in the frequency domain.

In the asymptotic limit of relatively high 𝒞I, η and ηω attain a near-constant limit away from the boundary, implying complete self-similarity,

(43) P y ( ω ) ω 2 m - 3 ,

where the proportionality coefficient is (U/K)2g3(m). Equation (43) predicts an exponent of the power spectral density which is independent of δ and has a power-law decay. The condition of high 𝒞I, η, and ηω is expected in systems that are dominated by erosion (high 𝒞I), far enough from the boundary (high η), and within small enough scales (high ηω). Again, a parallel with the Kolmogorov spectrum of turbulence (Pope2000; Katul et al.2019) can be drawn with the channelization index playing the role of the Reynolds number. Moreover, at a sufficient distance from the domain boundary, it may be assumed that the information regarding the domain geometry and direction is lost, becoming statistically isotropic. This is similar to the local isotropy at small scales (or eddies detached from the boundary) of fully developed turbulent flow.

Figure 8a shows the exponent of the power fits to PSDs for simulations with 𝒞I≥105, denoted by α, for different m values. This finding agrees with the relation α=2m-3 in the intermediate range of m and supports the validity of the assumption of complete self-similarity with respect to ηω. The inset of Fig. 8b also shows the function g3(m) from numerical simulations. The spectral scaling was also confirmed in the laboratory experiments and the real topography at the Calhoun Critical Zone Observatory (Fig. 8c–f).

Figure 8(a) PSDs of elevation longitudinal series for different values of m at an intermediate distance from the boundary. (b) Slope of the power fit to the declining part of the PSD from simulations with 𝒞I≥105, denoted by α, as a function of the exponent m. The data from the XLE physical experiment are also shown. The black line is the relation α=2m-3 derived from dimensional and self-similarity arguments in Eq. (43). The inset shows the function g3(m) from numerical simulations. Panel (c) shows an example of the XLE experimental landscape (Hooshyar et al.2019), of which panel (d) shows the PSD. (e) Portion of the landscape at the Calhoun Critical Zone Observatory in South Carolina, USA, from 1 m resolution lidar data. (f) The PSD is computed from one-dimensional transects within the area contained by the two parallel lines that are located at an intermediate distance from the main channel showing two distinct power-law scaling regions with similar exponents. After Hooshyar et al. (2021).

From a geomorphological point of view, the connection between the exponent α and the parameter m in the erosion term provides a useful link to landscape processes, since steep landscapes with debris‐flow‐dominated channels have been associated with smaller m, while flatter fluvial landscapes are characterized by larger values of m (Montgomery and Foufoula-Georgiou1993). The values of the spectral exponent are also interesting from a more theoretical point of view in relation to the role of nonlinearities as a function of scale. Our analysis has shown consistently values of α-2.5, while several studies previously reported exponents near −2 (Newman and Turcotte1990; Huang and Turcotte1989; Passalacqua et al.2006); the latter corresponds to a fractal dimension Dm=1.5 (Huang and Turcotte1989; Voss1985) and to the presence of an underlying fractional Brownian noise (Turcotte1987; Bell1975) with a Lorentzian spectrum and an exponential decay of the elevation autocorrelation. These models with α=-2 would therefore imply linear stochastic dynamics at odds with the known presence of nonlinear terms responsible for the very formation of the channel network. On the one hand, α=-2.5 would mean an erosion exponent m=0.25, therefore preserving the nonlinearity of the dynamics also for macroscopic scaling relationships, like the PSD scaling. On the other hand, the fact that the observed value is not far from α=-2 might mean that averages taken over complex landscape patterns cause an effective reduction of nonlinearity, whereby the activation of many degrees of freedom at high channelization regimes causes a statistical regularity which muffles the small-scale nonlinearities. This is certainly an interesting topic that deserves further investigation.

Finally, from a practical point of view, the spectral scaling of landscape elevation could be profitably utilized in developing efficient numerical simulations of landscape evolution (Passalacqua et al.2006). Such numerical schemes would potentially resemble large eddy simulation methods used in fluid turbulence (Pope2000), where the unsolved dynamics at finer scales are approximated by extrapolating the PSD. The improved speed would be a great asset for large-scale simulations of landscape evolution under different scenarios.

4 Conclusions

We have presented several examples of applications of the Π theorem and self-similar scaling to the partitioning of water and sediments driven by the terrestrial water cycle. It is time to draw to a close and ask ourselves whether by using dimensional analysis we have learned anything useful regarding these problems. An answer in the affirmative is suggested by considering that dimensional analysis helped reveal the dominant controls of the dryness index and storage index in the long-term rainfall partitioning, while in the weathering analysis it allowed us to extract the key nonlinear control of the dryness index, making order in the vast amount of information contained in global datasets from different experiments and environmental conditions. Finally, the analysis of the complex geometries obtained from landscape evolution models allowed us to discover new macroscopic relations among macroscopic variables in the partitioning of sediments, the elevation profile, and the spectral scaling. The analogy between landscape elevation and turbulent velocity fluctuations has also been fruitful and promises further results.

If these observations confirm the utility of dimensional analysis, it should also be clear that these methods are not a foolproof set of rules to achieve miraculous solutions. Rather, they are an iterative procedure to sharpen our hypotheses on physical processes. Thinking of dimensional arguments as a form of modeling allows an “explication of the role abstraction and multiple realisability; not as compatibility with other possible worlds but as compatibility with different fictional descriptions of our own world” (Pexton2014). The emergence of scaling laws is often only asymptotic (Koutsoyiannis et al.2018), as in the self-similarity approximations of Sect. 2.2.1, while in several instances it is either difficult to reduce the number of variables in the starting physical law or the problems are hard to formalize in terms of standard dimensions (e.g., in problems interfacing with the social sciences and economics). These challenges notwithstanding, we hope however that these considerations will help breathe new life into Strahler's view (Strahler1958) that dimensional analysis will become increasingly useful in hydrology, geomorphology, and beyond.

Appendix A: The Rayleigh–Riabouchinsky controversy: thermodynamic limit as complete self-similarity

In this Appendix we discuss in more detail the Rayleigh–Riabouchinsky controversy related to augmented dimensional analysis. In particular, after a brief historic background, we show how the thermodynamic approach of Rayleigh corresponds to a self-similar solution of the first kind with respect to the dimensionless group obtained from the Boltzmann constant, which serves as the dimensional unifier when augmenting the dimensional analysis from mechanics to thermodynamics.

As mentioned in the main text, Rayleigh (1915a)2 considered the problem of the heat flux h, the governed quantity, between a body of characteristic dimension  and a stream of incompressible and inviscid fluid moving with velocity v. Besides  and v, the other governing quantities are the specific heat capacity c and the thermal conductivity κ. Rayleigh formulated a physical law with n=5 governing variables and adopted length, time, energy, and temperature as primary dimensions (i.e., kRay=4), as typical in thermodynamics (of course one could also use LMT plus temperature, using mass instead of energy). As a result, he obtained one (n-kRay=5-4=1) governing dimensionless group, which he used successfully to describe the problem.

In a very short comment published soon after Rayleigh's paper, Riabouchinsky (1915) objected that, adopting the more advanced and detailed point of view of statistical mechanics, according to which temperature is the mean molecular kinetic energy, one could more parsimoniously limit the primary dimensions to length, time, and energy without involving temperature. As a result, he obtained two (n-kRia=5-3=2) dimensionless groups. In his response, Rayleigh rejected Riabouchinsky's alternative, commenting that

It would indeed be a paradox if the further knowledge of the nature of heat afforded by molecular theory put us in a worse position than before in dealing with a particular problem. The solution would seem to be that the Fourier equations embody something as to the nature of heat and temperature which is ignored in the alternative argument (Rayleigh1915b).

Rayleigh's reply was authoritative and sensible but not completely satisfactory. As pointed out by Buckingham (1915): “since he does not pursue the subject further and the reader may feel as if left in mid-air, it seems worth while that the point raised by M. Riabouchinsky should be somewhat further elucidated.” Since then the debate has been the subject of continued discussions. Here, taking a cue from previous analyses of this controversy (e.g., Gibbings1980; Butterfield2001) and adopting the use of a dimensional unifier (Panton2006), we show that Rayleigh's result may be seen as a complete self-similarity solution of the augmented physical law obtained by including the Boltzmann constant as a dimensional unifier to link thermodynamics to mechanics. Starting from a thermodynamic approach to the problem is tantamount to taking for granted the existence of this limit.

Thus, following Rayleigh, we specifically choose length (L), time (T), energy (heat, H), and temperature (Θ) as primary dimensions. With [h]=HT-1, [ℓ]=L, [θ]=Θ, [c]=HL-3Θ-1, [κ]=HL-1T-1Θ-1, [v]=LT-1, and [kB]=HΘ-1, we obtain the dimension matrix

(A1) h θ c κ v k B L 0 1 0 - 3 - 1 1 0 T - 1 0 0 0 - 1 - 1 0 H 1 0 0 1 1 0 1 Θ 0 0 1 - 1 - 1 0 - 1

of rank 4. As a result, the Π theorem gives 7-4=3 dimensionless groups.

Choosing the heat flux h as the governed quantity and selecting , θ, c, and κ as dimensionally independent (repeating) variables among the governing quantities leads to the physical law

(A2) h = f , θ , c , κ ; v , k B ,

from which the Π theorem then gives

(A3) h d θ κ = φ v c κ , k B 3 c .

With a system of units for which kB=1 and inverting the second Π group, one obtains the result of Riabouchinsky (1915),

(A4) h θ κ = φ Ria 1 3 c .

This however would imply the use of very small units, because in the usual SI system, kB=1.380649×1023 J K−1. As a result, apart from systems at the nanoscale, in normal conditions kB/(3c)0. Our everyday experience, on which thermodynamic concepts are based, shows that we can neglect this term in Eq. (A3) to arrive at the result, originally obtained by Rayleigh,

(A5) h θ κ = φ Ray v c κ .

Thus, looking at Eq. (7), it becomes clear that the thermodynamic limit, which allows us to go from statistical mechanics to continuum mechanics and thermodynamics, corresponds to a complete self-similar solution of Eq. (A3) with respect to kB/(3c). It should be stressed that this is not a mathematical fact based solely on the smallness of the Π group, but an empirical ascertainment, as for other self-similar asymptotic regimes (Barenblatt1996). Alternatively, the same result is obtained by straightforward application of the Π theorem to the restricted physical law

(A6) h = f Ray ( , θ , c , κ ; v ) ,

as originally done by Rayleigh.

Appendix B: Derivation of Eq. (28)

We follow here Bonetti et al. (2020) in deriving the equations for the minimalist LEM. The time evolution of the surface elevation z(x,y,t) is described by the sediment continuity equation

(B1) z t = U + f ,

where t is time, U is the uplift rate, and f is the total volumetric sediment flux, given by the sum of fluxes related to runoff erosion (fc) and soil creep processes (fd). The latter is assumed to be proportional to the topographic gradient, fd=δz, δ being a diffusion coefficient. In the so-called detachment-limited conditions all eroded material is considered transported outside the model domain. Thus, with no sediment redeposition, the runoff erosion becomes a sink term, fcKaqm|z|n, where Ka is an erosion coefficient, q is the specific discharge (per unit width of contour line), and m and n are model parameters. As a result, Eq. (B1) becomes

(B2) z t = U - K a q m | z | n + δ Δ z .

The surface water flux q is linked to the continuity equation

(B3) h t = Q + ( q n ) = Q + ( h v n ) ,

where h is the water height, n is the direction of the flow, Q is the part of the rainfall rate that has not infiltrated and that effectively contributes to runoff production, and v is the water velocity. Equation (B3) can be simplified assuming steady-state conditions with constant, representative runoff rate, QQ0, and speed of water flow v0 in the direction opposite to the landscape gradient (i.e., n=-z/|z|). In such conditions, with the specific discharge being proportional to the water height, q=v0h, Eq. (B3) becomes

(B4) h t = 0 = Q 0 + v 0 h z | z | .

Dividing it by Q0, Eq. (B4) becomes the equation for the specific catchment area, Eq. (28), with a=v0h/Q0. Finally, since q=aQ0, one can write Kaq=Kaam, where Ka=KaQ0m, thus obtaining Eq. (27).

Appendix C: Analogy with the velocity profile in turbulent flows

Here we report the equation for the streamwise velocity u for the flow between parallel plates to better understand the analogy between the statistics of near-wall turbulence and those of landscape evolution, discussed in Sect. 3.3. We begin with the Navier–Stokes equations for incompressible fluids (Munson et al.2013),

(C1) ρ D v D t = - p + μ Δ v


(C2) v = 0 .

For the flow between infinite parallel plates, the previous equations give, for the streamwise component u,

(C3) u t = p x - u u x - v u y - w u z + μ Δ u .

Indicating with G the constant pressure gradient and with q={v,w} the normal velocity

(C4) u t = G - q u + μ Δ u

where the Laplacian is only in y and z, the equation becomes similar to the LEM equation (Eq. B2) with m=n=1. Since here the continuity equation is

(C5) v y + w z = 0 ,

we can write q in terms of a single streamfunction, which automatically satisfies it, i.e.,

(C6) v = ψ z and w = - ψ y .

We thus have a direct correspondence between the streamwise velocity and the landscape elevation, uz, and an indirect one (this is where the equations and the physics obviously differ) between the streamfunction and the specific contributing area ψa. It is also interesting to note that, in these conditions of steady state and x independence, the equations for the other two velocity components and thus also for the streamfunction become decoupled from the equation of u, which means that, unlike the elevation profile, the complex velocity profile in turbulence is shaped by nonstationary terms and streamwise inhomogeneities.

Data availability

No data sets were used in this article.

Competing interests

The contact author has declared that there are no competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


Much of this work has been the result of invaluable collaborations with Edoardo Daly, Jun Yin, Xue Feng, Salvatore Calabrese, Sara Bonetti, Milad Hooshyar, and Shashank Anand. The author is also grateful to Paolo D'Odorico, Gaby Katul, Luca Ridolfi, and Ignacio Rodriguez-Iturbe for continued friendship, support, and advice. We thank Stefan Hergarten, Demetris Koutsoyiannis, Lamberto Rondoni, and two anonymous reviewers for very valuable suggestions.

Financial support

This research has been supported by the US National Science Foundation (NSF) grant nos. EAR-1331846 and EAR-1338694, the BP through the Carbon Mitigation Initiative (CMI) at Princeton University, and the Moore Foundation.

Review statement

This paper was edited by Nunzio Romano and reviewed by Demetris Koutsoyiannis, Stefan Hergarten, and two anonymous referees.


Anand, S. K., Hooshyar, M., and Porporato, A.: Linear layout of multiple flow-direction networks for landscape-evolution simulations, Environ. Model. Softw., 133, 104804,, 2020. a

Aronson, D. G. and Graveleau, J.: A selfsimilar solution to the focusing problem for the porous medium equation, Eur. J. Appl. Math., 4, 65–81, 1993. a

Barenblatt, G. I.: Scaling, self-similarity, and intermediate asymptotics: dimensional analysis and intermediate asymptotics, 14, Cambridge University Press,, 1996. a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u

Barenblatt, G. I., Chorin, A. J., and Prostokishin, V. M.: Scaling laws for fully developed turbulent flow in pipes: Discussion of experimental data, P. Natl. Acad. Sci. USA, 94, 773–776,, 1997. a

Barr, D. I.: Consolidation of basics of dimensional analysis, J. Eng. Mech., 110, 1357–1376, 1984. a

Bell Jr., T.: Statistical features of sea-floor topography, in: Deep Sea Research and Oceanographic Abstracts, vol. 22, Elsevier, 883–892,, 1975. a

Bhaskar, R. and Nigam, A.: Qualitative physics using dimensional analysis, Artific. Intel., 45, 73–111, 1990. a

Bluman, G. W. and Cole, J. D.: Similarity methods for differential equations, in: vol. 13, Springer Science & Business Media, ISBN 978-0-387-90107-7, 2012. a

Bolster, D., Hershberger, R. E., and Donnelly, R. J.: Dynamic similarity, the dimensionless science, Phys. Today, 64, 42–47, 2011. a

Bonetti, S., Bragg, A., and Porporato, A.: On the theory of drainage area for regular and non-regular points, P. Roy. Soc. A, 474, 20170693,, 2018. a, b

Bonetti, S., Hooshyar, M., Camporeale, C., and Porporato, A.: Channelization cascade in landscape evolution, P. Natl. Acad. Sci. USA, 117, 1375–1382, 2020. a, b, c, d, e

Bridgman, P. W.: Dimensional analysis, Yale University Press, ISBN 10:0548597286, 1922. a, b

Buckingham, E.: The principle of similitude, Nature, 96, 396–397, 1915. a

Burton, J. and Taborek, P.: Two-dimensional inviscid pinch-off: An example of self-similarity of the second kind, Phys. Fluids, 19, 102109,, 2007. a

Butterfield, R.: Dimensional analysis revisited, P. Inst. Mech. Eng. Pt. C, 215, 1365–1375, 2001. a, b

Calabrese, S. and Porporato, A.: Wetness controls on global chemical weathering, Environ. Res. Commun., 2, 085005,, 2020. a, b, c, d

Chen, A., Darbon, J., and Morel, J.-M.: Landscape evolution models: A review of their fundamental equations, Geomorphology, 219, 68–86, 2014. a, b

Daly, E. and Porporato, A.: A note on groundwater flow along a hillslope, Water Resour. Res., 40, W01601,, 2004a. a

Daly, E. and Porporato, A.: Similarity solutions of nonlinear diffusion problems related to mathematical hydraulics and the Fokker-Planck equation, Phys. Rev. E, 70, 056303,, 2004b. a, b, c

Daly, E., Calabrese, S., Yin, J., and Porporato, A.: Hydrological Spaces of Long-Term Catchment Water Balance, Water Resour. Res., 55, 10747–10764, 2019. a, b, c, d, e, f

Davy, P. and Lague, D.: Fluvial erosion/transport equation of landscape evolution models revisited, J. Geophys. Res.-Earth, 114, F03007,, 2009. a

Dietrich, W. E. and Perron, J. T.: The search for a topographic signature of life, Nature, 439, 411–418, 2006. a

Dimitrakopoulos, E. G. and DeJong, M. J.: Revisiting the rocking block: closed-form solutions and similarity laws, P. Roy. Soc. A, 468, 2294–2318, 2012. a

Dooge, J. C.: Sensitivity of runoff to climate change: A Hortonian approach, B. Am. Meteorol. Soc., 73, 2013–2024, 1992. a

Efthimiou, C. J. and Llewellyn, R. A.: Cinema, Fermi problems and general education, Phys. Educ., 42, 253–261, 2007. a

Eggers, J. and Fontelos, M. A.: The role of self-similarity in singularities of partial differential equations, Nonlinearity, 22, R1–R44,, 2008. a

Eyink, G. L.: Dissipative anomalies in singular Euler flows, Physica D, 237, 1956–1968, 2008. a

Feng, X., Vico, G., and Porporato, A.: On the effects of seasonality on soil water balance and plant growth, Water Resour. Res., 48, W05543,, 2012. a, b

Feng, X., Porporato, A., and Rodriguez-Iturbe, I.: Stochastic soil water balance under seasonal climates, P. Roy. Soc. A, 471, 20140623,, 2015. a, b

Feng, X., Ackerly, D. D., Dawson, T. E., Manzoni, S., Skelton, R. P., Vico, G., and Thompson, S. E.: The ecohydrological context of drought and classification of plant responses, Ecol. Lett., 21, 1723–1736, 2018. a

Gagnon, J.-S., Lovejoy, S., and Schertzer, D.: Multifractal earth topography, Nonlin. Processes Geophys., 13, 541–570,, 2006. a

Galilei, G.: Dialogues concerning two new sciences, Dover, 1914. a

Garrels, R. M.: The carbonate-silicate geochemical cycle and its effect on atmospheric carbon dioxide over the past 100 million years, Am. J. Sci., 283, 641–683, 1983. a

Gibbings, J.: On dimensional analysis, J. Phys. A, 13, 75–89,, 1980. a, b

Gibbings, J.: A logic of dimensional analysis, J. Phys. A, 15, 1991–2002,, 1982. a

Gilmore, R.: Lie groups, Lie algebras, and some of their applications, Courier Corporation, ISBN 0-486-44529-1, 2012. a

Goldenfeld, N.: Lectures on phase transitions and the renormalization group, CRC Press, ISBN 9780429962042, 2018. a

Gratton, J. and Minotti, F.: Self-similar viscous gravity currents: phase-plane formalism, J. Fluid Mech., 210, 155–182, 1990. a, b

Hankey, A. and Stanley, H. E.: Systematic application of generalized homogeneous functions to static scaling, dynamic scaling, and universality, Phys. Rev. B, 6, 3515,, 1972. a, b

Hills, C. P. and Moffatt, H.: Rotary honing: a variant of the Taylor paint-scraper problem, J. Fluid Mech., 418, 119–135, 2000. a

Hooshyar, M. and Porporato, A.: Mean Dynamics and Elevation-Contributing Area Covariance in Landscape Evolution Models, Water Resour. Res., 57, e2021WR029727,, 2021. a

Hooshyar, M., Singh, A., Wang, D., and Foufoula-Georgiou, E.: Climatic Controls on Landscape Dissection and Network Structure in the Absence of Vegetation, Geophys. Res. Lett., 46, 3216–3224, 2019. a, b

Hooshyar, M., Bonetti, S., Singh, A., Foufoula-Georgiou, E., and Porporato, A.: From turbulence to landscapes: Logarithmic mean profiles in bounded complex systems, Phys. Rev. E, 102, 033107,, 2020. a, b, c, d, e, f, g

Hooshyar, M., Katul, G., and Porporato, A.: Spectral Signature of Landscape Channelization, Geophys. Res. Lett., 48, e2020GL091015,, 2021. a, b, c, d, e

Howard, A. D.: A detachment-limited model of drainage basin evolution, Water Resour. Res., 30, 2261–2285, 1994. a

Huang, J. and Turcotte, D.: Fractal mapping of digitized images: application to the topography of Arizona and comparisons with synthetic images, J. Geophys. Res.-Solid, 94, 7491–7495, 1989. a, b

Hunt, A.: Soil formation, vegetation growth, and water balance: A theory for Budyko, in: Hydrogeology, Chemical Weathering, and Soil Formation, edited by: Hunt, A. G., Egli, M., and Faybishenko, B., American Geophysical Union, 67–80,, 2021. a

Huntley, H. E.: Dimensional analysis, Dover Publications, Open Library OL14610272MLC Control Number 67017978 Library Thing 2164960, 1967. a

Kader, B. and Yaglom, A.: Mean fields and fluctuation moments in unstably stratified turbulent boundary layers, J. Fluid Mech., 212, 637–662, 1990. a, b

Katul, G., Li, D., and Manes, C.: A primer on turbulence in hydrology and hydraulics: The power of dimensional analysis, Wiley Interdisciplin. Rev.: Water, 6, e1336,, 2019. a, b, c, d, e, f

Kolmogorov, A. N.: The local structure of turbulence in incompressible viscous fluid for very large Reynolds numbers, CR Acad. Sci. URSS, 30, 301–305, 1941. a

Koutsoyiannis, D., Dimitriadis, P., Lombardo, F., and Stevens, S.: From fractals to stochastics: Seeking theoretical consistency in analysis of geophysical data, in: Advances in nonlinear geosciences, Springer, 237–278,, 2018. a

Kramer, M. G. and Chadwick, O. A.: Climate-driven thresholds in reactive mineral retention of soil carbon at the global scale, Nat. Clim. Change, 8, 1104–1108, 2018. a

Lodge, A.: The Multiplication And Division Of Concrete Quantities, General Report, Association for the Improvement of Geometrical Teaching, 14, 47–70, 1888. a

Logan, J. D.: Applied mathematics, John Wiley & Sons, ISBN 13:978-1118475805, ISBN 10:1118475801, 2013. a, b, c

Mahajan, S.: Street-fighting mathematics: the art of educated guessing and opportunistic problem solving, The MIT Press, ISBN 9780262514293, 2010. a

Maher, K. and Chamberlain, C.: Hydrologic regulation of chemical weathering and the geologic carbon cycle, Science, 343, 1502–1504, 2014. a

Montgomery, D. R. and Foufoula-Georgiou, E.: Channel network source representation using digital elevation models, Water Resour. Res., 29, 3925–3934, 1993. a

Moon, P. and Spencer, D. E.: A modern approach to “dimensions”, J. Franklin Inst., 248, 495–521, 1949. a

Moran, M. and Marshek, K.: Some matrix aspects of generalized dimensional analysis, J. Eng. Math., 6, 291–303, 1972. a, b

Munson, B. R., Okiishi, T. H., Huebsch, W. W., and Rothmayer, A. P.: Fluid mechanics, Wiley, Singapore, ISBN 13:978-0470262849, ISBN 10:0470262842, 2013. a, b, c

Newman, W. I. and Turcotte, D. L.: Cascade model for fluvial geomorphology, Geophys. J. Int., 100, 433–439, 1990. a

Or, D. and Lehmann, P.: Surface evaporative capacitance: How soil type and rainfall characteristics affect global-scale surface evaporation, Water Resour. Res., 55, 519–539, 2019. a

Panton, R. L.: Incompressible flow, John Wiley & Sons, ISBN 13:978-111801343, 2006. a, b, c, d

Passalacqua, P., Porté-Agel, F., Foufoula-Georgiou, E., and Paola, C.: Application of dynamic subgrid-scale concepts from large-eddy simulation to modeling landscape evolution, Water Resour. Res., 42, W06D11,, 2006. a, b

Pelletier, J. D.: Fluvial and slope-wash erosion of soil-mantled landscapes: detachment-or transport-limited?, Earth Surf. Proc. Land., 37, 37–51, 2012. a

Persico, E.: Commemoration of Enrico Fermi, Enrico Fermi: His Work and Legacy, in: Enrico Fermi: his work and legacy, edited by: Bernardini, C. and Bonolis, L., Springer, p. 36, ISBN 978-3-642-06053-3, 2004. a

Pexton, M.: How dimensional analysis can explain, Synthese, 191, 2333–2351, 2014. a

Pope, S. B.: Turbulent flows, Cambridge University Press, ISBN 0521598869, 9780521598866, 2000. a, b

Porporato, A. and Yin, J.: Ecohydrology, Cambridge University Press, Cambridge, ISBN 978-1-108-84054-5, 2022. a

Porporato, A., Daly, E., and Rodriguez-Iturbe, I.: Soil water balance and ecosystem response to climate change, Am. Nat., 164, 625–632, 2004. a, b, c, d

Rayleigh, L.: The principle of similitude, Nature, 95, 66–68, 1915a. a, b

Rayleigh, L.: Reply, Nature, 95, 644, 1915b. a

Riabouchinsky, D.: The principle of similitude, Nature, 95, 591–591, 1915. a, b, c

Richter, D. D. and Billings, S. A.: `One physical system': Tansley's ecosystem as Earth's critical zone, New Phytol., 206, 900–912, 2015. a

Riebe, C. S., Kirchner, J. W., and Finkel, R. C.: Erosional and climatic effects on long-term chemical weathering rates in granitic landscapes spanning diverse climate regimes, Earth Planet. Sc. Lett., 224, 547–562, 2004. a, b

Rodríguez-Iturbe, I. and Rinaldo, A.: Fractal river basins: chance and self-organization, Cambridge University Press, ISBN 13:978-0521004053, ISBN 10:0521004055, 2001. a, b

Rodríguez-Iturbe, I. and Porporato, A.: Ecohydrology of water-controlled ecosystems: soil moisture and plant dynamics, Cambridge University Press, ISBN 9780511535727,, 2005. a

Shen, W., Davis, T., Lin, D. K., and Nachtsheim, C. J.: Dimensional analysis and its applications in statistics, J. Qual. Technol., 46, 185–198, 2014. a

Siano, D. B.: Orientational analysis – a supplement to dimensional analysis – I, J. Franklin Inst.e, 320, 267–283, 1985. a

Smits, A. J., McKeon, B. J., and Marusic, I.: High–Reynolds number wall turbulence, Annu. Rev. Fluid Mech., 43, 353–375, 2011. a

Sornette, D.: Critical phenomena in natural sciences: chaos, fractals, selforganization and disorder: concepts and tools, Springer Science & Business Media, ISBN 13:978-3540308829, ISBN 10:3540308822, 2006.  a, b

Spagnoli, A.: Self-similarity and fractals in the Paris range of fatigue crack growth, Mech. Mater., 37, 519–529, 2005. a

Sposito, G.: Scale dependence and scale invariance in hydrology, Cambridge University Press, ISBN 13:978-0521571258, ISBN 10:0521571251, 2008. a

Strahler, A. N.: Dimensional analysis applied to fluvially eroded landforms, Geol. Soc. Am. Bull., 69, 279–300, 1958. a, b

Sun, B.-H.: Scaling law for the propagation speed of domino toppling, AIP Adv., 10, 095124,, 2020. a

Szirtes, T.: Applied dimensional analysis and modeling, Butterworth-Heinemann, ISBN 13:978-0123706201, ISBN 10:0123706203, 2007. a, b

Taylor, G. I.: The formation of a blast wave by a very intense explosion I. Theoretical discussion, P. Roy. Soc. Lond. A, 201, 159–174, 1950a. a

Taylor, G. I.: The formation of a blast wave by a very intense explosion. – II. The atomic explosion of 1945, P. Roy. Soc. Lond A, 201, 175–186, 1950b. a

Theodoratos, N., Seybold, H., and Kirchner, J. W.: Scaling and similarity of a stream-power incision and linear diffusion landscape evolution model, Earth Surf. Dynam., 6, 779–808,, 2018. a

Turcotte, D. L.: A fractal interpretation of topography and geoid spectra on the Earth, Moon, Venus, and Mars, J. Geophys. Res.-Solid, 92, E597–E601, 1987. a

Voss, R. F.: Random fractal forgeries, in: Fundamental algorithms for computer graphics, Springer, 805–835,, 1985. a

Widom, B.: Scaling laws, revision #91749, Scholarpedia, 4, 9054,, 2009. a

Williams, W.: On the relation of the dimensions of physical quantities to directions in space, P. Phys. Soc. Lond., 11, 357–398, 1890. a

Xue, N. and Stone, H. A.: Self-Similar Draining near a Vertical Edge, Phys. Rev. Lett., 125, 064502,, 2020. a

Yin, J., Calabrese, S., Daly, E., and Porporato, A.: The energy side of Budyko: Surface-energy partitioning from hydrological observations, Geophys. Res. Lett., 46, 7456–7463, 2019. a, b

Zheng, Z., Christov, I. C., and Stone, H. A.: Influence of heterogeneity on second-kind self-similar solutions for viscous gravity currents, J. Fluid Mech., 747, 218–246,, 2014. a


For the purposes of this article, it may be useful to note how in mathematics one almost always assumes to deal with dimensionless quantities. For example, in Eqs. (1) and (2), the argument of the function f is assumed dimensionless (we know that if the function is transcendental, its argument should be dimensionless; see, e.g., Barenblatt (1996)). When thinking dimensionally, λ as a scale factor should have the dimension of the inverse of x. Thus, either one considers the scaling in Eq. (1) with regard to dimensionless quantities (then λ, x, and f are all dimensionless) or assumes the presence of unit factors to take care of the units: f(λ×x)=(1×a)nf(1×x), where the different “1” have different dimensions. Similarly, with x dimensional, log (x) only makes sense if one implies log(x/1), and the same is true when we develop a function in power series.


This Nature paper contains the famous quote of Lord Rayleigh: “I have often been impressed by the scanty attention paid even by original workers in physics to the great principle of similitude. It happens not infrequently that results in the form of `laws' are put forward as novelties on the basis of elaborate experiments, which might have been predicted a priori after a few minutes' consideration.”

Short summary
Applying dimensional analysis to the partitioning of water and soil on terrestrial landscapes reveals their dominant environmental controls. We discuss how the dryness index and the storage index affect the long-term rainfall partitioning, the key nonlinear control of the dryness index in global datasets of weathering rates, and the existence of new macroscopic relations among average variables in landscape evolution statistics with tantalizing analogies with turbulent fluctuations.