Interactions between propagating rifts and linear weaknesses in the lower crust Open Access

The continental lithosphere is characterized by rheological heterogeneities and mechanical anisotropies inherited from previous tectonic events (e.g., Audet and Bürgmann, 2011). Rheological heterogeneities are commonly associated with past orogenic or intraplate processes (i.e., collisional belts, shear zones, rifts, thermal weakening) that result in areas of the lithosphere with reduced strength (e.g., Vauchez et al., 1997; Corti, 2008; Gueydan et al., 2008). Mechanical anisotropies arise as a consequence of tectonic fabrics (i.e., foliations, lineations, crystallographic preferred orientations, reduced grain size) that define a mechanically weaker region than the surrounding lithosphere. It has long been recognized that they control the pattern, propagation and overall structural and sedimentological evolution of continental rifts (e.g., Rosendahl, 1987; Smith and Mosley, 1993; Vauchez et al., 1998; Will and Frimmel, 2017). The nature of the rheological heterogeneities and/or mechanical anisotropies, whether discrete or pervasive, has a first-order influence on the distribution of deformation (Morley, 1999). Many of these pre-existing features lead to linear zones of anomalously weak material within the lithosphere. Heterogeneities and anisotropies can therefore be treated as equivalents and are referred to as weak seeds that can be incorporated in analogue models (e.g., Zwaan and Schreurs, 2017; Molnar et al., 2017, 2018). Factors that play a major role on rift evolution include the location and obliquity of pre-existing structures with respect to the direction of extension (e.g., Acocella et al., 1999; Ebinger et al., 2000; Agostini et al., 2009) and temporal variations of the regional stress field (e.g., Ebinger et al., 2000; Morley, 2010). Observations from many continental rifts have greatly improved our understanding of rift evolution as a function of the obliquity of pre-existing weaknesses (e.g., Ring, 1994; Morley, 1999; Corti et al., 2003; Morley et al., 2004; Autin et al., 2010). Natural examples show a wide variety of geometric relationships, covering the entire range between rift-subparallel weaknesses (e.g., Rukwa Rift, East Africa; Morley, 2010), to approximately rift-perpendicular weaknesses (e.g., Main Ethiopian Rift, East Africa; Bonini et al., 2005, Corti et al., 2018a, 2018b) (Fig. 1). Understanding the temporal evolution of rifts relies on the preservation of timing and/or kinematic constraints within the geological record. When such evidence is lacking or when it is obscured by the superimposition of tectonic events, deciphering the temporal evolution of rifts can become a challenge.

Complementary analogue (e.g., McClay and White, 1995; Acocella et al., 1999; Corti, 2008; Agostini et al., 2009; Erbello et al., 2016) and numerical experiments (e.g., Brune, 2014; Brune et al., 2017; Naliboff et al., 2017; Mondy et al., 2018) have provided important insights on the development and interaction of extensional faults, which help to constrain the overall temporal evolution of rifts. Temporal changes in the extension direction are a common feature of rifts in nature (e.g., Morley, 2010; Nemčok et al., 2016; Whittaker et al., 2016; Pongwapee et al., 2018). Such changes have typically not been included in geodynamic models of continental extension. Previous modeling studies that focused on the role of multi-phase rifts have successfully recreated different extension directions through time (e.g., Bonini et al., 1997; Keep and McClay, 1997; Henza et al., 2010; Erbello et al., 2016; Duffy et al., 2017), but achieved this by imposing laterally homogeneous extensional boundary conditions during each phase. Conversely, the inclusion of a rotational kinematic boundary condition in experiments (e.g., Tron and Brun, 1991; Benes and Scott, 1996; Molnar et al., 2017) has effectively simulated progressive changes in extension direction over time. Such changes may induce a laterally heterogeneous strain field that results in a gradient in the amount of extension along the rift axis (e.g., Gofa Province of the Main Ethiopian Rift; Philippon et al., 2014). When these conditions are met, rift branches usually propagate laterally toward a pole of rotation. How they propagate across transversal pre-existing crustal heterogeneities remains poorly understood.

Although both pre-existing heterogeneities and rotational extension are often identified as important factors in ancient and modern rifts, their combined influence has rarely been tested by means of analogue or numerical modeling (e.g., Benes and Scott, 1996; Molnar et al., 2017). In this study, we present a series of three-dimensional laboratory experiments addressing this topic. By measuring surface deformation in the experiments at high spatio-temporal resolution, we characterize the detailed structural evolution of continental rifts as a function of the obliquity of crustal heterogeneities. We compare our results to selected natural examples in order to provide new insights on rift migration and propagation across discrete pre-existing crustal heterogeneities.

All laboratory experiments presented here are composed of a three-layer, brittle-ductile model lithosphere that floats isostatically on a low-viscosity model asthenosphere, contained within an acrylic tank (Fig. 2). All models have initial dimensions of 44 × 44 × 3 cm and are constructed to sit within two U-shaped confining walls (Fig. 2A). The experimental setup was designed to simulate continental rifts that propagate and interact with variably oriented crustal anisotropies. We therefore mimic propagating rifts by imposing a rotational kinematic boundary condition on the models following the approach used in Molnar et al. (2017). The model lithospheric plate is attached to a moving U-shaped wall and pulled by a linear actuator at a constant divergence rate (Fig. 2A). A hinge between the two confining walls acts as a fixed vertical pivot, creating a horizontal gradient in the imposed extensional strain field (Fig. 2A). The amount of extension in the models increases away from the pivot (Fig. 2C), which is analogous to creating an along-rift axis gradient in the amount of extension, as inferred from modern examples (e.g., Gofa Province; Philippon et al., 2014; Woodlark Basin, southwestern Pacific Ocean; Mondy et al., 2018).

