Modeling the Shape and Evolution of Normal-Fault Facets

Facets formed along the footwalls of active normal-fault blocks display a variety of longitudinal profile forms, with variations in gradient, shape, degree of soil cover, and presence or absence of a slope break at the fault trace. We show that a two-dimensional, process-oriented cellular automaton model of facet profile evolution can account for the observed morphologic diversity. The model uses two dimensionless parameters to represent fault slip, progressive rock weathering, and downslope colluvial-soil transport driven by gravity and stochastic disturbance events. The parameters represent rock weathering and soil disturbance rates, respectively, scaled by fault slip rate; both can be derived from field-estimated rate coefficients. In the model's transport-limited regime, slope gradient depends on the ratio of disturbance to slip rate, with a maximum that represents the angle of repose for colluvium. In this regime, facet evolution is consistent with nonlinear diffusion models of soil-mantled hillslope evolution. Under the weathering-limited regime, bedrock becomes partly exposed but microtopography helps trap some colluvium even when facet gradient exceeds the threshold angle. Whereas the model predicts a continuous gradient from footwall to colluvial wedge under transport-limited behavior, fully weathering-limited facets tend to develop a slope break between footwall and basal colluvium as a result of reduced transport efficiency on the rocky footwall slope. To the extent that the model provides a reasonable analogy for natural facets, its behavior suggests that facet profile morphology can provide useful constraints on relative potential rates of rock weathering, soil disturbance, and fault slip.


Introduction
Mountain fronts in extensional tectonic settings often display facets: steep, basin-facing hillslopes that follow the surface trace of the bounding fault and mark the transition from footwall to hanging wall ( Figure 1). In many cases, dissection of a footwall range by transverse streams creates triangular facets: facet slopes that are flanked by V-shaped transverse valleys, which present a triangular shape when viewed from the adjacent basin (see especially Figures 1a,1e,and 1i). In other settings, facets may be trapezoidal in profile or even compose a more or less continuous surface along a weakly dissected footwall range (e.g., Wallace, 1978) (Figures 1g and 1h). This paper explores the sculpting of facet cross-sectional profiles, using a process-based numerical model as an interpretive tool. We examine the extent to which the model can account for the diversity in facet morphology and in particular diversity in slope angle, regolith cover, profile shape, and presence or absence of a slope break across the range-bounding fault. We also use the model to frame testable predictions for the relationship between facet morphology, erosion rate, and fault slip rate.

