Impact of Inherited Foreland Relief on Retro‐Foreland Basin Architecture

We use a Landscape Evolution Model including flexural isostasy to investigate the influence of inherited foreland relief on the stratigraphic evolution of the retro‐foreland domain during mountain building. We show models with four different types of initial relief in the foreland domain: at sea level, elevated (+300 m), a 1 km‐deep and 100 km‐wide foreland basin associated with either a forebulge at sea level or elevated at +300 m. During the first 10 Myr of simulation, the landscape evolution of the foreland is significantly altered by its inherited bathymetry/topography. The impact is then smoothed out once the foreland slope has stabilized and develops a transverse drainage network. Models record a long‐term shallowing‐up mega‐sequence driven by the increase in sediment production rate in the uplifting range and the decrease in the rate of flexural accommodation space creation in the foreland basin. The initial relief of the foreland domain alters the timing of its transition from the under‐filled to the over‐filled phase. An initially deep foreland basin is twice as thick as an initially elevated foreland. It records deep marine deposits while a foreland initially at sea level records thin shallow marine and an elevated foreland records continental deposits. The forebulge is buried by continental deposits in an initially elevated foreland while it is buried by marine sediments in other models. Alluvial fans at the foot of the range are more elevated in initially elevated forelands. We discuss our results of modeled stratigraphic architecture in comparison with the Pyrenean, Alpine and Andean retro‐foreland basins.

Another parameter impacting foreland stratigraphic evolution is the inheritance from rifting before the collision. Mountain ranges often develop in previously rifted domains, as for instance in Tethyan orogenic systems such as the Pyrenees and the Alps (e.g., Beaumont et al., 2000;Desegaulx et al., 1991;Schlunegger et al., 1997;Stampfli & Hochard, 2009;Vacherat et al., 2017). Watts (1992), Stewart and Watts (1997) and Leever et al. (2006) have analyzed the impact of inheritance from a previous rifting event on the evolving EET and evolution of foreland basins. In the northern Pyrenees, Angrand et al. (2018), Desegaulx et al. (1991) and Desegaulx and Brunet (1990) showed that tectonic inheritance, and more specifically, the thermal and crustal structure has a strong impact on foreland basin geometry. Despite these studies, the effect of inherited topography or bathymetry in the proto-foreland domain on stratigraphic architecture of foreland basins have not been addressed yet. Indeed, the impact of structures inherited from rifting on orogenic deformation has been widely studied (e.g., Erdos et al., 2014;Wolf et al., 2021), but the implications for the stratigraphic architecture of foreland basins remains to be constrained. Remnants of a previous extensional phase in the relief of an initial foreland domain can potentially significantly impact its capacity to trap sediments produced in the mountain range and in doing so, the shape of the foreland basin, the duration of the underfilled/overfilled phases, as well as the associated paleo-environments. The aim of this work is to explore this effect and compare its magnitude to other controlling parameters such as mountain range uplift rate, EET of the foreland lithosphere, and erosion/transport/sedimentation efficiency.
To do this, we use a Landscape Evolution Model (LEM) taking into account flexural isostasy and both marine and continental sedimentary processes (FastScape S2S; Yuan, Braun, Guerit, Rouby, & Cordonnier, 2019;Yuan, Braun, Guerit, Simon, et al., 2019). The LEM allows assessing, in 3D, the relationships between flexural isostasy, landscape evolution and stratigraphic architecture of the foreland basin ( Figure 1). We focus on the stratigraphic architecture of the retro-wedge foredeep, between the frontal tip of the orogenic wedge and the forebulge (DeCelles & Giles, 1996). This allows us to approximate the syn-orogenic evolution using only vertical motion (uplift and flexural isostasy) as retro-wedge systems of small to intermediate size orogens are less affected by horizontal advection related to thrusting (Grool et al., 2018;Naylor & Sinclair, 2008;Willett et al., 1993;Wolf et al., 2021). We compare our simple cylindrical setups and generic approach to the Pyrenean, Andean, and Alpine retro-foreland systems, in order to discuss the potential effect of inherited topography and/or bathymetry on the evolution of these foreland basins.