We scale down forces, length, and time using dimensional scaling theory (Hubbert, 1937; Ramberg, 1967; Davy and Cobbold, 1991) in order to achieve an appropriately sized model that is kinematically and dynamically similar to the natural system. The scaling ratios (Table 1) for the experiments were set so that 1 cm in the models corresponds to 25 km in nature and 1 h in the experiment represents ∼0.8 m.y. in nature. Since all experiments were carried out under the normal field of gravity (1 g), the gravitational scaling ratio is 1:1. The velocity scaling ratio is derived from the time, length, gravity, and viscosity scales, such that the divergence rate of ∼5 mm/hr used in the experiments corresponds to ∼16 mm/yr in nature. Scaling factors and experimental parameters are detailed in Table 1.

The model lithosphere is composed of a mixture of quartz sand and hollow ceramic microspheres to replicate the brittle upper crust, polydimethylsiloxane (PDMS) to model the ductile lower crust, and a mixture of black Colorific Plasticine®, 3M® hollow glass microspheres and PDMS to replicate the lithospheric mantle (Figs. 2A and 2B). The upper crustal granular material mixture obeys the Mohr-Coulomb failure criterion, with an angle of internal friction ϕ< 38° and cohesion ∼9 Pa (Molnar et al., 2017). PDMS and the PDMS-based mixture used for the ductile layers are quasi-Newtonian fluids frequently used in analogue modeling studies (e.g., Weijermars, 1986; Cruden et al., 2006; Riller et al., 2012; Molnar et al., 2017). The model lithosphere overlies a model asthenosphere that contains a solution of Natrosol^® 250 HH plus sodium chloride in deionized water (Davaille, 1999; Boutelier et al., 2016). A complete description of the experimental setup, materials, scaling, and model construction can be found in Molnar et al. (2017).

We introduced discrete pre-existing crustal anisotropies (Morley, 1999) into the experiments following previous analogue modeling approaches (e.g., Corti, 2008; Le Calvez and Vendeville, 2002; Zwaan et al., 2016). Two ∼3-mm-wide and ∼1-mm-high linear semi-cylindrical seeds of PDMS were positioned on top of the homogenous lower crust prior to sifting the brittle upper crustal granular mixture. Since the brittle upper crust determines the maximum strength of the model lithosphere (Fig. 2B), the inclusion of these linear seeds locally reduces its thickness and hence the integrated resistance of the lithosphere to deformation. This technique may be used to simulate planar zones of weakness such as fault planes, shear zones or areas of previously thinned crust (Morley, 1999; Zwaan et al., 2016), although the thickness of the weak seeds (i.e., <10 km) would be more appropriate to compare with strike-slip fault zones or major ductile shear zones in crystalline basement that have been exhumed to upper crustal levels.

Rather than trying to reproduce a specific natural scenario, our aim is to characterize how rifts propagate across crustal weaknesses as a function of their obliquity (α) with respect to the initial extension-normal direction (Fig. 2D). We set a distance of 5 cm between the two discrete linear anisotropies in order to allow enough space for rift branches to develop freely between them and to study their propagation and temporal evolution. For ease of analysis and discussion, we subdivide the resulting fault pattern into proximal, central, and distal domains, which are delimited by the position of the linear anisotropies (Fig. 2E).

We employed a particle image velocimetry (PIV) system to monitor deformation and surface topography during the experiments. The PIV system consists of two stereoscopic cameras and one top camera that are set to capture digital images of the experiments at 2 min intervals, which provides quasi-continuous monitoring of deformation throughout the entire 36 h duration of the experiments. The digital images are processed using a cross-correlation algorithm (Adam et al., 2005), which allows the extraction of high-resolution surface topography, differential and cumulative strain, and velocity data from the models. See Supplemental Data¹ for details on the methods used for PIV data calculations.

A total of 12 experiments were carried out. Each experiment had identical rheological layering and kinematic boundary conditions. The angle of obliquity of the crustal anisotropies was increased by 15° steps between each experiment. We present and describe the detailed results of six experiments in which the linear anisotropies were progressively rotated clockwise between experiments. Additional experiments with anticlockwise-rotated weak seeds were carried out for comparison and showed similar patterns and timing of deformation. Figure 2E and Table 2 summarize the initial conditions of the experiments. Results are divided into three stages based on experiment running time and the corresponding amount of displacement imposed at the free side of the model. The pole of rotation is designated as the north side for reference; hence, the amount of extension increases from north to south (Fig. 2C). Surface topography and differential strain for each time-step are used to display the general evolution of propagating rifts in each experiment. Rose diagrams of fault distributions and fault activity graphs are presented for selected experiments. Fault activity graphs (panel D in Figs. 3–8) illustrate the experimental time (X axis) during which each structure had been accommodating extension, arranged by the distance to the pole of rotation (Y axis). In this way, it is possible to visualize the timing of activity/inactivity for the main extensional faults and depict the temporal evolution of the main rift faults for each of the previously defined proximal, central, and distal domains.

