Modeling the formation of toma hills based on fluid dynamics with a modified Voellmy rheology

Toma hills are perhaps the most enigmatic morphological feature found in rock avalanche deposits. While it has been proposed that toma hills might emerge from the fluid-like behavior of rock avalanches, there still seems to be no consistent explanation for their occurrence. This paper presents numerical results based on a modified version of the Voellmy rheology, which was recently developed to explain the long runout of rock avalanches. In contrast to the widely used original version, the modified Voellmy rheology defines distinct regimes of Coulomb friction at low velocities and velocity-dependent friction at high velocities. When movement slows down, returning to Coulomb friction may cause a sudden increase in friction. Material accumulates in the region upstream of the point where this happens. In turn, high velocities may persist for some time in the downstream and lateral ranges, ultimately resulting in a thin deposit layer. In combination, both processes generate more or less isolated hills with shapes and sizes similar to those of toma hills found in real rock avalanche deposits. Thus, the modified Voellmy rheology suggests a simple mechanism for the formation of toma hills.

1 Introduction

Toma hills are cone-, pyramid-, or roof-shaped hills consisting mainly of landslide material (Abele, 1974). While toma hills seem to be more or less isolated, they are typically clustered. In the European Alps, such structures have been documented in the Rhine Valley close to Domat/Ems (Abele, 1974), at Fern Pass (Prager et al., 2006), in Almtal (van Husen et al., 2007), in Obernberg (Ostermann et al., 2012), at Eibsee (Ostermann and Prager, 2016), and at Pragser Wildsee (Ostermann et al., 2020).

The formation of toma hills is still enigmatic. Results from dating (e.g., Ivy-Ochs et al., 2009; Ostermann et al., 2012) refute any effects of glaciation. Since the deposits of rock avalanches are typically hummocky, it makes sense to assume that toma hills are already formed during emplacement and are not the result of later modification. Ostermann et al. (2012) have discussed the occurrence of transverse ridges and toma hills in the valley of Obernberg in the context of the fluid-like behavior of rock avalanches. Besides the potential occurrence of waves and surges, the collapse of this fluid-like behavior below a minimum required kinetic energy level also played a part in this discussion.

However, there is still no direct evidence that toma hills emerge from the fluid-like behavior of rock avalanches. Several laboratory experiments on granular flow have been conducted, with a focus on rock avalanches (e.g., Pouliquen et al., 1997; Shea and van Wyk de Vries, 2008; Paguican et al., 2014; Valderrama et al., 2018), and have been partly able to predict the occurrence of hummocky topographies and ridges. Cone- or pyramid-shaped hills, however, have not been found. The same holds for particle-based numerical simulations (e.g., Campbell et al., 1995; Thompson et al., 2010; Mead and Cleary, 2015; Johnson et al., 2016) and for continuum simulations based on the fundamental theory proposed by Savage and Hutter (1989).

The lack of direct evidence from granular dynamics has motivated alternative approaches. More and Wolkersdorfer (2019) revisited the idea of a long-term modification of the topography via internal erosion. Knapp et al. (2022) proposed a concept involving an already-existing lake, which may have been formed by a previous rock avalanche. This idea is, however, still on a qualitative level and appears to be site-specific.

This study takes up the idea that toma hills emerge from the complex, fluid-like behavior of rock avalanches. It starts with the modification of the Voellmy rheology (Voellmy, 1955; Salm, 1993) proposed by Hergarten (2024c). The Voellmy rheology was originally developed for snow avalanches, but nowadays, it is also applied widely to rock avalanches. In its original form, however, it is not able to reproduce the long runout of large rock avalanches without assuming an artificially low coefficient of friction. The modified rheology reproduces the long runout reasonably well without assuming specific processes that reduce friction, such as frictional heating (Erismann, 1979; de Blasio and Elverhøi, 2008; Lucas et al., 2014) or acoustic fluidization (Johnson et al., 2016).