10.1029/2022JB024967
3 of 20 to a distal open marine domain (200 × 400 km; Figure 1; Table 1). The foreland domain includes the foreland basin and the forebulge controlled by flexural isostasy (Figure 1) which does not migrate. The distal open marine domain is a boundary condition inspired by the Pyrenean system (Bernard & Sinclair, 2022;Ortiz et al., 2020) although we did not simulate the non-cylindrical configuration of the natural case. The distal marine domain has been included to avoid border effects which would affect the dynamics of continental deposition in our models. This distal domain is included, but not discussed. We present four models with varying initial topography and bathymetry in the foreland ( The initial bathymetry at the foot of the range in M3 and M4 represents pre-existing rift related topography that can be encountered in natural orogenic systems such as the Pyrenees (e.g., Desegaulx et al., 1991). The associated elevated forebulge is consistent with the paleogeographic reconstruction of the Northern Pyrenean system ∼55 Myr ago by Vacherat et al. (2017). The initially elevated foreland domain in models M2 and M4 represents stable Phanerozoic continents that have an average elevation of ∼400 ± 400 m (e.g., Theunissen et al., 2022). Model duration (25 Myr), is consistent with the time span between the main phase of Pyrenean topography emergence (∼55 Myr) and the syn-to post-orogenic transition ∼27 Myr (e.g., Curry et al., 2019;Vacherat et al., 2017). Although inspired by the Pyrenees, this simple configuration can also be compared to other orogenic systems such as the Alpes or the Andes to discuss their retro-forelands sedimentary dynamics. To initiate rivers grading toward the foreland domain, we impose a gentle initial tilt of the uplifted domain (α = 0.076°; Figure 2).
For the analysis of the stratigraphic architecture, we use the definition of Catuneanu (2004Catuneanu ( , 2017 and Sinclair and Allen (1992) of the underfilled, filled, and overfilled stages in the evolution of the foreland basin in which depositional processes are dominated by deep marine, shallow marine, or fluvial sedimentation, respectively.

Reference Model M1
In reference model M1, the mountain range grows gradually to an average topography of 1.7 km elevation at 25 Myr (Figures 3, 4a, and 5a). The basement of the flexural foreland basin subsides progressively under the increasing load of the mountain range topography and of the sediments to a maximum depth of 2.6 km (Figures 3 and 5c). Initially, the depositional environments are shallow marine and the forebulge is partly submerged (Figures 3a and 4a; Movie S1). Part of the sediments produced by erosion of the mountain belt fills the flexural foreland basin while the remainder is exported to the marine domain (Figures 3a and 4a). After a 3.7 Myr, alluvial deposits and initially isolated and progressively coalescing alluvial fans (i.e., sediment deposited at a slope >0.4°; Bull, 1964;Milana & Ruzycki, 1999) form at the foot of the mountain range (Figures 3b and 4a). This transition is slightly diachronous laterally, due to relief heterogeneities depending on the local relief of the mountain range and previously deposited sediments (Movie S1). Afterward, the shoreline and the continental deposits progressively migrate from the foot of the mountain range toward the forebulge as the foreland basin evolves toward marine to continental environments (continentalization; Figures 3b and 4a). At the same time, alluvial fans at the foot of the range migrate away and toward the mountain range. These oscillations are driven by the competition between local erosion and the space available for deposition, which is controlled by the deposition of previous fans, local relief, and individual drainage dynamics (Movie S1). These short-term oscillations do not alter the general migration of continental deposits across the foreland basin. At 15 Myr, continental deposits reach the forebulge and at 25 Myr, the foreland basin and the forebulge are entirely continentalized (Figures 3c, 3d, and 4a; Movie S1).

Models M2 to M4 With Inherited Topography/Bathymetry in the Foreland Domain
Models with inherited topography and/or bathymetry in the foreland domain (M2 to M4) follow a first order trend common with the reference model M1: initial building of the mountain range topography, development of the flexural foreland basin and forebulge, formation and coalescence of alluvial fans at the foot of the mountain range, and progressive migration of continental deposits across of the foreland domain ( Figure 4, and Figures S1-S3 in Supporting Information S1). Inherited topography and/or bathymetry in the foreland domain do nevertheless have a significant impact on the timing and depositional environments of this trend ( Figure 4, and Figures S1-S3 in Supporting Information S1).
The initially elevated foreland domain of model M2 (+300 m; Figure 2b) records only continental sediments (alluvial deposits and fans), allows for less material to be trapped at the foot of the range and is rapidly incised by a regressive erosion that connects the mountain range to the open marine domain (within 1 Myr; Figure 4b and Figure S1 in Supporting Information S1; Movie S2). These river networks, not only export sediments produced in the mountain range, but also remobilize sediments previously stored in the foreland basin. In contrast with the reference model M1, this open marine domain is entirely filled and continentalized after 25 Myr suggesting that more sediments were exported (Figure 4b and Figure S1 in Supporting Information S1).
In model M3, sediments produced in the mountain range are trapped in the initially water filled basin at the foot of the range (1,000 m deep; Figure 2c) in deep marine depositional environments during the first 7 Myr (Figure 4c and Figure S2 in Supporting Information S1; Movie S3). The foreland basin progressively fills up to shallow marine deposits and the forebulge is progressively submerged and buried by shallow marine deposits. First fluvial deposits emplaced at the foot of the range at 3 Myr, followed by alluvial fans that start accumulating by 6 Myr (Figure 4c and Figure S2 in Supporting Information S1; Movie S3). The shoreline and continental deposits migrate across the foreland to reach the forebulge that is starting to be buried by continental sediments by 13 Myr (Figure 4c and Figure S2 in Supporting Information S1; Movie S3).
In model M4, combining an initially deep water filled basin and an elevated foreland domain, the foreland also traps deep marine sediments, and the initially elevated foreland is rapidly incised (within 1 Myr; Figure 4d and Figure S3 in Supporting Information S1; Movie S4). Sediments produced by erosion of the forebulge area contribute to the filling of the foreland basin in addition to the sediments produced by erosion of the range. Similar to model M2, once the initial bathymetry is filled to shallow marine depositional environments and then  Stock and Montgomery (1999). h Parameters from Braun and Willett (2013). i Parameters from Jordan and Flemings (1991). j Parameters from Rouby et al. (2013). k Parameters from Yuan, Braun, Guerit, Simon, et al. (2019). l Parameters from Simon et al. (2022). continentalized (∼4.6 Myr), the whole foreland domain is incised by regressive erosion remobilizing previously deposited sediments.

General Characteristics of Mountain Range and Foreland Basin Evolution
We evaluate the evolution of mean elevation, mean erosion rate, foreland basin depth, and sediment volume in the models ( Figure 5). The evolution of the mean elevation of the uplifted domain is very similar in the four models and shows a progressive build-up to ∼1.7 km after 25 Myr, without reaching steady state ( Figure 5a). Associated mean erosion rates in the mountain range follow a similar overall increase to 3.5-4.0 × 10 −4 m/yr at 25 Myr ( Figure 5b). During this increase however, all models undergo abrupt drops in mean erosion rates (ca. two-fold decrease; Figure 5b). The timing of the drops in erosion rate varies from one model to the other (5.2, 4.9, and 3 Myr in Models 1, 2, and 3 respectively). Model M4 shows a more complex behavior with a first drop at 4.6 Myr and a second one at 11.9 Myr associated with a few oscillations. After the drops, all models return to a trend of increasing mean erosion rates over time ( Figure 5b). This particular behavior is further discussed below. Figure S4 in Supporting Information S1 provides a top view of the erosion and deposition rates above sea level through time.
The maximum basement depths of the foreland basins of M1-M4 exhibit similar deepening trends but reach different final depths at 25 Myr (e.g., 2.6 km for M1, 2.3 km for M2, 3.6 km for M3 and 3.4 km for M4; Figure 5c). The total volume of sediments produced in the mountain range is similar in the four models (4-4.5 × 10 14 m 3 ). However, the volume of sediment accumulated in the foreland basin is quite different between the models (1.30 × 10 14 , 0.95 × 10 14 , 2.20 × 10 14 and 1.80 × 10 14 m 3 for M1, M2, M3, and M4 respectively; Figure 5d). This is mirrored by different proportions of sediments exported to the open marine domain ( Figure 5d).
Models M1-M4 are based on identical uplift rates of the mountain range, erodibility, and EET. To assess the impact of the initial bathymetry/topography of the foreland in comparison to one of these parameters, we also performed a sensitivity analysis of model M1 to varying uplift rate, erodibility, and EET. Supplementary model SM1 shows that increasing uplift rates (0.1-1 mm/yr) increases the mean topography of the mountain range (100-1,970 m respectively), the topographic load and the flexural isostatic response of the foreland and ultimately the thickness the foreland basin (maximum depth from 0 to 2,500 m respectively; Figure S9 and Table S1 in Supporting Information S1). Supplementary model SM4 shows that increasing erodibility (0.5-9 × 10 −5 m 0.2 /yr), decreases the mountain range mean topography (2,620-360 m respectively) and the associated flexure in the foreland basin (maximum depth from 3,850 to 450 m respectively; Figure  S10 and Table S1 in Supporting Information S1). Finally, supplementary model SM5 shows that increasing the EET (5-25 km), decreases the amplitude of the foreland basin flexure, and the thickness of the foreland basin infill (maximum thickness from 3,850 to 450 m respectively; Figure S11 and Table S1 in Supporting Information S1). This analysis shows that the impact of the initial bathymetry/topography of the foreland geometry is as significant as that of uplift rate, erodibility, and EET.

Foreland Basin Stratigraphic Architecture
For each model, we show the stratigraphic architecture of the foreland basin along a cross-section as well as the corresponding Wheeler diagram of the depositional bathymetry/elevation through time ( Figure 6; Sections location   in Figure 4). We highlight the transition from alluvial fan to fluvial plain deposits for sediments with a depositional slope >0.4° ( Figure 6 and Figure S5 in Supporting Information S1; Bull, 1964;Milana & Ruzycki, 1999).
The foreland basin in reference model M1 has a maximum thickness of 2.6 km at the mountain front ( Figure 6a). It shows continuous accumulation and is progressively filled up. Thus, it shows a shallowing trend, with shallow marine depositional environments at the base (water depth <100 m) evolving to continental depositional environments that record progressively increasing elevations ( The foreland basin in model M2 is thinner than in reference model M1 (2 km) and includes only continental deposits (Figures 6b and 6f). It is filled by continental deposits at progressively higher elevation (up to 700 m). The continental sediments migrate across the foreland basin, onlap the forebulge and bury it by 11 Myr (i.e., slightly before model M1; Figures 6b and 6f). The continental foreland domain shows several local incisions, as fluvial incision (channels of a few kilometers) or wider erosion areas (∼80 km) develop, particularly above the buried forebulge, and remobilize previously deposited sediments. Alluvial fans propagate back and forth up to ∼90 km within the foreland basin ( Figure 6j).
The foreland basin in model M3 is ultimately significantly thicker than in reference model M1 (4.5 km) and includes deeper marine deposits (water depths >300 m; Figure 6c). As in model M1, it is progressively filled up in a shallowing trend, but from deep marine depositional environments at the base evolving to shallow marine and continental depositional environments that record progressively increasing elevations (up to 500 m; Figures 6c  and 6g). The deep marine sediments progressively onlap the forebulge and shallow marine deposits bury it by 10 Myr (i.e., earlier than in model M1). Subsequently, similarly to model M1, the shoreline propagates across the foreland domain although alluvial fans propagate back and forth further within the foreland basin (up to ∼80 km; Figure 6k).
The foreland basin in model M4 is as thick as in model M3 (4.4 km) and includes deep marine deposits as well (water depths >300 m; Figure 6d). As in model M3, it is progressively filled up in a shallowing trend, with deep marine depositional environments at the base evolving to shallow marine and continental deposits that however record ultimately higher elevations (>700 m; Figures 6d and 6h). As in model M1, the shoreline propagates across the foreland domain. However, the forebulge remains above sea level throughout the foreland basin infill, even feeding it with sediments resulting from its erosion. It is buried by continental deposits by 12 Myr. The continental foreland domain shows several local incisions, as fluvial incision develops and remobilizes previously deposited sediments. Alluvial fans propagate back and forth by more than 100 km in the foreland basin ( Figure 6l).

Erosion and Accumulation Dynamics
Our models show specific erosion and sediment accumulation features. In reference model M1, erosion rates in the mountain range reduce sharply at 5.2 Myr and increase again steadily afterward (Figures 5b and 7). For M1, this drop of erosion rate is coeval with the coalescence of alluvial fans at the foot of the mountain range (Figures 7b and 8a). Models M2-M4 exhibit similar behavior with one or more drops in erosion rate that are also correlated to changes in depositional environments (transition from marine to continental fluvial plain deposits) or alluvial fan dynamics (transition from continental fluvial plain to alluvial fan deposits; Figure 5b, Figures  S6-S8 in Supporting Information S1). Figure 8 shows the co-evolution of depositional environments and bathymetry/elevation at the foot of the mountain range. In models M1, M3, and M4, transition from marine to continental depositional environments occurs between 3.0 and 4.6 Myr. In model M1, M2, and M3, alluvial fan build-up occurs between 4.9 and 6.0 Myr. Model M4 presents specific features in comparison to other models. Transition from marine to continental depositional environments corresponds to a first alluvial fan build-up (i.e., without preceding fluvial plain deposits) and a second phase of alluvial fan build-up occurs at 11.9 Myr (Figure 8d). The shoreline migration rates across the foreland are 23 and 17 km/Myr for M1 and M3, respectively (Figures 8a and 8c). For all models, the maximum elevation of the alluvial fan varies from 600 to 800 m at 25 Myr ( Figure 8). It is higher for M2 and M4 (elevated foreland) than M1 and M3 (foreland at sea level). Interestingly, drops in erosion rates in the mountain range are coeval with the transition from marine to continental depositional environments or with alluvial fan coalescence in the foreland basin (Figures 5b, 7, and 8; Figures S6-S8 in Supporting Information S1).

Stratigraphic Trends of the Foreland Domain Common to All Models
In terms of stratigraphic architecture, all models record a long-term prograding mega-sequence that is characteristic of foreland basin stratigraphic architectures (Figure 9a; e.g., DeCelles & Giles, 1996). During the initial stage of the simulation, the topographic load of the rising mountain range creates accommodation in the foreland basin by flexural isostasy and allows the storage of sediments at the foot of the mountain range. For a foreland domain initially at sea level (model M1), these sediments are initially deposited under shallow marine environments (Figure 9a). After 3.6 Myr, alluvial and alluvial fan deposits at the foot of the range (Figure 8a) mark the beginning of the prograding mega-sequence and of the transition of the foreland basin to the continentalization phase. Then, the whole depositional profile (shoreline, alluvial, alluvial fans) migrates away from the mountain range recording the ongoing prograding (coarsening-up) mega-sequence (Figures 6 and 9). All models (M2-M4, SM1-SM4) display this long-term trend. In all experiments, the sediment supply [S] produced by erosion of the uplifting mountain range is increasing over the 25 Myr of the simulation (Figure 5d). In parallel, the subsidence of the foreland basin basement is slowing down, hence reducing the accommodation space [A] creation in the foreland basin (Figure 5c). This reduction in accommodation space creation is driven by the slowing down of topographic loading over the time (Figure 5a). The drainage development in the mountain range and the associated erosion are progressively increasing in efficiency (as indicated by the volume of sediment produced; Figure 5d), reducing the rate of the mean topographic rise toward steady state (e.g., Babault et al., 2005;Carretier & Lucazeau, 2005;Tucker & van Der Beek, 2013), although it is not reached in 25 Myr of simulation ( Figure 5a). As a result of these two coeval trends, the volume of sediments initially produced by erosion of the uplifting range is not sufficient to fill the space available for sedimentation. Then, progressively, the volume of sediments becomes equal and higher than the accommodation creation (overfilled phase). Although the proportion of sediments by-passing the foreland basin and exported to the open marine domain increases through time, the volume of sediment preserved in the foreland basin is increasing as well ( Figure 5d) and drives the long-term prograding mega-sequence.
Another feature common to all models is that, once the foreland domain is fully continentalized (i.e., the long term progradation reaches the static forebulge), the sedimentary load is not limited to the foreland basin but is distributed over the entire accumulation area (foreland basin and forebulge; Figures 3, 6, and 9). Consequently, the impact of the sedimentary load in terms of differential subsidence/uplift between the foreland basin and the forebulge decreases after the continentalization and the burial of the forebulge (between 10 and 13 Myr in the models shown here). The accommodation in the foreland domain is then mostly controlled by the rise of the mountain and the distal marine base level. Over that period, the foreland basin progressively widens (Movies S1, S3, and S4). Leever et al. (2006) proposed that the migration of the orogenic load by the propagation of the thrust front away from the range could widen the flexural basin. We show here that, with a static orogenic load and thrust front, the sediment distribution can also induce a widening of the foreland basin even with a steady uplift of the range.

Influence of the Initial Elevation on the Stratigraphic Trend of the Foreland Domain
During the first 10 Myr, the landscape evolution of models M1-M4 is significantly different as a result of their inherited foreland domain bathymetry/ topography ( Figure 4). Interestingly, the influence of the initial relief largely disappears in the landscape of all models after ∼10-13 Myr (Figure 4), after the foreland domain has been continentalized and its slope stabilized, that is, after the long-term prograding trend reaches the static forebulge. Afterward, all models show very similar landscape evolution with a continental foreland domain developing a transverse hydrographic network (Figure 4). Nevertheless, the initially different landscape evolution results in major differences in the stratigraphic architecture of the foreland basin.
First, the ultimate thickness of the foreland basin is very different. The foreland domain initially at sea level preserves up to 2.6 km of sediments, while the case with an elevated foreland domain preserves only 2 km and the case with an initially deep foreland domain preserves up to 4.5 km (Figures 6 and 9). As all models show a similar sediment production (Figure 5d), this difference results from the flexural response combined with the initial capacity of the foreland domain to store sediments (directly controlled by its initial elevation). The additional load of the sediments trapped in the initial deep basin in the foreland domain (models M3 and M4) amplifies the flexural subsidence to create a thicker basin with respect to reference model M1 (Figures 5c  and 5c). Conversely, the reduced sediment load in model M2 with an initially elevated foreland domain dampens the flexural subsidence and produces a thinner foreland basin with respect to model M1 (Figures 5a, 5c, and 5d). Flemings and Jordan (1989) have proposed that an increase in erosion or deposition efficiency, a higher EET or a slower thrust rate advance can result in a thicker foreland basin. We show here that an inherited bathymetry, such as a rift remnant, is another mechanism to produce a thicker foreland basin. Our models suggest that the deeper the initial rift remnant, the thicker the foreland basin deposits.
Second, the foreland domain initially at sea level preserves only 0.5 km of (shallow) marine sediments, while the elevated foreland domain preserves only continental deposits and the initially deep foreland domain up to 2.5 km of (deep) marine sediments (Figures 6 and 9). This is again the direct result of the foreland domain initial elevation impacting the depositional environments at which the long-term prograding mega-sequence initiates: deep marine, shallow marine, or continental. Our models show that, for the parameters used, an initially deep basin in the foreland domain is required to preserve a significant proportion of marine deposits in the foreland basin. Fourth, the forebulge is buried under continental sediments in models with an initially elevated foreland (with or without a deep inherited basin; models M2 and M4) while it is buried by marine sediments in cases with a foreland domain initially at sea level (with or without a deep inherited basin; models M1 and M3; Figures 6  and 9). The flexural subsidence in the foreland domain is sufficient to submerge the forebulge in experiments with a foreland domain at sea level but not in the experiment with an initially elevated foreland that remains above sea level throughout the experiment (Figures 4, 6 and 9). In addition, in an initially elevated foreland domain (with or without a deep inherited basin; models M2 and M4), continental and alluvial fan deposits reach ultimately higher elevations than in a foreland initially at sea level (with or without a deep inherited basin; models M1 and M3; Figures 6 and 9). Indeed, as the forebulge is acting as the base level once the foreland basin is filled up and continentalized, the more elevated the forebulge, the more elevated the upstream continental deposits (over 700 m for models M2 and M4 and about 500 m for models M1 and M3). The higher elevation of the alluvial fans also allows for their spreading further away from the mountain range (up to 111 km for model M2 and M4 and 40-80 km for models M1 and M2).
In summary, the occurrence of an inherited topography/bathymetry in the foreland domain does not alter the longterm prograding (shallowing-up) trend of the foreland basin. Indeed, for a constantly uplifting range, constant erodibility, for a given EET, the mean elevation of the range converges toward an equilibrium state and results in a coeval decrease in the rate of sediment production ( Figure 5b) and of flexural accommodation space creation in the foreland basin. In all experiments, the decay in topographic building rates (Figure 5a) causes the attenuation of flexural subsidence in the foreland, which is in favor of the sediment volume that progressively fills up the accommodation volume and produces a prograding mega-sequence ( Figure 9). Nonetheless, the decrease of the topographic building rates, associated with inherited topography bathymetry results in different rates of decay of sediment production in the range and accommodation creation in the foreland domain ( Figure 5). This produces either deep marine, shallow marine or continental initial depositional environments in the foreland domain. It also impacts rates of progradation in the timing of the transition from marine to continental conditions. In an initially elevated foreland, the transition from marine to continental conditions occurs earlier and the progradation is faster than in initially deep foreland (Figure 9).

Feedback of the Foreland Domain Dynamics on Erosion Rates in the Mountain Range
We show above how inherited topography and/or bathymetry in the foreland domain impacts its subsidence and accumulation history. Our models also show that landscape dynamics in the foreland domain provide a feedback affecting the erosion dynamics of the mountain range. The abrupt drops in erosion rates in the uplifted domain (Figures 5b and 7) are synchronous with changes in the depositional systems at the foot of the mountain range. They systematically correspond to a transition from marine to continental depositional environments or from fluvial to alluvial fan deposits (Figures 7 and 8, Figures S6-S8 in Supporting Information S1). The continentalization of the foreland domain and the build-up and coalescence of alluvial fans, are associated with a raise of the base-level at the foot of the mountain range. Indeed, deposition at the orogenic piedmont will shift the elevation of the drainage basins upward since relief denudates at a lower rate than the uplift rate (Babault et al., 2005). This results in the increase of the absolute elevation of the topography by an amount equal to the mean elevation of the alluvial fan apex, which defines the base level of the uplifting relief (Babault et al., 2005). This base level rise is responsible for reducing the erosive potential for upstream areas and is thus responsible for the transient drops of erosion rates observed in the mountain range (Figures 5b and 7; Babault et al., 2005;Carretier & Lucazeau, 2005). After a while, the hydrographic network returns to its previous base level and erosion rates in the mountain range gradually return to similar but lower trends (Figure 5b). This autogenic feedback has previously been documented using both analog (Babault et al., 2005) and numerical modeling studies (Carretier & Lucazeau, 2005). This suggests that erosion rates in mountain ranges driven by climatic or tectonic forcing can also be modulated by downstream sedimentary dynamics. For example, the endorheic phase of the Ebro basin in the Eo-Oligo-Miocene (∼37 and 16.5 Myr) resulted in substantial accumulation of continental deposits among others at high elevation at the foot of the southern Pyrenean range (e.g., Babault et al., 2005;Garcia-Castellanos et al., 2003). This period is characterized by a lowering in erosion rates in the mountain range that has been interpreted as resulting from the rise of the regional base-level (Babault et al., 2005;Garcia-Castellanos et al., 2003). The high-frequency oscillations in erosion rates observed in models M1 and M4 correspond to rapid coalescence and dispersal events of the alluvial fans at the foot of the range inducing transient rise and fall of the local base level (<500 kyr; e.g., Figure 5b). Indeed, the alluvial fans coalescence is not definitive as fans coalesce and disunite during a few time steps (Movies S1 and S4). However, these rapid oscillations do not impact the long-term erosion dynamics of the mountain range (Figure 5b).

Comparison With Natural Retro-Foreland Systems
We consider the implications of our model results for the northern retro-foreland system of the Pyrenees (Figure 10). The northern Pyrenees and the Aquitaine basin-Bay of Biscay system is a classic example of retro-wedge flexural foreland basin (Angrand et al., 2018;Bernard et al., 2019;Ortiz et al., 2020). Our set-up does not include several features of the Pyrenean system such as: horizontal displacement of thrusts, postrift thermal subsidence (Angrand et al., 2018;Vacherat et al., 2014), basement heterogeneities in the retro-foreland basin (Angrand et al., 2018), geological and geometric complexities during mountain building (Vacherat et al., 2017), lateral variations in exhumation and uplift of the mountain range Fillon & van der Beek, 2012;Fitzgerald et al., 1999), and the elbowed geometry of the North-Pyrenees-Aquitaine-Bay of Biscay system. Nevertheless, our simplified models exhibit first order features useful to understand the Pyrenean retro-foreland basins systems. The mean mountain range elevation after 25 Myr in the order of 1.5-2 km ( Figure 5a) is consistent with the reconstructed mean elevation of the Pyrenean mountain range at the end of the syn-orogenic phase (e.g., Curry et al., 2019;Huyghe et al., 2012). The northern Pyrenean retro-foreland basin records a prograding mega-sequence similar to our models (e.g., Ford et al., 2016;Ortiz et al., 2020;Rougier et al., 2016). The Pyrenees developed by inversion of an inherited rifted domain (e.g., Desegaulx et al., 1991;Vacherat et al., 2017). The inherited pre-orogenic rift formed an initially deep foreland domain in the western sector of the Aquitaine foreland similar to models including an initial bathymetry at the foot of the range (M3 and M4). The maximum total subsidence at the thickest part of the Pyrenean retro-foreland (Central Pyrenees; close to ECORS line; Roure et al., 1989) ranges from 4 to 5 km-depth . These basement depths are consistent with our models including an initial bathymetry at the foot of the range (models M3 and M4; Figures 6 and 10) while it is shallower in models without (M1 and M2; Figures 6a, 10b, and 10c). The basement of the model M3 foreland basin is about 1 km deeper than in the Pyrenean case ( Figure 10c), but its first-order stratigraphic architecture is consistent with the main trends of the present-day Pyrenean retro-foreland basin (Figures 10b and 10c). It shows deep initial depositional environments in the foreland basin similarly to the northern Pyrenean flysch, deposited in the late Cretaceous, during early convergence Puigdefabregas & Souquet, 1986). It includes a significant section of marine sedimentary deposits like the one preserved in the Pyrenean retro-foreland basin (Serrano et al., 2006; Figure 10b). It also includes marine sedimentary deposits onlapping and burying the forebulge as in the Pyrenean retro-foreland basin (Serrano et al., 2006; Figure 10b). Our models show that an inherited bathymetry in the foreland basin associated with a forebulge area initially at sea level is critical to preserve sediments deposited in deep-marine depositional environments, to produce a thick marine sedimentary section in the retro-foreland, and to create a forebulge onlapped and buried by marine sediments (Figures 6c, 9, and 10). The configuration with a deep rift remnant in the foreland basin area and a low elevation of the forebulge area in model M3 is consistent with the paleogeographic reconstruction of Vacherat et al. (2017).
The northern Andean sediment routing system (∼3,500 km for the Amazonian drainage area; Bajolet et al., 2022) is significantly longer than the northern Pyrenean one (∼700 km) and the Andes has a significantly higher elevation. It is, however, interesting that the only regions in the Andean retro-foreland basin where marine, and especially deep marine, sediments are  Angrand et al. (2018) and Serrano et al. (2006)). The inset shows Pyrenees and its associated Aquitaine retro-foreland basin location (red square) at the scale of Western Europe. (b) Cross-section of the Pyrenean retro-foreland stratigraphy (modified after Serrano et al., 2006). The transition from marine to continental depositional environment (Maastrichtian) is deduced from Rougier et al. (2016). Other marine deposits related to marine incursions later in the stratigraphic sequence (Ypresian) are not represented because of their limited thickness (<50 m; Rougier et al., 2016). (c) Cross-sections of foreland basin stratigraphic architectures and basement depth in models M1-M3. We plotted the present-day smoothed base of the Pyrenean retro-foreland from Angrand et al. (2018). AB: Aquitaine basin, outlined in black (Ortiz et al., 2020). NPFT: Northern Pyrenean Frontal Thrust. preserved (e.g., Llanos basin (Colombia), Oriente basin (Ecuador), Ucayali basin (Peru)) correspond to former retro-arc basins (Horton, 2018). These retro-arc basins could have been associated to deep foreland areas like the one tested in M3 and M4. Regions located in between these basins, in the northern Andes area present a classical coarsening-up prograding sequence from shallow marine to continental deposits (Horton, 2018) more comparable to M1 (Figure 6a). In areas where the Andean foreland domain was initially elevated (outcropping basement of South American shields), the forebulge is mainly buried by continental sediments (Bajolet et al., 2022) as predicted by models M2 and M4 (Figure 6d).
In the Alpine retro-foreland, Paleogene to Miocene sedimentation in the basin is controlled by preserved Mesozoic extension related geometries (Turrini et al., 2016). Where initial accommodation space has been preserved from the rifting phase, the foreland basin is thicker preserving deep marine sediments (Turrini et al., 2016) as predicted by M3 and M4 (Figure 6c).
Our cylindrical set-up does not include several features of these natural retro-foreland systems (e.g., longitudinal drainage, horizontal displacement of thrusts, post-rift thermal subsidence, basement heterogeneities, lateral variations in exhumation and uplift of the mountain range, among others). However, it provides insight into the impact of several usually underestimated factors specifically the initial bathymetry and relief in the foreland domain that can provide an alternative explanation for some of the first order observations of the stratigraphic architecture of retro-foreland basins (thickness, preservation of deep marine deposits in the foreland basin and forebulge area, timing of the marine to continental transition).

Model Limitations
Thrust front propagation affects the syn-orogenic dynamics of foreland basins (Simpson, 2006), in particular by causing foreland migration of facies belts remobilizing previously deposited sediments at the foot of the mountain range as well as by inducing retrogradation phases in the foreland basin at the onset of thrusting events (Flemings & Jordan, 1990). The effects of thrust propagation are significant in pro-foreland systems where thrust front migration can exceed 100 km as for instance in the southern Pyrenean pro-wedge (Grool et al., 2018). Our models do not include horizontal deformation and cannot be used as an analog for pro-wedge systems. However, they are useful for understanding retro-foreland systems in which the maximum propagation of the deformation front is limited and less than 100 km. The northern Pyrenees are characterized by shortening of about 20 km (Grool et al., 2018). Naylor and Sinclair (2008) suggested that, in retro-foreland basins, the stratigraphic architecture is mostly controlled by the load of mountain range topography and the associated flexural isostatic subsidence of the foreland, whereas horizontal thrust propagation plays a subordinate role.
Natural examples of mountain range-foreland systems may also display lateral variations in the degree of shortening, amount of erosion, associated sediment delivery to the foreland or the location of the distal open marine domain. In the Pyrenees, basement depth varies from 1 to 3 km in the east to >5 km in the west. These variations have been mainly related to variations in extensional inheritance in the foreland (Angrand et al., 2018). The diachronous onset of the exhumation and topographic build-up from east to west (Vacherat et al., 2017) is also responsible for along strike variations in sediment supply which impact foreland basin filling and stratigraphic architecture (Michael et al., 2014;Ortiz et al., 2022;Verges, 2007). Our cylindrical modeling setup does not allow testing for these lateral variations, which may be investigated in future work using a non-cylindrical model setup. Nonetheless, Fastscape S2S does simulate in three dimensions depositional systems and lateral variations of deltas or alluvial fans at kilometric scale (Figure 4; Movies S1-S4). Note that these local sediment migrations along strike do not affect the long-term trends in sedimentary filling and stratigraphic architecture.
For sake of simplicity, in our models, global sea level, precipitation rate, uplift and sediment transport coefficient (K f ) are constant through time and homogenous in space. Furthermore, we do not include a multi-grain size distribution of the marine deposition and marine diffusion of sediments (e.g., sand vs. silt; Rouby et al., 2013;Yuan, Braun, Guerit, Simon, et al., 2019). Investigation of tectonic or climate-driven variations of sediment supply and the detailed stratigraphic architecture of the open-marine domain is, however, beyond the scope of our study as we focus on the long-term stratigraphic architecture of the foreland basin.

Conclusions
We investigate the influence of inherited foreland relief on the stratigraphic evolution of the foreland domain during the building of a mountain range using a LEM that couples continental and marine surface processes with flexural isostasy. We show models with four different reliefs in the foreland domain: one initially at sea level, one initially at +300 m (continental foreland), one with a pre-existing 1 km-deep and 100 km-wide basin associated with either a forebulge area at sea level or elevated at +300 m.
Our models show that, during the first 10-13 Myr of model simulation, the landscape evolution of the foreland is significantly affected by inherited bathymetry/topography. However, the impact of the initial relief disappears in the landscape evolution once the foreland slope has stabilized and develops a transverse drainage network.
All models record a long-term prograding (coarsening-up and/or shallowing-up) trend of the foreland domain showing that, throughout the simulation, sediment production in the uplifting range is greater than the flexural creation of accommodation space in the foreland basin. Once the foreland domain is fully continentalized and the sedimentary load is distributed over the entire accumulation area (foreland basin and forebulge), the contribution of the differential subsidence between the foreland basin and the forebulge in the stratigraphic evolution decreases (after 10-13 Myr in the models shown here).
Models with different initial topography and/or bathymetry result in major differences in stratigraphic architectures of the foreland basin.
1. Initially deep basins lead to significantly thicker foreland deposits compared to a scenario with an initially elevated foreland. 2. An initially deep basin results in deposition and preservation of thick deep marine deposits while a foreland initially at sea level records only thin shallow marine deposits and the elevated foreland case only continental deposits. 3. The forebulge is buried under continental sediments in an initially elevated foreland (with or without a deep foreland basin) while it is buried by marine sediments in a foreland domain initially at sea level. 4. The elevation of alluvial fans at the foot of the range is higher (up to 200 m) in initially elevated foreland (with or without a deep foreland basin) than a foreland domain initially at sea level. 5. The initial topography and/or bathymetry of the foreland domain alters the timing of the transition from marine to continental phase: it occurs up to 5 Myr earlier in an initially elevated foreland compared with an initially deep foreland and the progradation rate is up to 35% faster.
All the models exhibit alluvial deposits and/or alluvial fan coalescence at the foot of the mountain belt that induces transient drops of erosion rates in the range by raising the local base level, showing how the dynamics of the depositional system at the foot of the mountain range may exert feedback on the erosion dynamics in the mountain range.
Comparison with the Pyrenean, Alpine and Andean retro-foreland basins shows that inherited bathymetry related to pre-orogenic rift structure, allows the deposition of a significant amount of syn-orogenic deep marine deposits and that a forebulge initially at sea level can be onlapped and buried under marine deposits. Although our cylindrical set-up does not include several features of these natural retro-forelands, it can provide usually underestimated factors related to the initial bathymetry and relief in the foreland domain to explain first order observations of the stratigraphic architecture of the retro-foreland basins (e.g., thickness, preservation of deep marine depositional environments in the foreland basin and forebulge area, timing of the marine to continental transition).