Strain initially localized in the distal domain as extension was accommodated by linear, normal faults that formed perpendicular to the extension direction (Fig. 3A). Increasing stretching resulted in the propagation of these major boundary faults toward the pole of rotation, forming northward pointing V-shaped basins (Figs. 3A–3B). Differential strain measurements indicate the incipient formation of new rift graben in the central domain, bounded by linear normal faults that also strike approximately perpendicular to the stretching direction (Fig. 3B). After 18–19 h (14.6–15.4 m.y.), activity on the distal domain progressively decreased, rift tip propagation slowed down and the central domain rift graben started to accommodate most of the extension (Fig. 3D). The bounding normal faults in the central domain formed V-shaped graben that propagated both north and south, with the eastern graben developing more rapidly and forming deeper basins (Figs. 3B–3C).

The overall deformation was gradually transferred toward the east and to the proximal domain, with faults that grew toward the pole of rotation due to the anticlockwise opening of the moving plate (Figs. 3A–3C). The main boundary faults in the distal domain were active for 26 h (i.e., 21.1 m.y.) (Fig. 3D), after which strain was transferred to the rift depressions (Fig. 3C). Normal faults in the central domain accommodated most extension between 18 and 27 h (i.e., 14.6 – 21.9 m.y.) (Fig. 3D). Finally, extension was predominantly accommodated in the proximal domain after 27 h (i.e., 21.9 m.y.). Comparison between fault orientations over the different time steps indicate that the boundary faults were oriented perpendicular to the extension direction during the entire experiment, and rift graben propagated both away from and toward the pole of rotation without being laterally offset as they propagated (Figs. 3A–3C).

Deformation became highly localized at early stages, forming two elongated rifts aligned parallel to the crustal anisotropies. The presence of a linear weak seed produced anomalously long and straight bounding normal faults that reach their near-final lengths early in the experiments (Fig. 4A). Differential strain measurements showed a more evolved stage in the eastern rift graben, associated with its proximity to the moving plate (Fig. 4A). An incipient N-S–oriented V-shaped central rift branch, located between the two main rift compartments, was detected by differential strain measurements before developing a visible feature on the surface (Fig. 4A). Rift branches that propagated northward beyond the crustal anisotropies showed a more irregular pattern, but orientation analysis indicates that the fault trends remained perpendicular to the extension direction (Fig. 4A).

Deformation evolution followed the behavior observed at early stages, with the widening of the two main rift branches accommodating most of the extension (Fig. 4B). Compared to the early stages, the differential strain field map shows maximum strain values closer to the pole of rotation, indicating a progressive northward migration of deformation (cf. Figs. 4A and 4B). Rifts beyond the northern tips of the linear anisotropies continued their northward propagation as short arcuate segments but with an overall N-S trend. Fault development in the proximal domain was not controlled by pre-existing anisotropies, resulting in smaller, slightly laterally offset rift basins that subsequently merged into a single graben (Figs. 4A–4C). Although deformation in the eastern main branch initially propagated farther north, higher values of differential strain were recorded in the western main branch at advanced stages (cf. Figs. 4A and 4C). Secondary rift branches with faults trending perpendicular to the extension direction developed on both sides of the linear anisotropy (Fig. 4C). During the final stages, maximum values of strain were observed in the proximal domain, near the pole of rotation, while deformation was more distributed and limited to the basin floors in the distal domain (Fig. 4C).

Deformation became strongly localized early in the experiment due to the presence of the N30°E oriented crustal anisotropies, as observed in Model B (α = 0°) (Fig. 5A). Extension was accommodated almost entirely by elongated, N30°E oriented sinistral oblique normal faults that dislocated the distal and central domains (Fig. 5A). The only structural features that developed outside the linear weak area were two short, N-S–oriented faults that delimited a rift depression in the distal domain of the model (Fig. 5A). With increasing deformation, this rift branch propagated northward and interacted with the southernmost anisotropy. This interaction resulted in a locally heterogeneous strain field with relatively high differential strain in the south of the anisotropy, while extension continued to be accommodated by the early-formed N30°E oriented elongated normal faults toward the north of the anisotropy (Fig. 5A). When deformation propagated beyond the anisotropy and into the central domain, the new rift bounding normal faults developed at high angles with respect to the stretching direction (Figs. 5B–5C).

Differential strain measurements indicate a progressive north-eastward propagation of strain along the southernmost anisotropy (Figs. 5A–5C). Deformation northeast of this anisotropy created a fault system with a horsetail geometry (Fig. 5B), in which extension was accommodated by three curved splays that propagated toward the pole of rotation (Figs. 5B–5C). The westernmost rift branches of the horsetail splay were bounded by normal faults trending approximately parallel to the central domain rift branch (Fig. 5C). Strain was mainly localized along the southernmost crustal anisotropy throughout the entire experiment (cf. Figs. 5A–5C). Consequently, only minor sinistral oblique-slip extensional kinematics developed along the northern crustal anisotropy and faults were completely absent in the proximal domain (Fig. 5C). This indicates that the low obliquity angle of the crustal anisotropy favored dislocation of the model lithosphere into kinematically independent blocks, resulting in low relative strain within the central and proximal domains (Figs. 5A–5C).

Early stages of deformation were characterized by partitioning of deformation into two N-S–oriented rift graben in the distal domain and two other V-shaped rift basins bounded by N45°E oriented linear normal faults localized along the crustal anisotropies (Fig. 6A). These narrow structures have sinistral oblique-slip components, indicating that the weaker crust favored dislocation of the model lithosphere into blocks with different displacement rates as in Model C (α = 30°) (Figs. 6A–6C).