The existence of two distinct regimes of granular flow is the main idea behind the modified Voellmy rheology. The ideas of Voellmy (1955) are adopted for high velocities, resulting in an effective friction that is proportional to the square of the velocity. In turn, Coulomb friction is assumed at low velocities. The velocity at which the transition occurs was derived from reinterpreting the concept of random kinetic energy (Buser and Bartelt, 2009; Bartelt and Buser, 2010). As a key point, friction may be much lower in the fast regime than in the slow regime. In this case, friction may increase suddenly when movement slows down and returns to Coulomb friction. In this sense, the modified Voellmy rheology is a perfect starting point for further elaborating on ideas about a fluid dynamical origin of toma hills, as discussed by Ostermann et al. (2012).

2 Approach

There was no specific strategy for selecting a suitable study site except for keeping it simple. As a first step, the model for rockslide disposition proposed by Hergarten (2012) was applied to a given topography – specifically, a digital terrain model of Tyrol with a mesh width of 5 m (Land Tirol, 2023). The model was only used for detecting locations where a large volume could be detached along a sufficiently steep failure surface. This model is a statistical model of progressive failure and involves only two parameters, s_min and s_max. These two values define the range of slopes over which the probability of failure increases from 0 to 1, whereby the values used by Argentin et al. (2021) (s_min=1; s_max=5) were implemented in the present study. The most recent version of the model was used, which also takes into account the local orientation of the failure surface during progressive failure. This modification prevents the failure of entire mountains in multiple directions and thus reduces the number of obviously unrealistic detached volumes.

The valley of Obernberg was selected as a study site by screening the spatial distribution of 10⁶ potential events visually. Figure 1 shows the topography, including the predicted area of detachment (red), where the largest predicted volume (about 8×10⁷ m³) was chosen. It should be kept in mind that this is only the initial state of the simulation of the runout performed in this study and not a result. Since the model does not take into account the existence and orientation of faults or any other properties of the rocks, we cannot tell whether such a rockslide is geologically possible at this location or even assign a probability to it.

Figure 1 The study site, located in the valley of Obernberg. The red region depicts the assumed detachment area. The black line illustrates the scarp of the real rock avalanche (redrawn from Ostermann et al., 2012). Yellow patches indicate removed local maxima. Coordinates refer to EPSG:31254 (MGI/Austria GK West).

While the chosen initial state is somewhat arbitrary, the valley of Obernberg is particularly interesting because it was subject to a rock avalanche in the past, which was investigated in detail by Ostermann et al. (2012). Since the valley already contains toma hills, its morphology should, in principle, be suitable for generating toma hills. The assumed detachment area is in a different branch of the valley than the original event. An advantage of considering a rockslide in the different branch is that the existing deposits can interfere with the simulated rock avalanche only in the lower part of the valley, whereas the original branch was reshaped considerably by the deposits of the real rockslide. While the lower part of the valley contains only a small fraction of the deposited volume, it hosts the toma hills. Therefore, the existing hills were removed before simulating the rock avalanche. For this purpose, closed contour lines in the valley were computed at 0.25 m intervals. All points within these contour lines with an elevation higher than the contour line were replaced by minimal surfaces (yellow areas in Fig. 1). In combination with the closed contour lines, this procedure is almost the same as clipping the areas inside the contour lines to the elevations of the respective contour lines, but it yields a slightly smoother transition at the edges.

Version 2 of the model MinVoellmy (Hergarten, 2024e) was used for simulating the rock avalanche. As a modification of version 1, described by Hergarten (2024d), it uses thickness-weighted central difference quotients for the gradient of the surface. One-sided difference quotients toward the steepest descent are only used at local maxima. This version was also included in the model intercomparison by Wirbel et al. (2024).

The model MinVoellmy implements the modification of the Voellmy rheology proposed by Hergarten (2024c). Friction is represented by shear stress, expressed as