Background
Fault scarps and facets have intrigued geologists since at least the late nineteenth century, when mapping expeditions in western North American brought surveyors and geologists to the spectacular terrain of the Basin and Range physiographic providence. Among the first published remarks and illustrations on facet geomorphology in the western United States were those of Gilbert (1875Gilbert ( , 1928 and Davis (1903Davis ( , 1909. Both viewed facets as exhumed fault planes, with only minor modification by erosion. Later workers, however, noted a discrepancy between the dip angle of facets and of the fault planes beneath. Where facets often dip between 20 • and 35 • (Anderson, 1977;Blackwelder, 1928;Davis, 1903Davis, , 1909Fuller, 1931;Gilluly, 1928;  Menges, 1990;Pack, 1926;Petit et al., 2009;Wallace, 1978;Wilkinson et al., 2015), the bedrock fault planes below commonly form angles of 50 • to 70 • with respect to the horizontal (Blackwelder, 1928;Gilluly, 1928;Fuller, 1931;Pack, 1926;Schneider, 1925;Wallace, 1978;Wilkinson et al., 2015). The difference in dip between a normal fault plane and the facets above it implies that facets, despite their often strikingly planar form, are erosionally modified features (Gilluly, 1928;Pack, 1926). Gilluly (1928), in his work on the Oquirrh Range (Utah, USA), pointed out an interesting implication of this erosional modification: "as the dip of the fault averages more than 60 degrees along this part of the range front, a wedge having an apical angle of 30 to 40 degrees has evidently been removed from each facet." Because the tip of a facet gets exposed to erosion earlier than the base, it undergoes more cumulative erosion. In this sense, the tip of a facet may be considered geomorphically "older" than the base (Gilbert, 1928;Menges, 1990;Wallace, 1978).
An unusually well-documented example of the contrast between fault plane and facet morphology comes from a study of the Campo Felice fault in the Italian central Apennines by Wilkinson et al. (2015). Detailed maps of the Campo Felice fault plane obtained from terrestrial laser scans of the bedrock fault scarp, together with ground-penetrating radar images of the fault plane in the subsurface, revealed a fault dipping at 57 ± 4 • , whereas the facet surface above the exposed fault plane dips at 40 ± 5 • and the debris below the fault trace dips at 36 ± 3 • . Along the Campo Felice, therefore, Gilully's "removed wedge" would have an apical angle of about 17 • .
Facets show a wide diversity in morphology. Although the dip angles of many facets range between 20 • and 35 • , facets have been reported to have dips as low as several degrees (e.g., Menges, 1990) or as high as ≥40 • (e.g., Wilkinson et al., 2015) (Figure 2). Their apices may vary in height from tens to hundreds of meters above the fault trace. Some facets are more or less continuously mantled in soil, as, for example, those along a portion of the Sangre de Cristo Range in New Mexico, USA, studied by Menges (1990), and on some facets along the eastern margin of the American Basin and Range (Figures 1d, 1f, and 1i). Others, including facets developed on carbonate rocks in the Italian central Apennines, are rocky, with a shallow, discontinuous colluvium (Tucker et al., 2011) (Figures 1g and 1h). The longitudinal profiles of facets may be planar, slightly convex upward, or slightly concave upward (Figure 2). Some faceted mountain fronts display  Table S2). a clear break in slope between the base of the facets and an adjacent colluvial apron (Figures 1g, 3b, and 3c). For example, facets and adjacent colluvial wedges surveyed by Bubeck et al. (2015) in the Apennines showed a distinct slope break. Along other mountain fronts, the basal colluvium dips at a similar angle to the facet above ( Figure 3a), with the contact between the two sometimes marked by a fault scarp (Figure 1h), and sometimes obscured (Figure 1i).
The diversity in facet morphology raises the question of whether facets may encode useful information about tectonic processes, as several studies have suggested. Hamblin (1976) and Anderson (1977) identified flights of facet-like surfaces along the Wasatch Front (Utah, USA) separated by bench-like spurs and interpreted these as reflecting alternating episodes of rapid slip and tectonic quiescence. Menges (1990) noted that facets along the southern Sangre de Cristo Range tend to be steeper and taller toward the middle of fault segments, as opposed to zones of overlap between adjacent segments. DePolo and Anderson (2000) compiled morphologic data on 45 normal faults with independent slip-rate estimates in the arid to semiarid environment of the Great Basin, USA. Faults with a slip rate in excess of a few tens of microns per year were associated with facets, and among these, DePolo and Anderson (2000) demonstrated a correlation between slip rate and facet height. In a study of facets along four segments of the Wasatch fault system, Zuchiewicz and McCalpin (2000) noted multiple potential controls on facet geometry, including lithology, but they considered slip rate to be the primary control. In laboratory experiments by Strak et al. (2011), facet angle increased with slip rate up to a limiting threshold angle. By contrast, Densmore et al. (1998) and Ellis et al. (1999) suggested on the basis of numerical model experiments that facet erosion might be controlled chiefly by bedrock landsliding, such that facet angle represents a threshold angle for stability that does not correlate with slip rate.
If facets take shape through a collaboration between tectonics and erosion, then it stands to reason that their morphology might also encode useful information about rates of geomorphic processes. Menges (1990) noted, for example, that the degree of soil development generally increases upslope on facets along the southern Sangre de Cristo Range. Using a simple geometric model, Tucker et al. (2011) noted that the difference in dip angle between a facet and its basal fault-in other words, the apical angle of Gilluly's wedge of missing rock ( in Figure 4)-could be related quantitatively to the ratio of the rates of fault slip and facet erosion. Expressed in terms of surface-normal erosion rate E n (as opposed to the arc-wise vector used in the original paper), the relation is as follows: (b) Facet and colluvial wedge profiles for which the slope of both is less than 34 • , and fractional regolith cover nears 100%, but where the facet and colluvial wedge exhibit different slopes. (c) Facet and colluvial wedge profiles for which the slope of the facet is greater than 34 • , and the surface of the facet appears rocky (fractional regolith cover <100%).
(d-f) Airborne lidar hillshade examples of profiles in relation to topography: (d) gentle and smooth slopes with uniform regolith cover; (e) intermediate slope with uniform regolith cover; and (f) steep slope with bedrock exposed. Line segments correspond to plotted profiles in (a), (b), and (c). Fractional soil cover estimated from lidar hillshade and field observations (profile coordinates listed in Table S3).
where is the fault dip angle, is the facet dip angle, and V is the fault slip rate (these and other symbols are listed in Table S1 in the supporting information). The concept is illustrated in Figure 4. One implication of the geometry illustrated in Figure 4 is that if one knew the fault slip rate and the fault dip, one could estimate the erosion rate (averaged over the age of the facet). Conversely, independent knowledge of the erosion rate would allow estimation of the slip rate.
The geometric view of facets as surfaces that erode as they emerge from below ground leads to the question of what factors determine the erosion rate on a normal-fault facet-and this in turn requires us to understand the governing geomorphic processes acting on the surface. Wallace (1978) hypothesized that fault scarps created by tectonic offset should relax relatively quickly from the initial fault dip to a more stable angle of 30 • to 37 • and thereafter lay back much more slowly. As noted earlier, Densmore et al. (1998) and Ellis et al. (1999) proposed, on the basis of numerical model experiments, that many facets form by bedrock landsliding and that the facet surfaces essentially represent failure planes. Petit et al. (2009) questioned this interpretation, in part because of the scarcity of features such as head scarps and debris lobes that one would expect to be associated with bedrock landsliding. Our own observations of facets in the Italian Apennines and the Wasatch fault system, USA, also lead us to support the view of Petit et al. (2009) that facets often undergo progressive weathering and erosion, as opposed to deep-seated landsliding (though landsliding undoubtedly does occur occasionally on some faceted mountain fronts). Menges (1990) suggested that facets might be effectively transport limited, yet the fact that some facets have extensive bedrock outcrops suggests that this is not always the case. Whereas Menges (1990) noted some evidence for slope wash, Gilbert (1928) argued that "rains accomplish little in the way of erosion … [due to] absorption and retardation by porous talus and in conditions unfavorable to concentration of flow." Finally, variation in facet morphology and scale with lithology implies that facet materials differ in their susceptibility to weathering and transport (Menges, 1990;Zuchiewicz & McCalpin, 2000).
Numerical models of extensional mountain range evolution can reproduce classic landforms such as facets, spurs, and wineglass-shaped valleys (wide in the headwaters and narrowing downstream; e.g., Leeder & Jackson, 1993), but models differ in their representation of the governing processes. Densmore et al. (1998) and Ellis et al. (1999) introduced a model that included rock weathering, regolith creep, and bedrock landsliding and explored a part of the parameter space in which facet erosion occurred primarily by landsliding. A model developed by Petit et al. (2009) represented hillslope erosion using a diffusion formulation, with a higher transport coefficient applied to slopes steeper than 40 • . In their study of facets in the Great Basin, DePolo and Anderson (2000) showed that the observed relationship between slip rate and facet angle was consistent with a nonlinear diffusion model. The planform landscape evolution model studied by Petit et al. (2009) showed little correlation between facet slope and fault dip angle, whereas the geometrical analysis of Tucker et al. (2011) implies that fault dip should be a primary control. In summary, the community has developed several published models of landscape evolution on extensional footwalls, but the implications of these models for facet evolution differ depending on the assumed process rules. To make further progress, we need models that can account for the observed diversity in facet morphology, including variations in slope angle, regolith cover, and shape, and that make field-testable predictions about the relationship between morphology, erosion rate, and slip rate.

Approach and Scope
We view normal-fault facets as unique natural experiments: slopes that are born as steep, seismotectonic fault scarps, and undergo progressive weathering and erosion as they are translated upward and away from the fault trace. We use a process-oriented cellular automaton model of facet cross-section evolution as an interpretive tool with which to address the following questions: Can a model that combines rock weathering with disturbance-driven soil creep account for the observed range in facet angle, regolith cover fraction, shape, and colluvial wedge angle? If so, what are the primary controlling factors? Does the model imply a systematic relationship between facet angle, erosion rate, and fault slip rate, and if so, what does that relationship look like? The answer to this last question is especially important, because it sets up a testable prediction: Facet angle is easy to measure, and erosion rate can in principle be obtained either through the geometric method summarized in Figure 4 (if fault slip rate is known) or by techniques such as cosmogenic nuclide analysis.
We focus on the cross-sectional geometry of facets, rather than their full three-dimensional (3-D) form; in other words, our interest is not in explaining why facets are often triangular (reflecting the geometry of transverse canyons) but rather in understanding what sets their gradient, profile shape, and regolith cover. That said, the influence of facet length-the distance between the fault trace and the facet's upper edge at a particular cross-sectional position-is considered by comparing models with varying domain length, as described below. We do not consider controls on the height of the facet edge, as this clearly depends in part on the spacing between transverse channels and on the topography of valley side slopes along them. We further restrict our consideration to unchanneled facet surfaces, where incision by concentrated flow plays little or no role. Finally, our main focus is on normal-fault facets, rather than facet-like features that form by other means, such as folding or rapid river incision (e.g., Cotton, 1950), though there are obvious similarities, and steep slopes formed by rapid vertical base-level fall have been studied theoretically using a version of the same model that we apply here (Tucker et al., 2018).
Note that we use the terms "soil," "regolith," and "colluvium" interchangeably to refer to loose granular material on a hillslope. We do not distinguish among the properties of different types of soil or other granular material but instead consider the primary material contrast to be the difference between bedrock and the disaggregated, granular material that results when that rock has been sufficiently weathered.

Cellular Hillslope Evolution Model
We model facet profile evolution using a continuous-time stochastic cellular automaton method. With this method, the spatial domain is represented using a discrete lattice of cells, while the dynamics are modeled by cell-state transition events that occur at random intervals in time, according to specified rate coefficients (Narteau et al., 2001(Narteau et al., , 2009). This approach offers several important advantages for modeling facet slope evolution. First, it accounts for a continuous range of slope forms, from smooth, soil-mantled slopes to steep, irregular, and rocky ones (Tucker et al., 2018). It provides a direct treatment of stochastic disturbances, rather than lumping their effects into a "black-box" rate coefficient, and these disturbances can trigger both short-range and long-range (nonlocal) ravel-like sediment motion, depending on topography. Nonetheless, the model's parameters can be translated to field-estimated rate coefficients for weathering and sediment transport (Tucker et al., 2018). Furthermore, because the cellular framework allows a fully two-dimensional representation (as opposed to the more common profile representation, in which surface elevation is a function of one independent spatial dimension), it provides a natural way to treat combined vertical and horizontal tectonic offset. Finally, the cellular approach described below (and in greater detail by Tucker et al., 2016Tucker et al., , 2018 honors the occurrence of a patchy or otherwise incomplete regolith cover, as is observed on some facet surfaces. The facet profile evolution model builds on the "Grain Hill" cellular automaton framework (Tucker et al., 2018), with the addition of a 60 • dipping normal fault. The model domain consists of a lattice of hexagonal cells that represents a vertical cross section through a hypothetical facet and its adjacent colluvial wedge ( Figure 5). The height of each hexagonal cell (defined as twice the apothem, which is the distance from hexagon's center to the midpoint of one of the sides) is denoted by . The lower left corner of the model represents a longitudinal stream that removes any debris delivered to it, and whose elevation can either be fixed to the hanging wall block or allowed to rise over time at a prescribed rate. In either case, the hanging wall serves as the reference frame for the model. The right side of the model represents the upper edge of the facet at that particular cross section; increasing or decreasing the width of the domain broadly equates to moving the cross section toward or away from the facet tip ( Figure 5).
Each hexagonal model cell represents one of three types of material: air, rock, or regolith ( Figure 6). Regolith cells may be stationary or may be in a state of motion in one of the six lattice directions. Collectively, these materials and motion directions are represented by assigning one of nine integer state codes to each cell in the domain (one each for air, rock, and stationary regolith, plus one for each of the six motion directions). Stochastic, pairwise transitions represent the processes of rock weathering, regolith disturbance, and ensuing regolith motion (Tucker et al., 2016(Tucker et al., , 2018. For example, the rock cell in a rock-air pair has a user-specified probability per unit time of transitioning to a regolith cell, representing weathering. Instead of clicking through a series of time steps of fixed duration, the model iterates over a sequence of these stochastic-in-time transition events, in which one or both cells in an adjacent pair change state. The algorithm works by scheduling each potential transition event at a randomly generated future time, using an exponential probability distribution function of interevent waiting times. The program then iterates in chronological order through these scheduled events. As the domain evolves, new transitions are scheduled, and some previously scheduled ones are invalidated. Tucker et al. (2016) provide a comprehensive description of the continuous-time stochastic framework and the algorithms that implement it. Tucker et al. (2018) present the rule set for the Grain Hill model. Here, we briefly summarize these rules and describe two additions to the original version: the implementation of steady slip on a 60 • dipping normal fault, and the addition of a rule that represents dissolution of bedrock. Adding the dissolution rule permits comparison of the numerical model with a simple analytical expression and also acknowledges the occurrence of facets in soluble rock, such as the platform carbonates common to the Mediterranean region.
To represent normal-fault slip, a 60 • dipping fault crosses the grid lattice at a user-specified location ( Figure 5, bottom). Fault slip is treated as quasi steady, with small displacement events occurring at regular intervals. During each slip event, cells in the footwall block are shifted up and to the right on a 60 • angle, with a displacement distance of √ 3 lattice units per time interval s ( √ 3 corresponds to shifting cells to the right by one column and upward by one and a half cell widths). We will also make use of a related quantity, = s ∕ √ 3, which is the average time interval between unit slip events. The slip rate is therefore V = ∕ [L/T], where is the width of a grid cell; slip rate is controlled in the model by setting s . Note that the square brackets here and below are used to indicate the dimensions of certain quantities, with L denoting length and T denoting time (so for example "[L/T]" should be read as "this variable has units that represent length per time," such as meters per second, fathoms per fortnight, or potrzebies per kovac; Knuth, 1955).
Production of regolith from bedrock is represented by a transition from a rock-air pair to a regolith-air pair (Figures 6a and 6b). The rate constant w [1/T] represents the average transition frequency. We denote the corresponding characteristic weathering velocity as W = w [L/T]. Given a cell width of , the expected bare-bedrock weathering rate is 2aW, where a is a surface-roughness coefficient (described below); the factor of 2 reflects the fact that for a planar surface, the hex lattice geometry exposes an average of two faces per cell. Note that the model's treatment of weathering means that it effectively occurs normal to the air-rock interface, whatever the orientation of that interface may be.
In order to explore the case of completely weathering-limited slopes, we introduce a second weathering rule to represent dissolution. When this rule is invoked, rock-air cell pairs transition to air-air pairs-representing rock dissolution-with an average rate s [1/T], and a charateristic dissolution velocity W s = s. The expected bare rock dissolution rate [L/T] is therefore 2aW s .
We assume that regolith transport can occur by two means: (1) displacement by a disturbance event and subsequent motion or (2) spontaneous gravitational failure, when the local angle of repose is exceeded. The model's disturbance rule represents the action of processes like animal burrowing, frost heave, tree throw, and other mechanisms that tend to displace regolith outward from the surface. The disturbance transition rule applies to locations where a resting regolith cell lies adjacent to an air cell. When the transition occurs, the regolith and air trade places, and the regolith state switches from resting to moving in the direction from which disturbance originated (Figures 6b and 6c). The disturbance rate parameter d [1/T] represents the average disturbance frequency and functionally equates to the disturbance frequency parameter N a in the probabilistic theory of soil creep developed by Furbish et al. (2009). The corresponding characteristic disturbance velocity is denoted by D = d. Tucker et al. (2018) show that for relatively gentle slopes, the disturbance rate relates directly to the commonly used soil transport efficiency factor ("hillslope diffusivity"), D s , according to For regolith-mantled slopes steeper than about 15 • , the effective transport efficiency increases progressively with slope angle, diverging at a 30 • effective angle of repose (Tucker et al., 2018, their Figure 10).
Moving particles follow a set of transition rules that mimic the kinematics of inelastic grain motion in a gravitational field. Although these motion and collision rules are necessarily heuristic, they effectively capture the settling motion of disturbed particles (Furbish et al., 2009). The motion, gravitational, and inelastic collision rules are illustrated in Figure 6 and described in greater detail by Tucker et al. (2016Tucker et al. ( , 2018. One rule to note in particular is that a regolith cell lying above and adjacent to an air cell can transition to a moving state, with a high transition rate parameter, representing spontaneous gravitational destabilization. Because of the lattice geometry, this rule imposes an effective 30 • angle of repose. Regolith cells on slopes steeper than this will tend to undergo spontaneous downslope motion without requiring a disturbance event to mobilize them. This treatment provides a way to represent ravel transport. Sensitivity experiments indicate that the exact nature of the motion and settling rules is not especially important (see the supporting information). What matters more is that there exists a timescale separation between disturbance (with intervals on the order of years) and settling (timescale on the order of seconds or less). The rule set described above can capture a range of slope forms, including regolith-mantled and convex-upward forms, planar angle of repose, and partially mantled "rocky" slopes (Tucker et al., 2018). This generality makes the model an appropriate one for exploring bedrock fault scarps and facets, which emerge during earthquakes as steep, rock slopes and can subsequently evolve into partially or fully regolith mantled erosional slopes. In the following section, we present a systematic parameter exploration, beginning with the simple case of a purely weathering-limited slope.

Weathering Limited Case: Facet Dissolution
We start with a simple test: If a facet erodes at a steady, uniform rate, its evolution should follow the geometry illustrated in Figure 4. The profile should be linear, and the dip angle should relate to the rates of erosion and slip according to equation (1). To perform this test, we run the model with dissolution activated and without any rock-to-regolith conversion. On a planar surface, with an average of two faces per cell, the expected average rate of erosion by dissolution would be 2 s. Surface roughness further increases the effective dissolution rate by a factor a ≈ 1.8 (supporting information section S2).
We define a dimensionless dissolution efficiency as the ratio of expected dissolution rate (2a s) to fault slip rate: Here 2a s [L/T] is the expected slope-normal erosion rate, and S ′ therefore represents the ratio of slope-normal erosion rate to fault slip rate. From equation (1), the predicted facet dip angle is where and are both in radians. Figure 7 presents simulated profiles for facets eroded by dissolution, under different values of dimensionless weathering rate S ′ . The figure also compares shorter and longer facets (top and bottom rows, respectively). The simulated facets show a linear relation between angle and dissolution rate, consistent with the analytical expectation (equation (4)) ( Figure 8). As expected, there is no apparent relation between slope angle and the width of the cross section.

Facets With Regolith
We next consider the case in which rock weathers isovolumetrically to regolith rather than dissolving, with regolith motion driven by stochastic disturbance events ( Figure 6). Two dimensionless parameters  The first of these represents the characteristic weathering velocity (W) relative to fault slip rate (V), whereas the second represents characteristic disturbance velocity (D) relative to fault slip rate.
Figures 9 and 10 illustrate the role of these two parameters for the case of a facet of fixed length. The model shows three behavior regimes. When w ′ ≫ d ′ , gradient depends on the disturbance frequency regardless of weathering rate. In this transport-limited regime, solutions with d ′ ≤ 1 produce angle-of-repose slopes (recall that the model has a 30 • effective angle of repose) (Figure 9, panels in upper left quadrant). Although these angle-of-repose solutions are not inevitable, they occupy a large part of the model's parameter space and could be thought of as an attractor state. When w ′ < 1, slope angle depends primarily on w ′ . But even in this weathering-controlled regime, disturbance rate continues to have some influence, except in the extreme dissolution-limited case when there is no regolith to move (Figure 10, solid curve and circles).
Facet shape depends on both w ′ and d ′ . With sufficiently effective weathering and transport, represented by the parameter space w ′ > 1 and d ′ > 1, the model produces convex-upward facets (Figure 9, upper right). If either of these parameters is less than unity, the model generates planar or nearly planar facets. The regolith cover that forms on these surfaces remains relatively thin (one to a few thick) because the nearly complete cover shields the bedrock from further weathering. Figure 11. Modeled regolith cover proportion for facets in quasi steady state, as a function of weathering and disturbance rate parameters. Scatter around the sigmoidal curve reflects stochastic variability (each data point represents one snapshot in time rather than a temporal average).
Facet length has little influence on slope angle (cf. Figure 10 left and right), except in the case of d ′ ≫ 1 (star and triangle symbols in Figure 10). Here, longer slopes are systematically steeper, because in this transport-limited, subthreshold mode, increased downslope regolith flux necessitates a steeper gradient.
The fractional regolith cover depends mainly on w ′ and secondarily on d ′ (Figure 11). The relation follows a sigmoid-like curve, with the steepest segment of the curve corresponding to w ′ ≈ 1 and a fractional cover of approximately 50%. The calculations show >80% cover when w ′ > 10 and less than 50% cover when w ′ < 1. For a given w ′ , a facet with a higher disturbance rate will tend have a thinner cover, all else equal. Values of d ′ ≪ 1 are associated with significantly enhanced stochastic variability in cover percentage in both space and time.

Colluvial Wedges
We explore the formation of colluvial wedges with a series of model runs similar to those in Figure 9 but with the fault location shifted to the right ( Figure 12). This geometry places the base-level outboard of the fault trace, thereby allowing the creation of a colluvial wedge between the fault and base level. In the angle-of-repose regime, with w ′ > 1 and d ′ ≤ 1, the wedge shares the same slope gradient as the facet itself, so that the fault trace appears partway up the facet slope ( Figure 12, upper left). Solutions with w ′ < 1 and d ′ < 1 generate slope breaks because the bedrock facet rises more steeply than the colluvial wedge, which can be no steeper than the angle of repose ( Figure 12, lower left). With w ′ < 1 and d ′ > 1, sediment transport becomes efficient enough that the wedge slope is less than the angle of repose ( Figure 12, lower right). Solutions with w ′ ≈ 1 and d ′ > 1 produce a relatively muted slope break in which both the facet and the wedge lie below the stability angle.

Base-Level Change
Up to now we have assumed that the base-level elevation remains fixed to the rocks in the hanging wall. In real extensional systems, however, this will rarely be the case. With the hanging wall rock as a reference frame, base level may rise as the hanging wall basin undergoes aggradation: a common situation in extensional systems. On the other hand, if the entire extensional system rises relative to an external base level such as eustatic sea level, the hanging wall blocks may themselves be incised, such that the local base level within a given hanging wall falls. Hanging wall incision has occurred, for example, in parts of the central Apennines extensional province in Italy (e.g., D' Agostino et al., 2003).
To explore how base-level rise might influence facet evolution, we performed a series of models runs identical to those in the previous section but with the addition of a numerical "sill" along the left boundary ( Figure 13). The sill's elevation rises through time at 25% of the vertical speed of footwall motion. The experiment is imperfect, as there is no easy way to keep the facet length constant. Nonetheless, the results illustrate potential impacts of hanging wall aggradation. In the reference frame of hanging wall bedrock, the fault trace rises and migrates away from the basin. Aggradation reduces the gradient of both the bedrock facet and the colluvial wedge. When w ′ and facet length are sufficiently low, the rate of accommodation space creation outpaces the sediment flux from the eroding facet, and all the incoming sediment is trapped ( Figure 13). In these cases, sediment from an external source, such as an axial river or a nearby fan fed by a transverse stream, might be expected to fill the extra accommodation space, creating a sharp contact between the facet and the (relatively flat) basin floor, with the fault trace essentially at the slope base, as is sometimes observed (e.g., Figure 1f).
When the weathering and (especially) disturbance rates are high relative to slip rate, the model realizations produce nearly flat footwall and hanging wall topography (Figure 13, top three rows of rightmost column). Effectively, the footwall is planed off as it rises, with the debris accumulating in the adjacent basin.

Relation Between Facet Angle and Erosion Rate
To examine the predicted relationship between erosion rate and facet dip angle, we define a dimensionless slope-averaged vertical erosion rate as The denominator represents the maximum rate of regolith removal by disturbance: It is the rate one would obtain if each disturbance triggers a ravel event that transports the disturbed material to the base of the slope. In other words, E ′ v = 1 can be thought of as representing "perfect" nonlocal transport (e.g., DiBiase et al., 2017;Doane et al., 2018;Tucker & Bradley, 2010). Cases with E ′ v < 1 arise either when some displaced regolith remains on the slope (indicating local, diffusive-like transport) or when the regolith cover is incomplete. Cases with E ′ v > 1 indicate some degree of direct gravitational failure of regolith material. To understand how E ′ v > 1 is possible, recall that the model allows for spontaneous mobilization of stationary regolith wherever a regolith cell lies diagonally above an adjacent air cell (Tucker et al., 2016, their Figures 14 and 15). This rule represents the gravitational failure that results when grains are perched in a configuration that exceeds their friction angle (imagine, e.g., Figure 14. Relationship between dimensionless erosion rate and steady-state facet angle, generated from a set of model runs with varying w (from 10 −5 yr −1 to 10 −2 yr −1 ), varying s (from 10 2 to 10 5 yr, corresponding to V from 0.005 to 5 mm/yr), and fixed d (10 −4 yr −1 ). Each dashed line connects a series of runs with identical w but varying s (decreasing from left to right). Color in (b) represents the ratio of actual slope-normal erosion rate, estimated from the run's last snapshot, to the theoretical maximum. For very steep cases, the estimation of erosion rate relies on just a small number of eroded cells and gives rise to variability in the erosion-rate ratio when angle is greater than about 50 • . Color scales shown by bars to right of each figure. placing a pebble on an angle-of-repose slope: the pebble rolls downhill simply because the gravitational force exceeds the frictional resistance, without the need for a disturbance event to mobilize it).
The predicted relationship between facet dip angle and E ′ v is shown in Figure 14. The figure depicts data from a set of model runs with varying weathering efficiency (w) and slip rate (V). Dotted lines connect runs with the same weathering efficiency but varying V; for each such set of points, slip rate rises progressively from left to right, accompanied by an increase in facet angle. Color shading in (a) indicates the proportional regolith cover at the end of each run, while shading in (b) shows the ratio of actual slope-normal erosion rate to the theoretical maximum that would apply if regolith were instantly removed as soon as it was created.
The results illustrate two regimes of behavior. Regolith-mantled facets show a strongly nonlinear relationship between angle and erosion rate, with a rapid acceleration in erosion rate near the model's 30 • angle of repose (yellow [light-colored] points in Figure 14a). Rocky facets, on the other hand, have a roughly log-linear relationship (blue [dark] points in Figure 14a).
The first regime represents transport-limited behavior. As slope angle increases, the model predicts a nonlinear increase in average erosion and transport rates and a transition in morphology from convex upward to planar (Tucker et al., 2018). For fully regolith-mantled slopes, the model's 30 • effective angle of repose imposes an upper limit (yellow points in Figure 14a). This behavior is consistent both with the Andrews-Bucknam transport law (Andrews & Bucknam, 1987;Roering et al., 1999) and with experimental and field data on soil-mantled hillslopes (Binnie et al., 2007;DiBiase et al., 2012;Ouimet et al., 2009;Roering et al., 2001;Roering, 2008;Stock et al., 2009).

10.1029/2019JF005305
The second regime represents weathering-limited behavior. In this regime, for a given weathering efficiency (w), the erosion rate rises approximately exponentially with slope angle, as indicated by the roughly linear trends for facets with less than ∼80% regolith cover in Figure 14a (note the logarithmic axis). The relationship reflects increased regolith mobility, and decreased cover fraction, on steeper slopes. Because of surface roughness, the regolith is not infinitely mobile but rather can become temporarily trapped and stored in microdepressions. These microdepressions appear even on very steep simulated bedrock facets, like those along the bottom row of Figure 9. One consequence is that the slope-normal erosion rate remains below the theoretical maximum of E max = 2a w, which would represent a case of pure dissolution and no regolith production (symbol colors in Figure 14b illustrate the ratio of actual erosion rate to this theoretical maximum). This weathering-limited regime illustrates an interesting feedback: Weathering produces microtopographic roughness features that trap colluvium, preventing slopes above the angle of repose from become completely bare and (according to the rules of the model) shielding some of the bedrock from weathering. As the slip rate increases, the ability of microtopography to trap material is diminished, the fractional cover declines, and the erosion rate rises.

Understanding the Variability of Facet Form
Facets and their colluvial toes exhibit a fascinating variety of forms, with notable differences in slope angle, regolith cover, profile shape, and the presence or absence of a slope break at the fault trace. Although the process model developed here has only two dimensionless parameters (d ′ and either w ′ or s ′ ), it can reproduce much of the observed variety. The model manifests two general modes of behavior. Modeled slopes for which weathering efficiency is much greater than slip rate (w ′ ≫ 1) exhibit transport-limited behavior, in which slope gradient is largely insensitive to weathering efficiency ( Figure 10, points with w ′ > 1). In this mode, the equilibrium dip angle for modeled facets with a slip rate greater than the characteristic rate of soil disturbance (d ′ ≤ 1) is the angle of repose (30 • in the model) (Figure 10, small circles and pluses at w ′ > 1). When disturbance frequency exceeds the slip rate, facet slopes act as diffusion-like hillslopes, with convex-upward longitudinal profiles, and average dip angles below the angle of repose (illustrated by the upper right panels in Figure 9 and by the stars and triangles at w ′ > 1 in Figure 10). In either case, the facet's profile shape and dip angle are largely insensitive to weathering efficiency.
By contrast, facet slopes with a slip rate faster than the potential weathering efficiency (w ′ < 1) operate in a weathering-limited mode. In this mode, dip angle depends on weathering efficiency, slip rate, and disturbance frequency (Figure 10, points for w ′ < 1). In this mode, facets are planar and may be steeper or gentler than the angle of repose. Although weathering limits the rate of erosion in this mode, the processes of soil disturbance and downslope transport still exert an important influence; all else equal, facets with lower disturbance frequency are predicted to be steeper, and vice versa. The morphologic signature that most clearly distinguishes weathering-limited and transport-limited modes is the proportional regolith cover, which depends primarily on weathering efficiency and secondarily on disturbance frequency (Figure 11).
The model also predicts that a break in slope between the facet surface and the colluvial apron below will occur either when the facet dip angle exceeds the angle of repose or when the lack of regolith limits the transport and erosion rate on the facet (Figure 12). Both cases occur in the weathering-limited regime and are marked by an incomplete regolith cover on the facet. The slope break reflects a contrast in erosion and transport efficiency between the rocky facet surface and the colluvial wedge. For a system with a relatively steady base level, presence of a slope break suggests a weathering-limited regime. The converse is not quite true; model realizations in the transitional regime of w ′ ≈ d ′ ≈ 1 show only a subtle slope break that might be difficult to detect in the field, especially under base level rise (Figures 9 and 12).
A terrestrial laser scanning study of the Campo Felice fault in the Italian central Apennines provides an interesting case study of an active fault with a slope break between facet and colluvial wedge. A set of 25 topographic profiles collected by Wilkinson et al. (2015) reveal a footwall slope angle of 40.3 • ± 4.6 • , with the colluvial wedge below dipping at 35.7 • ± 3.2 • . The facet surface on the footwall is rocky, with a colluvial cover fraction that we estimate from a reconnaissance survey to be ≈80%. Both the slope break and the incomplete colluvial cover imply a weathering-limited regime.
Similarly, field observations of facets along the Wasatch fault zone are also consistent with the model. Where slopes are below the angle of repose, bedrock breaks the surface only rarely, and no slope break is seen between facet and wedge (Figures 3a and 3d). Where facet angles reach and begin to exceed likely angles of repose, bedrock exposure becomes much more obvious and also highly spatially variable. This is accompanied by the appearance of a distinct break in slope between the facet and the wedge. These observations suggest a transition from the w ′ > 1 to w ′ < 1 domains, presumably driven by spatial variation in fault slip rate on different segments of the fault.

Tectonics From Process and Topography
The two-parameter cellular model implies a systematic relationship between facet angle, regolith cover, and erosion rate. For fully regolith-covered slopes, the model's behavior is broadly consistent with the Andrews-Bucknam transport law, which expresses a nonlinear relation between sediment flux per width, q s , and slope gradient, S: where D s [L 2 ∕T] is a rate coefficient and S c is a threshold slope. The key difference is that whereas equation (5) diverges at S = S c (and therefore only applies when S < S c ), the cellular model implies the possibility that S can exceed S c . The transition from below-threshold to above-threshold form coincides with a shift from soil-mantled to partly rocky and represents a limitation on regolith transport as the cover thins. The limitation has two elements. First, the reduction in cover limits the degree of mobilization by disturbance. Second, microtopography limits the transport length of disturbed particles and effectively traps regolith even when the average slope angle exceeds the angle of repose. Both effects must occur on real hillslopes. Capturing them with a geomorphic transport law may require a nonlocal formulation that replaces the divergence in the Andrews-Bucknam law with a reduction in transport rate as regolith cover shrinks. One potential approach would be an explicit representation of characteristic transport length, which would vary both with S∕S c and with the ratio of mean regolith thickness to bedrock roughness.
Another way to envision the slope-erosion relation is to consider a set of facets with identical lithology and climate (and hence fixed W and D) but varying slip rate ( Figure 15). Facets with a lower rate of slip would tend to be fully regolith mantled (unless W is especially low), and for this subset, increasing slip rate would be met with an increase in facet angle (Figure 15, points). Once the angle comes close to the threshold angle for colluvium, any further increases in slip rate would be accompanied by an increase in gradient above the threshold angle and a reduction in colluvial cover fraction (Figure 15, dark-colored points).
The cellular model suggests an interesting role for topographic roughness. The stochastic nature of weathering and regolith-disturbance events creates surface roughness, and especially so when a substantial fraction of rock is exposed at the surface (Figures 7 and 9). The development of roughness can accelerate rock weathering because it increases the exposed surface area. At the same time, as noted above, a rough surface can trap and hold regolith even when the average slope angle is well above the angle of repose, which allows even very steep slopes to retain some degree of regolith cover. On vegetated hillslopes, vegetation also contributes to sediment trapping (DiBiase & Lamb, 2013;Doane et al., 2018). In the model, trapping of regolith by roughness elements retards weathering by shielding the rock. In reality, it is possible that a thin regolith cover could actually accentuate rock weathering by trapping chemically reactive water and/or enhancing plant growth and vegetation-related weathering, as hypothesized by Gilbert (1877). The consequences of such an effect for facet evolution might be explored numerically with a modified version of the cellular model that allows faster weathering for rock in contact with near-surface regolith.
One indicator of relatively efficient weathering is the presence of a more or less complete soil mantle. Such soil-mantled facets clearly exist, as documented, for example, by Menges (1990). Given that slip rates on faceted mountain fronts are commonly 0.1-2 mm/yr (DePolo & Anderson, 2000), an implication is that maximum rates of rock weathering on facets can be at least this high. Rates on that order are broadly consistent with the findings of Heimsath et al. (2012), who demonstrated that regolith production tends to be faster on steeper slopes.

10.1029/2019JF005305
To test the idea that facet angle and soil cover depend on rates of weathering and fault slip, one place to look would be extensional systems in a hyperarid environment. In such an environment, faults slipping at millimeters per year would be expected to produce rocky facets steeper than the threshold angle for granular material, all else equal (and assuming the rock is strong enough to resist deep-seated landsliding). The same could be said of facets formed in highly resistant lithologies such as quartzite.
One limitation to the scenarios considered here regards the potential role of partial dissolution in soluble lithologies. For example, carbonates underlie many of the prominent facets in the Mediterranean region, and facets we have surveyed in central Italy often host a thin (order 10 cm) and discontinuous soil cover (Tucker et al., 2011). Copious clastic debris in fans and colluvial wedges implies that mechanical weathering plays a major role, but it is also possible that partial dissolution of hillslope colluvium results in a soil cover that is thinner than it otherwise would be (it is also conceivable that the soils have thinned over historic time as a result of heavy grazing).
The behavior of the cellular model is consistent with a geometric model that relates facet angle to fault dip and the ratio of slope-normal erosion rate to slip rate (equation (1)). An important caveat in applying this model is that the erosion rate can be expected to vary systematically with slope angle, even on weathering-limited slopes. The cellular model predicts that the relationship between erosion rate and slope angle depends on the degree of regolith cover, with transport-limited slopes showing a nonlinear relationship in which erosion rate rises rapidly as the threshold angle is approached, and weatheringlimited slopes exhibiting a gentler, exponential-like relationship ( Figure 14). Hillslopes can transition from transport-limited to weathering-limited under changes in slip rate, as might be seen along strike in major normal fault systems, or with variation in weathering and transport efficiency through time under a changing climate. For all but the most strongly weathering-limited systems, this transition creates a notable kink in erosion rates close to the angle of repose, but whose exact form will be challenging to constrain from field observations alone. This kink also in part explains the tendency of real facets to show dip angles clustered in the range of 32-40 • : In this range, changes of only a few degrees can accommodate large changes in erosion rate, creating a geometric attractor state in this window. The occurrence of facets considerably steeper than 40 • requires a weathering efficiency (maximum weathering rate) less than the rate of fault slip (in other words, w ′ < 1). However, where this does occur, it is likely that facets offer records of past slip rates.

Testable Predictions
The prediction that erosion rates exhibit significant and highly nonlinear variations with facet slope could be tested by obtaining erosion-rate estimates from facets above faults with an independently known slip rate, perhaps using cosmogenic radionuclide methods. The model suggests that the form of this relationship will be challenging to predict a priori for a given environment from field observations of slope and bedrock exposure alone (e.g., Figure 14). However, given independent constraints on actual erosion rates across the full range of slopes in a given facet array, the model suggests the potential to quantify fault slip rate from facet form. Methods that assume either a constant erosion rate or an Andrews-Bucknam style slope threshold will both result in poor calibration of slip rate from facet slope.
The predicted correlation between facet angle and fractional rock exposure should also be testable. Recent work using lidar topography data has identified correlations between bedrock exposure and variables such as drainage density, erosion rate, surface roughness, fracture density, and mean slope gradient (DiBiase et al., , 2018Milodowski et al., 2015;Rossi et al., 2020). One challenge is that bedrock exposure metrics tend to rely on local slope gradient as a proxy, which carries the risk of spurious correlation. Overcoming this risk would require a slope-independent bedrock exposure proxy, such as one based on photographic imagery rather than topography.

Model Applications and Limitations
The cellular model we have explored has the advantage of simplicity (just two dimensionless parameters, w ′ and d ′ ), as well as a treatment of soil transport that honors its episodic and stochastic nature. Nonetheless, there are plenty of limitations. In terms of geomorphic process, the model assumes that transport and erosion by overland flow are negligible. This assumption may be justified for the planar, unincised cases that we address here but is clearly not applicable to facets that are sliced by gully networks. We have also not considered the potential role of deep-seated landsliding in limiting the relief of facets on fast-slipping faults, though as noted earlier, this process does not seem to be especially common. In terms of tectonic processes, we have ignored both fault-plane rotation, and the possibility that the near-surface fault plane dip might not be identical to the long-term slip direction (e.g., in softer material, the near-surface fault dip may effectively be steeper than the slip orientation; e.g., McCalpin, 2009). We have also ignored the possibility of time-varying fault motion, which, as Hamblin (1976) suggested, might leave an imprint on mountain-front morphology. Similarly, we have not addressed transient facet growth, fault breccia (or other forms of lithologic heterogeneity), varying fault dip, the occurrence of multiple parallel faults, or the effect of a basinward shift in the position of the active range-bounding fault. These and other effects could be explored theoretically using the cellular framework, which could in turn point toward additional field-testable predictions.
Questions concerning the formation and evolution of bedrock fault scarps could be addressed by incorporating discrete slip events, as opposed to the quasi-continuous slip considered here. Used in that fashion, the model could provide a useful tool for analyzing and interpreting cosmogenic radionuclide samples on bedrock fault scarps. The numerical model could also provide a template for developing an improved continuum theory of soil transport on steep slopes: one that overcomes the limited range of applicability of the Andrews-Bucknam law and accounts for the role of microtopography in trapping sediment, as well as the role of bedrock exposure in limiting downslope sediment flux.

Conclusions
The longitudinal profiles of normal-fault facets show a wide diversity in form. Facets can vary in slope angle from several degrees to over 40 • . Soil entirely blankets some of them, while others reveal bare rock with a patchy colluvial cover. Some show a distinct slope break between the facet and the colluvial wedge below, while others exhibit a continuous, uniform gradient across the fault trace from footwall to basal colluvium. Many have gradients close to a threshold angle for colluvium, but many others dip more gently.
A stochastic cellular automaton model of hillslope evolution with two dimensionless parameters can account for the observed diversity in longitudinal profile form. The model addresses the morphotectonic evolution of the near-surface portion of the footwall of a 60 • dipping active normal fault. It uses transition rules to represent bedrock weathering and the intermittent disturbance and transport of regolith. As shown in previous work, the model's parameters have a physical meaning and can be related to the more common parameters of geomorphic transport laws for regolith production and downslope transport. The two dimensionless process parameters represent characteristic velocities for weathering and regolith disturbance, respectively, normalized by the fault slip rate.
The model exhibits two modes of behavior: transport limited and weathering limited. Under the former regime, which represents relatively large weathering efficiency, facet angle depends on disturbance frequency and fault slip rate. Regolith fully blankets the modeled facets, and the threshold angle for granular-material failure imposes an upper limit to gradient. This soil-mantled, threshold behavior emerges when the potential weathering rate exceeds the fault slip rate, but the slip rate outpaces the disturbance rate. In this regime, modeled slopes show a continuous gradient across the fault trace, from colluvial wedge to footwall.
Weathering-limited behavior occurs when the fault-slip rate exceeds the maximum regolith production rate. Modeled facets in this regime show intermittent regolith cover, with the cover percent depending mainly on the ratio of (potential) weathering and slip rates and secondarily on disturbance frequency. Gradient may exceed the threshold angle for colluvium; in this case, a partial colluvial cover remains held in place by microtopography. If weathered rock dissolves instead of transitioning to mobile regolith, the model recaptures a geometric analytical solution for facet angle as a function of the ratio of erosion to slip rate. Under weathering-limited behavior, modeled facets can exhibit a slope break between footwall and colluvial wedge when either the footwall slope exceeds the threshold angle, or disturbance frequency is sufficiently high (or slip rate sufficiently low) to transport eroded material at a lower-than-threshold gradient. Despite the term, the model's weathering-limited mode involves important interactions between weathering and regolith transport processes. For example, the rate of erosion in this mode depends on regolith disturbance frequency, which helps set the degree of rock exposure.
The model predicts a nonlinear relation between erosion rate and gradient. The nature of the relation depends on the mode. A soil-mantled, transport-limited facet slope shows a nonlinear curve that resembles the Andrews-Bucknam transport law, with a nonlinear approach to an asymptote at the threshold angle for stability. But where the Andrews-Bucknam law diverges at the threshold angle and says nothing about 10.1029/2019JF005305 what happens above it, the cellular model implies a transition to weathering-limited behavior when gradient exceeds the threshold. For rocky, weathering-limited slopes, the model predicts an exponential-like relation between erosion rate and gradient.
The degree to which the model provides an accurate representation of facet profile evolution remains to be tested. Estimating the erosion rates on facets of varying angle, lithology, and soil cover using cosmogenic nuclide analysis could provide a test of the predicted controls on erosion rate. Direct collection and measurement of downslope sediment flux could provide a shorter-term test that might also illuminate the governing regolith production and transport processes. Ultimately, such data, in combination with the process-based model presented here (or a suitably modified version thereof), could provide the basis for a geomorphic transport law that captures the transition in behavior across the threshold angle for soil stability and the related transition from soil-mantled to rocky slopes.