Boundary faults of the distal domain rift branches propagated northward, gradually changing their orientation as they propagated toward the southernmost linear anisotropy (Fig. 6A). A series of three new rift branches with bounding faults trending ∼N30°W formed in the central domain after 12 h (∼9.7 m.y.) (Figs. 6B and 6D). These basins formed at slightly obtuse angles with respect to the extension direction, and at ∼80° with respect to the weakness-related faults (Figs. 6B and 6E). The orientation of the central domain graben may be associated with a local modification of the stress field caused by horizontal displacement along the two linear anisotropies. This produced a zig-zag pattern and consequently resulted in a rhombic basin geometry (Figs. 6B–6C). Rift branches that developed in the central domain were laterally offset with respect to those formed previously in the early formed distal domain, except for the westernmost graben (Fig. 6B).

Differential strain measurements indicate that strain along the oblique-slip fault zones propagated in a northeast direction (Figs. 6A–6C). Figure 6B shows that strain along the northernmost oblique-slip fault zone stopped propagating after ∼16.5 h (∼13.4 m.y.), when strain was relayed into the N30°W trending proximal domain rifts (Fig. 6B). Rift branches in the proximal domain developed with a right-lateral offset with respect to those in the central domain (Figs. 6B–6C). Orientation analysis of the proximal domain normal faults show a progressive rotation toward the north due to their proximity to the pole of rotation (Figs. 6C and 6E). Distal and central domain structures kept accommodating deformation during the final stages, but to a lesser extent than those in the proximal domain (Figs. 6C and 6E). The highest differential strain values during the late stages were detected at the intersection between the northward propagating central domain graben and the northernmost crustal anisotropy. These areas, in turn, coincide with greater fault throws and hence formed the deepest basins recorded in the experiment (Fig. 6C).

Extension at early stages was predominantly accommodated by N-S–trending faults that delimit two main rift branches that propagated northward, toward the pole of rotation (Fig. 7A). Differential strain measurements show minor localization along the crustal anisotropies, which results in short, N60°E oriented sinistral oblique normal faults (Fig. 7A). The northward propagation of the main rift branches was laterally deflected eastward by the presence of the anisotropies, as evidenced by the sinusoidal pattern of the differential strain field (Fig. 7A).

Most extension during the experiment was accommodated by N-S rift branches (Figs. 7A–7C). Maximum computed values of differential strain at early stages were detected in the western rift branch, however, the faster northward propagation of the eastern rift indicates a progressive migration of deformation toward the moving plate (cf. Figs. 7A–7C). Advanced stages show distributed deformation in the distal domain, with low strains localized in the rift depression, indicative of activity on internal rift faults or hyperextension in the ductile layers of the lithosphere (e.g., Manatschal, 2004) (Fig. 7C). Oblique-slip extensional faults associated with the crustal anisotropies decreased their activity with ongoing deformation and became completely inactive after 23 h (∼18.6 m.y.) (Fig. 7C).

Two extension-normal rift branches developed early in the experiment and propagated northward defining narrow V-shaped basins (Fig. 8A). The orientation of the bounding normal faults was not affected by the crustal anisotropies, but a pulsed northward progression of strain was detected in the computed differential strain field maps (Fig. 8A). A decrease in the amount of strain toward the north coincided with the location of the rift-perpendicular linear anisotropies (Fig. 8A). This suggests that deformation migration stalled in the proximity of the anisotropies and rift branches continued to propagate only after additional strain was accumulated (Figs. 8A–8B).

Early in the experiment the eastern rift branch displayed a more advanced deformation, having propagated beyond the southernmost linear weak area and subsequently slowing down in the proximity of the northern anisotropy (Fig. 8A). However, with ongoing deformation, the western branch propagated northward faster, and beyond the eastern branch (Figs. 8A–8C), showing a similar progression to Model B (cf. Figs. 4 and 8).

Maximum computed strain values were located closer to the northern end of the model, indicating a progressive migration of deformation toward the pole of rotation (Fig. 8C). Differential strain measurements at the final stages do not exhibit compartmentalization of deformation, indicating that the propagation-inhibiting role of the crustal anisotropies was mostly active during the incipient rupture of the crust and not at later stages of rift deepening and propagation (cf. Figs. 8A–8C). Nevertheless, digital elevation maps display an along-axis variation in rift basin depth, with deeper sections located between the crustal anisotropies (Fig. 8C).

A complementary set of experiments was carried out with an anticlockwise orientation of the crustal heterogeneities and the same obliquity angles used above (models C2, D2, E2, G2, H2; Table 2; Fig. 9). Similar first-order structural patterns and temporal evolution of deformation were observed between experiments with the same angle of obliquity but opposite anisotropy orientation, with minor differences in the location of extensional faults outside the anisotropic linear zones (Fig. 9). While the position and relative orientation of the heterogeneities with respect to the stretching direction has a significant impact, fault orientations and the influence of the crustal heterogeneities on northward propagating rifts were comparable between experiments with clockwise and anticlockwise oriented anisotropies.

Comparison between differential displacement field maps of selected clockwise (Models C, D, and E) and anticlockwise (Models C2, D2, and E2) experiments provides insights on how the anisotropy orientation affects the sense of motion on oblique-slip fault zones. Sinistral oblique-slip motion was detected early in the experiments when the linear anisotropies were oriented clockwise with respect to the rift propagation direction. Conversely, anticlockwise oriented linear anisotropies favored dextral oblique-slip reactivation. The influence of the anisotropy orientation on the oblique-slip component of motion explains the differences observed in the lateral offset of propagating rift branches between displaced domains (Fig. 9).