$$\begin{array}{}\text{(1)}& \mathit{\tau }=\left\{\begin{array}{lcl}\mathit{\mu }\mathit{\sigma }& & v<v_{\mathrm{c}}\\ & \mathrm{for}& \\ \frac{\mathit{\rho }g}{\mathit{\xi }}{v}^{\mathrm{2}}& & v\ge v_{\mathrm{c}},\end{array}\right.\end{array}$$

at the bed. The upper row describes Coulomb friction with the coefficient μ, where σ represents the normal stress. The lower row describes velocity-dependent friction, as originally proposed by Voellmy (1955) for snow, where ρ, g, and v represent density, gravity, and vertically averaged velocity, respectively. The parameter ξ refers to the roughness of the bed.

In its original form, the Voellmy rheology assumes that the shear stress is the sum of the two terms in Eq. (1). Defining two distinct regimes of friction separated by a crossover velocity (v_c), which depends on the thickness (h) of the mobile layer, was the central idea behind the modified Voellmy rheology. Hergarten (2024c) derived a particular relation,

$$\begin{array}{}\text{(2)}& v_{\mathrm{c}}\propto \sqrt[\mathrm{3}]{\mathit{\xi }h},\end{array}$$

by reinterpreting the random-kinetic-energy model (Buser and Bartelt, 2009; Bartelt and Buser, 2010), which describes the supply of kinetic energy from random particle motion and its consumption. Equation (2) requires a factor of proportionality, which is defined by v_c at h=1 m (for a given value of ξ) in MinVoellmy. All occurrences of ρ finally cancel out since driving gravitational forces and normal stress (σ) are proportional to ρ. Thus, MinVoellmy requires three physical parameters (μ, ξ, and v_c at h=1 m).

Several simulations were performed with different values of ξ and v_c, while μ=0.75 was kept constant as a typical value for the Coulomb friction of rocks. Constant time increments of δt=0.01 s were used, ensuring that the transport distance at each step was much smaller than the grid spacing and thus maintaining numerical stability according to the Courant–Friedrichs–Lewy criterion. As a further technical setting, it was assumed that layers thinner than 1 mm cannot move.

3 Results and discussion

3.1 Morphology of the deposits

The parameter combination of ξ=500 m s⁻² and v_c=4 m s⁻¹ was selected for an in-depth analysis. The influence of these two parameters will be investigated in Sect. 3.4. Figure 2 shows the final deposits, which extend about 1 km further downstream than the deposits of the real rock avalanche described by Ostermann et al. (2012). Figure 3 shows a close-up of the deposits in the valley, with contour lines of peaks with more than 3 m of prominence. The respective outermost closed contour line defines the base for the height, and the inner contour lines are plotted at intervals of 2 m.

Figure 2 Final deposit thickness. The black outline encloses the region affected by the rock avalanche.

Figure 3 Morphology of the deposits in the valley. The colored lines represent closed elevation contour lines drawn at 2 m intervals around hills with at least 3 m of prominence, with the outermost contour line defining the base level. The red regions indicate the thickness of the deposits, as shown in Fig. 2. The black and gray lines, along with the labels, correspond to Figs. 5–7.

Figure 4 shows the results of the same analysis for the real topography. The biggest hills are obviously smaller than those obtained in the simulation. While the maximum prominence is about 19 m in the simulated deposits, it is only about 12 m in the real topography. However, a two-sample Kolmogorov–Smirnov test does not reveal any differences between the statistical distributions of the prominence values, with a p value of 56 %. Thus, it cannot be concluded whether the simulation simply yields more hills compared to the real topography (42 vs. 17) or whether the hills are also systematically larger.

Figure 4 Hill structures in the real topography. The colored lines represent closed contour lines drawn at 2 m intervals, with the outermost contour line defining the base level. The black lines correspond to Fig. 5.

Comparing the shapes in plan view is, however, difficult beyond the visual impression of some similarity. To investigate the cross-sectional geometry, the prominence of each hill is considered in relation to two different radii. One of them is the minimum radius, defined by the shortest distance from the peak to the outermost contour line, illustrated by black lines in Figs. 3 and 4. The second radius is the mean radius, which is obtained by converting the outermost contour line into a circle of equal area. Two-sample Kolmogorov–Smirnov tests do not reveal any differences between the statistical distributions of the simulated and real radii, with p values of 33 % (minimum radius) and 75 % (mean radius).

Figure 5 shows the prominence of all hills vs. the two radii. The straight lines were obtained by fitting linear functions through the origin. While individual hills scatter strongly around these lines, the lines themselves agree strikingly well between simulated and real hills. This means that the simulated hills have practically the same steepness as the real hills on average. The respective mean slopes are 0.46 in the steepest direction (minimum radius) and 0.23 on average (mean radius).

Figure 5 Prominence of hills in relation to the minimum radius (black lines in Figs. 3 and 4) and the mean radius. The straight lines represent fitted linear functions through the origin.

Figure 6 shows profiles of the numbered hills from Fig. 3 in the steepest direction. As a first observation, some profiles (hills 3, 4, and 6) contain a more or less straight segment that is considerably longer than the minimum radius, illustrated by the black lines. This property is due to the irregular shape of the hills and emphasizes the role of prominence as a lower bound for the height of the respective hill. For the two steepest hills (hills 3 and 5), the slope of the straight segment is close to the coefficient of Coulomb friction, μ=0.75. This value defines the maximum stable slope. Thus, some flanks of the hills are at the limit of stability. The majority, however, are less steep, although some are also partly straight.

Figure 6 Profiles across the 12 numbered hills from Fig. 3 in the direction of the shortest radius. The black lines refer to the black lines shown in Fig. 3. Negative values of the radius point in the opposite direction. Filled areas illustrate the thickness of the deposits.

3.2 Are the hills isolated?

As a second characteristic property beyond surface morphology, Abele (1974) described toma hills as isolated objects. However, the question of what this means exactly is not trivial. While toma hills typically look like islands standing out from an alluvial plain, some clustering is often found. Ostermann et al. (2012) even found connected hills forming ridges in the valley of Obernberg, which can partly be recognized in Fig. 4. These observations raise the question of whether toma hills are typically isolated patches of rockslide material or outcrops of a continuous layer in the subsurface. Resistivity measurements obtained in the Rhine Valley by Knapp et al. (2022) could not answer this question uniquely but suggest that there may be some rockslide material beneath the alluvial deposits.

As already visible in Fig. 3, the hills generated by the model are not just hummocks at the surface of a thick layer of rockslide material. There are even areas between the hills with very little thickness. However, several hills are not clearly separated. Some of the connections are like narrow bridges. Revisiting Fig. 6 confirms that the thickness between the hills is irregular. Hills 4 and 8 approach a very small deposit thickness in their steepest direction. For the majority of the hills, however, a thickness of several meters persists at the foot.

Figure 7 shows a longitudinal valley profile drawn through hills 1 and 10. Hill 4, which is crossed at a saddle between two local maxima, looks quite isolated in this direction. Hills 5 and 10 are not isolated in the profile direction but are much higher than the minimum deposit thicknesses between the respective hills and their neighbors. In turn, the lower part of the profile is characterized by a more continuous deposit layer, which is partly related to the road.

Figure 7 A longitudinal valley profile drawn through hills 1 and 10. For better identification with Fig. 3, easting is plotted instead of length along the profile. Elevation is exaggerated by a factor of 4.

Formally, however, the hills are not separated clearly in terms of deposit thickness. Figure 8 shows the minimum thickness around each hill, defined by the lowest thickness contour line around each hill that does not include another hill. This means that the respective hill would be isolated completely if we removed this thickness. The lowest thickness is 2.9 m, and only 3 out of the 42 hills are separated by a thickness less than 5 m. This separation is weaker than what can be expected from the previous considerations. It is, however, related to the complex pattern of the deposits already recognized in Fig. 3. Some hills are quite close to each other and are thus not separated clearly, but narrow, ridge-like connections also contribute to the weak separation.

Figure 8 Minimum deposit thickness around each hill. Colors indicate prominence, and marker sizes represent the area of prominence, defined by the area within the respective outermost contour line shown in Fig. 3.

Overall, the model predicts more or less isolated hills. This finding is not inconsistent with previous knowledge about toma hills. Taking into account the fact that the uppermost rockslide deposits may be reworked by fluvial processes or even anthropogenically, validating the results of the model quantitatively may be challenging.

3.3 The mechanism of formation

The formation of hills is closely related to the discontinuity in friction at the velocity (v_c) assumed in the model. If the fluid-like behavior collapses at any point, friction increases suddenly, which causes a rapid decrease in velocity and, ultimately, some kind of traffic jam upstream of the respective point.

Figure 9 illustrates at which time the fluid-like behavior collapses in a section of the valley. Overall, this collapse proceeds downstream due to the low friction when v≥v_c. If friction were high, particles would decelerate strongly when moving downstream, which could cause the collapse to start at the front and propagate upstream. Then, the entire mass could accumulate in a small area.

Figure 9 Time of collapse of the fluid-like behavior. The black lines represent the outlines of the hills shown in Fig. 3.

Since this is not the case here, the formation of hills hinges on local variations in the time of collapse. As expected, the hills in Fig. 9 are primarily formed in zones of early collapse or slightly upstream. Local minima of the time of collapse are particularly relevant here. The yellow points in Fig. 10 represent such local minima, with a (negative) prominence of at least 1 s. Thus, the collapse spreads from these points in all directions for at least 1 s.

Figure 10 Local minima of the time of collapse of the fluid-like behavior, with a (negative) prominence of at least 1 s. Contour lines are drawn around the minima at intervals of 1 s. The black lines represent the outlines of the hills shown in Fig. 3, and the black dots indicate the respective nearest peaks in the topography. The labels correspond to those in Fig. 11.

The solid lines in Fig. 11 show the velocities and thicknesses of the five seed points labeled in Fig. 10. The collapse of the fluid-like behavior takes place during a phase of decreasing velocity. However, this does not imply that the particles are already significantly decelerated since the plot does not follow moving particles but instead refers to the particles at a given position (i.e., a Eulerian description). Friction remains low as long as v≥v_c, and the decrease in velocity arises from the ceasing supply of material from the upstream region.

Figure 11 Velocity and thickness at the five hills marked in Fig. 10. The solid lines refer to the seed points (yellow points in Fig. 10), and the dashed lines refer to the respective peaks in elevation (black points in Fig. 10). The dots depict the collapse of the fluid-like behavior.

After returning to Coulomb friction, however, the particles strongly decelerate. Assuming a horizontal bed and a horizontal surface, this deceleration corresponds to μg≈7.3 m s⁻², causing the particles to come to a complete stop in less than 1.4 s from v=10 m s⁻¹ (the typical order of magnitude shown in Fig. 11a), with a traveling distance of less than 7 m. In turn, the change in thickness at the seed points is not unique. It may increase or decrease, but there is typically little change.

The dashed lines in Fig. 11 show the points that will become the respective peaks in elevation (black points in Fig. 10). By definition, these points return to Coulomb friction later than the respective seed points. Due to their location upstream of the seed points, material accumulates here, causing an increase in thickness until the material comes to rest. While the thickness increases by almost a factor of 2 at peaks 2 and 4, there is almost no increase at hill 5. The majority of the peaks not shown here are also in the range between almost no increase and a factor of 2. Thus, the accumulation of material alone cannot explain the occurrence of more or less isolated hills. There must also be a decrease in thickness in the surrounding area.

Figure 12 illustrates the velocity field around seed point 1. Time is measured relative to the time at which the seed point returns to Coulomb friction (t₀=157.48 s). Flow continues with little change in velocity for several seconds in the region upstream of the seed point, resulting in the accumulation described above. Some slight changes in the direction of velocity indicate that particles try to bypass the obstacle, but this effect is weak. Most of the incoming flux still moves toward the obstacle and causes it to grow rapidly. Downstream, the opposite effect is observed. Particles that returned to Coulomb friction rapidly decelerate, while particles that are already ahead keep moving at a high velocity. This difference leads to a depletion of material downstream of the obstacle and thus to a region with a low deposit thickness.

Figure 12 Velocity field around seed point 1 (shown in Figs. 10 and 11). Time is measured relative to when the seed point returns to Coulomb friction (t₀=157.48 s). The domain is centered on the seed point and is 200×200 m in size. The blue arrows indicate velocity-dependent friction (v≥v_c), and the yellow arrows indicate Coulomb friction (v<v_c). The dots show the points considered in Fig. 13.

The behavior concerning the lateral direction is less systematic. Figure 13 illustrates the change in thickness for six points. The point that will become the peak (black line) is the first to return to Coulomb friction and experiences the strongest increase in thickness. The points shifted by $\mathit{\delta }y=\pm \mathrm{25}$ m follow closely behind but behave differently. One of them (δy=25 m) behaves similarly to the peak, but the southward point ($\mathit{\delta }y=-\mathrm{25}$ m) even decreases in thickness. It is recognized in Fig. 12 that material accumulates somewhat upstream of this point. The point at δy=75 m stays at v≥v_c for almost 10 s longer. During this time span, thickness decreases from more than 8 m to less than 2 m and then further to about 1 m until movement stops completely. This depletion of material is not necessarily related to the dynamics at the hill. Here, it is the result of the ceasing supply from the upstream region, making the layer thinner until it finally comes to rest. In turn, laterally elongated hills or even transverse ridges may occur if the collapse of the flow propagates laterally before the supply ceases.

Figure 13 Thickness at the points marked in Fig. 12. The black line corresponds to the point that becomes the peak. The colored lines refer to points shifted in the north–south direction in steps of 25 m. The black dots indicate the time at which the fluid-like behavior collapses.

Thus, the formation of more or less isolated hills involves two components. The top of the hill is typically formed by the accumulation of material, which starts from a point that returns to Coulomb friction and propagates upstream. As a second contribution, material is depleted in the domain around the hill. This depletion is only partly related to the formation of the hill. The hill immediately suppresses the supply of material to its downstream range, leading to rapid depletion here. In turn, fast flow may continue for a longer time in the upstream and lateral ranges, leading to a depletion of material simply due to the cessation of the incoming flux. Since this contribution depends on what happens in the upstream range rather than on what happens at the hill, it seems to be responsible for the large variation observed in the shape of the hills.

3.4 Influence of the parameters

The results obtained so far were obtained from a combination of parameter values that allows for a relatively long runout in the valley. Figure 14 shows two scenarios that yield a shorter runout. Friction in the fluid-like regime was increased by a factor of 2 (ξ=250 m s⁻² instead of 500 m s⁻²), as shown in Fig. 14a. The other version (Fig. 14b) assumes the original friction but an earlier collapse of the fluid-like behavior (v_c=5 m s⁻¹ instead of 4 m s⁻¹ when h=1 m). The two parameter sets are related to each other through the minimum amount of random kinetic energy required to maintain the fluid-like behavior. According to the considerations of Hergarten (2024c), this amount is proportional to $\frac{v_{\mathrm{c}}^{\mathrm{3}}}{\mathit{\xi }}$. It is thus almost identical for the two parameter sets and twice as high as that for the original scenario.

Figure 14 Deposit thickness in a section of the valley obtained for two different parameter combinations. The value of v_c refers to a thickness of h=1 m. The black outline encloses the region affected by the rock avalanche. The straight lines correspond to the profile considered in Fig. 7.

Runout is more than 1 km shorter in both scenarios than in the original simulation (Fig. 2) and is also shorter than that of the real rock avalanche in the valley. The two scenarios are similar with regard to not only the runout but also the morphology of the deposits. Both versions produce hills, albeit over a shorter section of the valley due to the shorter runout. Some hills even seem to be located at similar positions in both scenarios.

The profiles shown in Fig. 15 confirm the similarity of the morphologies, although the version with the increased friction (Fig. 14a) produces thicker deposits along large parts of the profile compared to the version with an earlier collapse of the fluid-like behavior (Fig. 14b). In some sections of the profile, the morphology even looks similar to that of the original scenario (Fig. 3).

Figure 15 Comparison of the profiles along the straight lines shown in Fig. 14 with the original simulation (Figs. 3 and 7). Elevation is exaggerated by a factor of 4.

These results suggest that the minimum amount of random kinetic energy required to maintain the fluid-like behavior controls how far the rock avalanche travels along the valley. It therefore defines the potential to form toma hills. However, it seems that the formation of hills in the part of the valley covered by rock avalanche deposits is not very sensitive to the parameter values. The crucial condition seems to be that the friction in the fluid-like regime must be much lower than the Coulomb friction.

The spatial resolution also affects the results. Figure 16 was obtained with the original parameter values (ξ=500 m s⁻²; v_c=4 m s⁻¹) but with a grid spacing of 10 m instead of 5 m. The initial conditions were obtained by averaging over tiles of 2×2 points.

Figure 16 Deposit thickness in a section of the valley obtained with 10 m grid spacing instead of 5 m grid spacing. The black outline encloses the region affected by the rock avalanche.

Although the hills are clearly recognizable, their pattern is visually even less similar to the original pattern (Figs. 2 and 3) compared to the patterns shown in Fig. 14. Furthermore, the runout length has increased considerably, although this increase is partly due to a small amount of material traveling along a narrow section of the valley. Overall, the scenario of material flowing along a relatively flat valley seems to be sensitive not only to the model parameters but also to the spatial resolution. While this issue affects estimates of the runout length more than the formation of toma hills, it should be investigated more thoroughly in subsequent studies.

4 Conclusions and outlook

Based on the numerical simulation of a hypothetic rock avalanche, it was found that the modification of the widely used Voellmy rheology proposed by Hergarten (2024c) is able to reproduce the occurrence of distinct hills in the deposits. The key feature of the modified rheology is a sudden drop in friction when transitioning from Coulomb friction to a fluid-like behavior. While this property was originally introduced in order to explain the long runout of rock avalanches without specific assumptions, it turned out to be responsible for the formation of hills as well.

The mechanism of forming hills involves two components, starting with returning to Coulomb friction at one point. The sudden increase in friction decelerates the particles strongly, causing material to accumulate in the domain upstream of the seed point. In turn, particles downstream of the seed point continue to move at high velocities, causing a rapid depletion of material. However, the formation of distinct hills also requires depletion in the lateral region. This happens if flow at low friction levels continues until the supply ceases. In turn, hills become laterally elongated or even form transverse ridges if the fluid-like behavior collapses in the lateral region before the supply ceases.

Since the shape of the simulated hills is quite irregular, a quantitative comparison with real toma hills is difficult. For the valley considered, however, the steepness of the simulated hills is strikingly similar to the steepness of existing toma hills generated by a real rock avalanche. This similarity holds for the average steepness as well as for the slope in the steepest direction. For some of the simulated hills, the slope in the steepest direction is close to the maximum stable slope, which is defined by the coefficient of Coulomb friction. The majority of the hills, however, are less steep.

The simulated hills are more or less isolated in the sense that there are areas with minimal deposit thickness between them. In turn, there are also connections between hills with several meters of deposits. Thus, the hills are not clearly separated in terms of deposit thickness. Although real toma hills often look like islands on an alluvial plain, the partial connections formed by rock avalanche deposits are not inconsistent with observations. Knowledge about the subsurface structure around toma hills is still incomplete. Validating or refuting the model will require more simulations, a more detailed analysis, and more research on real toma hills.

The occurrence of hills that look much like toma hills does not seem to rely on specific parameter values. In terms of the parameters, it seems that the only important factor is that friction in the fluid-like regime is much weaker than Coulomb friction. Given that there is even some similarity in the spatial patterns of the hills for different parameter values, the shape of the valley and the morphology of the valley floor also seem to be important. Finding out which valley morphologies facilitate the formation of toma hills could be the subject of a subsequent study.

So far, the modified Voellmy rheology has been supported mainly by its ability to explain the long runout of rock avalanches and it dependence on volume (Hergarten, 2024c). The finding that it also predicts the formation of hills that look much like toma hills provides further support for this rheology.