The modeling results presented show how the presence and orientation of crustal anisotropies affects the propagation and evolution of continental rifts. The overall deformation pattern is broadly characterized by extensional structures that initiate at the southern end of the model and propagate northward toward the pole of rotation, in agreement with the general structural evolution previously described for propagating rifts (Courtillot, 1980; Hey et al., 1980; Courtillot, 1982, Martin, 1984). In the homogeneous lithosphere experiment, subparallel structures with a horst-graben style of deformation were distributed over a 10–11-cm-wide (250–275-km-wide) area in the south end of the model. This morphology resembles the classic wide rifting mode of extension (e.g., Benes and Davy, 1996) but with less distributed deformation. As rifts propagated toward the northern end, the deformed area narrowed (up to 100–125 km) and resulted in a style of deformation similar to a narrow rifting mode.

Differences in the spatio-temporal evolution and rift architecture between experiments with similar rheological and kinematic initial conditions were solely dependent on the angle of obliquity of the crustal anisotropies, α (Fig. 9). Extensional structures in models with α ≤ 30° (Models B, G, G2, C, and C2) were strongly localized along the crustal anisotropies at early stages. This localization results from the high angle between the anisotropies and the maximum instantaneous stretching direction and persisted throughout the experiment, although Model B (α = 0°) showed a wider distribution of deformation at advanced stages (Figs. 9 and 10). The overall fault development was greatly controlled by the linear weak sections; long, linear normal faults that aligned with the anisotropies accommodated most of the extension, while only minor structures, trending approximately perpendicular to the extension direction, developed outside or beyond the linear weak areas. Displacement field evolutionary maps (Figs. 10 and 11) illustrate how the reactivation of the linear anisotropies as normal faults with a left-lateral component was easier when the obliquity angle was between 15° and 45°. This resulted in differential relative displacement between the distal and central domains of the experiments, which was maximum for Model C (α = 30°, clockwise) (Fig. 10).

Experiments with a linear anisotropy oriented between 30° and 60° (Models D, D2, E, and E2; Fig. 11) were characterized by complex interactions between the northward propagating normal faults and the crustal anisotropies. Localization along the anisotropy was strongest at early stages, when normal fault sets with sinistral oblique-slip components developed, resulting in kinematic dislocation of the model lithosphere into discrete domains (Fig. 11). In both experiments, the oblique-slip activity gradually decreased with ongoing stretching and deformation was progressively transferred to extension-normal structures. This interaction resulted in rift branches that are laterally offset as they propagated northward across the anisotropies, resulting in rhombic fault-bounded domains (Figs. 9–11).

When α > 60° (Models H, H2, and F), the acute angle between the anisotropies and the stretching direction inhibited strain localization and extension was almost entirely accommodated by N-S–trending normal faults. The crustal anisotropies had minimal effects on fault trends during their northward propagation and induced minor lateral offsets (Figs. 9 and 10). A pulsed propagation of the extensional basins was observed only when the crustal heterogeneities were oriented parallel to the initial extension direction (Model F, Fig. 8). Deformation analyses show that strain accumulated in steps as it propagated along the main rift branches, and only penetrated through the pre-existing heterogeneities after a threshold strain was reached. This mechanism of propagation consequently led to deeper rifts in the vicinity of the crustal anisotropies.

The experimental results can be compared to selected natural examples of continental rifts that have been affected by heterogeneous extensional stress fields and/or that experienced changes in extension direction over time. Due to the intrinsic complexity of natural rift systems, we cannot replicate all of the variables that control the evolution of continental rifts. In addition, the dimension of the deformed regions in the selected examples is lower than that observed in the experiments. However, we argue that this does not change the extrapolation of the main findings of modeling to nature and that our experiments provide first-order insights on the spatio-temporal evolution of propagating rifts.

Here, we discuss comparisons between our experimental results and natural examples as a function of the angle of obliquity of the crustal anisotropies with respect to the extension-normal direction, α. For ease of visualization, some modeling results are mirrored and/or rotated in Figure 12 to compare with the natural prototypes.

The Tanganyika and Rukwa rift basins are located in the central segment of the Western Branch of the East African Rift System (EARS) (Fig. 12A). Present-day kinematics and recent fault activity is well constrained by geodetic, geophysical, and focal mechanism data (e.g., Morley et al., 1992; Delvaux and Barth, 2010; Saria et al., 2014), but the timing of events that led to the formation of these NW-SE–oriented narrow basins remains controversial (Morley, 2010 and references therein). It is widely accepted that a series of NW-SE rift basins of Karoo age (Permo–Triassic) (Morley et al., 1992; Morley, 1999; Delvaux, 2001) played an important role in controlling the Oligocene to present-day (Roberts et al., 2012) rifting of the Tanganyika-Rukwa area. The general direction of extension affecting the EARS is a contentious issue (e.g., Ebinger, 1989; Ring et al., 1992; Morley, 1999; Chorowicz, 2005; Calais et al., 2006; Bellahsen et al., 2013) (Fig. 12A), but recent studies support a general E-W regional present day extension (Corti et al., 2007; Delvaux and Barth, 2010; Heidbach et al., 2010; Morley, 2010). Additionally, two-dimensional numerical simulations indicate that larger rifts may capture rigid cratons, which can lead to local rotations (Koehn et al., 2008). These rotations accommodate the orthogonal motion of larger plates and are consistent with the present-day anticlockwise rotation of the Victoria plate—containing the Tanzania craton—with respect to the Nubian plate (Calais et al., 2006; Koehn et al., 2008; Stamps et al., 2008) (Fig. 12A). The Tanganyika and Rukwa rifts are good analogues for the experiments with low obliquity due to the characteristics of the boundary conditions; the Karoo-age pre-existing discontinuities are oriented almost perpendicular to the local direction of extension and there is a comparable rotational component in the extensional system Fig. 12 and Supplemental Figs. S2–S4 [^{footnote 1}]).

The Tanganyika-Rukwa rifts are narrow, deep basins bounded by a series of anomalously long, straight NW-SE–trending faults that follow pre-existing heterogeneities (Morley et al., 1992; Ring, 1994; Corti et al., 2007) (Fig. 12B). Pure dip-slip motion characterized these long, straight boundary faults at early stages (Foster and Jackson, 1998; Delvaux, 2001; Morley, 2010), indicating a local stretching sub-perpendicular to the pre-existing heterogeneities (Fig. 12B). Since the regional extension direction was roughly E-W, this local feature could be explained by a deflection of stresses by inherently weak basement fabrics (Morley, 2010). A similar behavior was observed in our low obliquity experiments, which developed long, linear, dip-slip normal faults (Fig. 9 and Supplemental Figs. S2-S4). The overall rift pattern in these models was reached at early stages and remained almost constant throughout the whole experiment (Figs. 4 and 5). This is also consistent with field observations in the Tunganyika-Rukwa rifts, which suggests that deformation pattern was relatively static since its initial formation (Morley, 2002).

The deflection of the Tanganyika rift from a NW-SE direction to a more N-S orientation can be explained by the closer proximity of the northward propagating rift to the relative pole of rotation of the Victoria plate (Stamps et al., 2008) (Fig. 12A). All of our low obliquity experiments showed strong initial localization along pre-existing heterogeneities, followed by later deflection of propagating rifts toward the pole of rotation (Figs. 4, 5, and 9). Models B, C2, G, and G2 have partial resemblance to this segment of the EAR. Only Lake Rukwa and the southern section of Lake Tanganyika were strongly influenced by pre-existing heterogeneities (Morley, 1999, 2010), which favored the formation of long linear faults near the anisotropy but led to more sigmoidal or arcuate normal faults in the central and northern Lake Tanganyika areas (Fig. 12B). Similar observations were replicated in Models B and G2 (Fig. 12C). The development of these sigmoidal faults probably occurred synchronously with the northward propagation of the rift (Macgregor et al., 2017). Smaller scale discrete fabrics in the upper crust also contribute to this pattern, as shown by previous analogue models of fabric reactivation (Corti et al., 2007). Additionally, an area in which a propagating rift branch ignored a pre-existing heterogeneity was observed in Model C2 (Fig. 12D), which could be analogous to the region west of southern Lake Tanganyika, in which no significant reactivation occurred.

Sediment thicknesses in the Rukwa rift are less than in the more extension-orthogonal central Tanganyika rift (Chorowicz, 2005; Corti et al., 2007; Misra and Mukherjee, 2015) (Fig. 12B). This is similar to the evolution of models G and G2, in which the extensional basins that formed perpendicular to the stretching direction were deeper (Fig. 12E and Supplemental Figs. S2–S3 [^{footnote 1}]). The two heterogeneity-related rift branches accommodated most of the extension early in the experiments, but soon after deformation was transferred to a single rift branch with larger fault throw and therefore more accommodation space for sedimentation (Fig. 12C). The experimental evolution of deformation also indicates strong uplift of the western rift shoulders, which compares to the higher elevations in the western areas of Lake Tanganyika and Lake Rukwa (cf. Figs. 12B and 12E), which were probably aided by magma underplating (Morley et al., 1992).

A series of Pleistocene NE-SW–trending extensional structures (e.g., Usangu rift and Lake Mweru) formed perpendicular to the main boundary faults of the Tanganyika-Rukwa rifts (Ebinger, 1989; Morley et al., 1992) (Fig. 12B). Field observations also suggest a relative change in the extension direction that reactivated earlier dip-slip faults as strike- and oblique-slip faults (Ring, 1994; Morley, 1999). The absence of normal faults with intermediate trends between the Tanganyika-Rukwa rifts and the Pleistocene NE-SW–oriented extensional structures indicate a rapid change in the extension direction. Although this change is supported by field observations, it may be only a local feature since some plate kinematic models (e.g., DeMets and Mekouriev, 2016) cannot explain a variation in the regional extension direction. Models B, C, C2, G, and G2 show a history of rift propagation that is partially analogous to the Oligocene–Pleistocene stages of the Tanganyika-Rukwa rifts. However, they do not explain the perpendicular structures that developed later, which implies a significant change in the extension direction that did not occur in our experiments. Such an abrupt change also explains dextral strike-slip reactivation of the main NW-SE–trending faults (Ring, 1994). Since strain tends to keep accumulating in the same locations (Morley, 2002), the pre-Pleistocene normal faults acted as pre-existing heterogeneities that facilitated their reactivation under the new stretching regime. At the same time, our modeling results show that the left-lateral motion between domains was strongest when the weaknesses where favorably oriented (i.e., 30°–45° clockwise, Fig. 11), and is consistent with a ∼60° to 90° change in the extension direction (Versfelt and Rosendahl, 1989; Delvaux et al., 1992; Ring et al., 1992) (Fig. 12B).

The Recôncavo-Tucano-Jatobá rift basins of northeastern Brazil formed during the Neocomian during the opening of the South Atlantic Ocean (Milani and Davison, 1988; Heine et al., 2013). Previous studies (e.g., Milani and Davison, 1988; Darros de Matos, 1992, 2000) have shown that complex structural framework of these rift basins was strongly controlled by pre-existing crustal heterogeneities associated with E-W to NE-SW–trending Pan-African shear zones (Milani, 1987; Magnavita, 1992; Fig. 13A). Several interpretations have proposed an anticlockwise rotation of the East Brazilian microplate (e.g., Milani and Davison, 1988; Szatmari and Milani, 1999) as a possible kinematic model to explain the spatio-temporal evolution of these basins, which opened from south to north (Fig. 13A). These conditions are comparable to experiments with intermediate anisotropy obliquity (Models D [α = 45°] and E [α = 60°]) (Fig. 13A). Fault trend analysis in the southernmost Recôncavo rift basin shows a dominant N25°E oriented normal fault network and a secondary N30°W trending strike-slip fault group (Milani and Davison, 1988). The final architecture of Model E2 shows a very similar fault distribution both in terms of both geometry and sense of motion (Fig. 13B).

Differential strain measurements of experiments with α = 60° show a top-view sigmoidal strain pattern as rifts propagated across the crustal anisotropies (Fig. 7 and Supplemental Fig. S5 [^{footnote 1}]), similar to the inferred deflecting role of pre-existing basement weaknesses in the Recôncavo rift (Milani and Davison, 1988). Such crustal anisotropies may have localized major transfer zones in the Recôncavo-Tucano-Jatobá system (Fig. 13B). Two of these shear zones dislocated the eastern section of the Recôncavo rift into three distinct compartments (Milani and Davison, 1988), similar to the displacement maps for our intermediate obliquity models (Fig. 11). The overall movement sense of reactivated shear zones and its change over time is difficult to determine in nature without corresponding kinematic indicators or crosscutting relationships. Analyzing the evolution of laboratory experiments therefore provides insights into the complex rifting history of natural examples such as northeastern Brazil.

Field observations and differences in the throw of boundary faults support a dextral strike-slip motion on the NW section of the Mata-Catu transfer fault (Milani and Davison, 1988) (Fig. 13B). However, there are few kinematic constrains for its SE section. A sinistral sense of motion has been interpreted, based on the present-day structural pattern (Milani and Davison, 1988). Our modeling provides additional constraints, since displacement field maps for model E2 (α = 60°, anticlockwise) show a similar structural pattern to the southeastern end of the Mata-Catu fault (Fig. 13B). The model results indicate constant dextral oblique- and strike-slip motion on transfer zones in the northernmost sections of the experiment (Fig. 13B), as inferred for the Itapicuru, Itaporanga, and north of the Mata-Catu transfer zones (Milani and Davison, 1988). However, the distal domain shows more distributed deformation and a more complex evolution. At later stages in the experiment there is a local increase in displacement rates southwest of the southernmost crustal heterogeneity (Fig. 13B). A possible explanation for this is progressive rotation of the moving plate. While initially unfavorably oriented, the pre-existing weakness was eventually reactivated with a dextral sense of motion (Fig. 13B). At later stages, a shift in the relative orientation between the stretching direction and the heterogeneity orientation may have caused a reactivation with a sinistral sense of motion (Fig. 13B), which correlates well with observations in the eastern end of the Recôncavo rift (Milani and Davison, 1988).

The main boundary faults north of the Jatobá Graben terminate abruptly against the Ibimirim fault, which is associated with the major ENE-WSW–trending Pernambuco Lineament (Fig. 13C). Similar behavior is observed in Model D (α = 45°). Rift propagation in the central domain halted when the northernmost crustal anisotropy was reached and strain was relayed onto the transfer zone (Fig. 6), as probably occurred in the Jatobá Graben-Ibimirim fault interaction zone (Milani and Davison, 1988). The North Tucano Graben has also been interpreted to cut through a basement weakness due to the latter’s unfavorable orientation for reactivation (Darros de Matos, 1992). This is analogous to models F (α = 90°) and H (α = 75°, anticlockwise), in which the crustal anisotropies did not affect the northward trend of rift propagation (Fig. 9).

The inherent complexity of the northeastern Brazil rift system makes it difficult to test in a single analogue experiment. However, the intermediate to high obliquity between the general extension direction and the pre-existing basement structures makes it a good reference to compare with several experimental results (i.e., Models D, E, F, and H). Detailed analysis of the experiments presented here provide plausible explanations for why planar basement weaknesses were either ignored or followed by the evolving rifts in the Recôncavo-Tucano-Jatobá system.

The NE-SW–trending Main Ethiopian Rift (MER) is located between the Afar depression and the Kenya rift in the Eastern Branch of the East African Rift System (Fig. 14A). The intricate network of extensional faults and volcanic history of the MER reveal a complex tectono-magmatic evolution (e.g., Ebinger et al., 2000; Bonini et al., 2005; Corti, 2009). The role of two major lineaments in the evolution and propagation of the MER is of particular interest to this study.

The approximately E-W–oriented Goba-Bonga (GB) and Yerer-Tullu Wellel (YTW) lineaments (Fig. 14A) are parallel to the general stretching direction (Abebe et al., 1998; Corti et al., 2018a, 2018b). The Somalian plate has been rotating clockwise with respect to the Nubian plate since the middle Miocene (Bonini et al., 2005; Stamps et al., 2008) (Fig. 14A). Recent work (Balestrieri et al., 2016; Boone et al., 2019) has illustrated a complex scenario of rift initiation in the MER; deformation propagated northward in the southern MER, while the northern MER propagated southward, into the central MER. The boundary conditions in model F (α = 90°) are therefore analogous to the tectonic setting of the northern MER (Fig. 14A).

Extensional structures developed in Southern Afar during the middle Miocene formation of the Afar triple junction but did not propagate beyond the YTW lineament at that time (Bonini et al., 2005) (Fig. 14B). This is consistent with the evolution of Model F, in which propagation of initial rift-parallel structures appeared to slow down when they approached the first crustal anisotropy (Fig. 14C). Extensional deformation continued to propagate southwards beyond the YTW lineament and into the central MER in the late Miocene as a result of clockwise Somalian plate rotation (Collet et al., 2000; Bonini et al., 2005) (Fig. 14B). At this time, rift propagation was toward the SW, perpendicular to the general direction of extension, as observed in the Model F (Fig. 14C). It has been suggested that the presence of the YTW led to a minor right lateral offset of the rift segments (Keranen and Klemperer, 2008; Corti, 2009), which was also observed in Model F (Fig. 14C). Southwestward rift propagation was again hindered during the Pliocene by the transversal GB lineament, as interpreted by Bonini et al. (2005) and similarly observed in Model F (Fig. 14B).

Earlier interpretations (Bonini et al., 2005) suggested that during this period the southern MER (SMER) underwent a tectonic pause, which is consistent with the evolution of our experiment. However, recent thermochronological data from the SMER (Balestriei et al., 2016) indicate that rifting in the SMER activated prior to extensional deformation in the central MER, which indicates an even more complex scenario. In spite of the differences observed between nature and experiments due to the high complexity of the Ethiopian Rift, the overall evolution of the northern MER shows a good match with Model F, in which progressive rotational stretching resulted in the propagation of rift branches across transversal crustal anisotropies (see Fig. 8). While the anisotropies appear to have little influence in the orientation of the extensional structures (Fig. 12), they did play an important role stalling rift propagation, resulting in a pulsed rift history (e.g., Bosworth, 2015). In the experiment, the final stages were characterized by the southward continuation of rift propagation across the linear crustal anisotropies (Fig. 14B). This scenario resembles earlier interpretations of the final stages of the MER evolution during the Pleistocene, characterized by the propagation of the rift across the GB lineament and subsequent reactivation of early Miocene structures in the SMER (Bonini et al., 2005) (Fig. 14B). Although the previously mentioned data from the SMER (Balestrieri et al., 2016) represents evidence against this order of segmented rift activation events, the significant role of crustal anisotropies as obstacles to rift propagation is still comparable. Northward propagation of the SMER, associated with deformation in the Kenya Rift (Bonini et al., 2005; Balestrieri et al., 2016), was temporarily halted when it reached the GB lineament during the early Miocene. We therefore consider crustal anisotropies as structures that can cause a partial impediment to rift propagation. The development of locked zones has been identified in oceanic rifts (e.g., Courtillot, 1982; van Wijk and Blackman, 2005) but their origin in continental rifts remains poorly understood.

Kinematic rotation in three-dimensional laboratory experiments simulates realistic initial conditions for natural extensional settings with laterally heterogeneous stress fields. The use of a stereoscopic PIV system to monitor experiments provides key insights on the temporal evolution of rift development and propagation across crustal heterogeneities. The results illustrate how the orientation and obliquity of rheological heterogeneities in the crust have a strong influence on fault trends and lengths, the sense of movement on reactivated oblique-slip zones, and on the role crustal anisotropies play in dislocating the crust into kinematically independent blocks.

The experiments show that extension-normal and low obliquity heterogeneities (α ≤ 30°) favor the formation of anomalously long normal faults, which reach their near-final length early in the evolution of rifts in comparison to cases with a homogeneous lithosphere. Low to intermediate obliquity angles (∼15° ≤ α ≤ 45°) oriented clockwise with respect to the extension direction favor sinistral oblique-slip reactivation of linear anisotropies and subsequent dislocation of the crust into blocks with different displacement rates. Such dislocation is most effective and long lived when α= 30°. Conversely, comparable obliquity angles with an anticlockwise orientation results in conditions that favor reactivation of heterogeneities as dip-slip or dextral oblique-slip normal faults. Nevertheless, the amount of oblique-slip or dip-slip deformation in all experiments varies over time and specific interpretations must be made for each particular scenario. Overall, our results demonstrate that the lateral sense of motion on faults is determined by the relative orientation of crustal anisotropies with respect to the stretching direction. Finally, very high obliquity crustal anisotropies (α > 60°) have little influence on the trends of border faults in propagating rifts. However, they may cause pulses in along-axis rift propagation, which can subsequently result in variations in subsidence rates, fault throws and depocenter depths along the rift.

N.E.M. was supported by a Monash Graduate Scholarship. A.R.C. was supported by the Australian Research Council Discovery grant DP110103387 and Monash University. The manuscript benefited from the valuable suggestions and constructive comments of reviewers Giacomo Corti and Chris Morley.