Modeling mobile shales under contraction: Critical analyses of new analog simulations of shale tectonics and comparison with salt-bearing systems

Where present in a basin, mobile substrates typically exert a primary control on deformation styles in all structural settings (e.g., Jackson et al., 1995; Wood, 2010; Jackson and Hudec, 2017; Rowan, 2020; Hudec et al., 2023). The presence of mobile layers profoundly influences fold and fault geometries, and diapiric bodies of mobile material localize deformation due to their weak nature, heavily influencing the structural grain especially under shortening (e.g., Nilsen et al., 1995; Rowan et al., 2004; Wood et al., 2004; Rowan and Vendeville, 2006; Dooley et al., 2009, 2015; Callot et al., 2012; Duffy et al., 2018; Restrepo-Pace, 2020; Rowan, 2020; Santolaria et al., 2021; Hudec et al., 2023; Spina and Mazzoli, 2023; Adineh et al., 2024). When discussing mobile substrates, the two most common are rock salt and mobile shales based on their relative ubiquity on earth (see maps of salt distributions and settings in Figures 1–5 of Hudec and Jackson [2007], and of shale distributions and settings in Figure 1a taken from Soto et al., 2021a). Although both of these materials are weak and ductile, for a long time they were viewed as similar in terms of deformation dynamics and kinematics, as mobile shales form many observable features on the surface and subsurface that are superficially similar to those seen in salt tectonics (Figure 1b–1e; Soto et al. [2021a, 2021b] and Hudec et al. [2023] for more discourse on this).

However, the mechanical behavior of salt and mobile shales are very different as evidenced by discussions in a variety of relatively recent publications, many with contrasting views on the mechanical behavior of “mobile” shales (e.g., Wood et al., 2004; Day-Stirrat et al., 2010; Maloney et al., 2010; Wood, 2010; Morley et al., 2017, 2018; Restrepo-Pace, 2020; Soto et al., 2021a, 2021b; Spina and Mazzoli, 2023; Hudec et al., 2023). Rock salt is essentially a crystalline solid that flows, by creep, in a ductile fashion under most geologic conditions because of its mineralogy (Figure 2; Urai et al., 1986; Jackson and Hudec, 2017). The viscosity of rock salt varies as a function of temperature, crystal size, and impurities among other factors (see Urai et al. [1986] and van Keken et al. [1993] for further details) but, fundamentally, rock salt is always weak compared with the country rocks that surround and lie above it. In contrast, mobile shale originates as a brittle solid and only becomes mobile when it reaches critical state (Figure 2; Brown, 1990; Soto et al., 2021a, 2021b; Hudec et al., 2023). Soto et al. (2021a) define mobile shales as “bodies of clay-rich sediment or sedimentary rock undergoing penetrative, (visco-) plastic deformation at the critical state”, and once it departs this critical state it returns to being a brittle solid. What determines this critical state, when the imposed shear stress exceeds the shear strength of the shales (Soto et al., 2021a), will depend on a variety of factors within the shale body, such as the magnitude of overpressure, diagenetic changes, hydrocarbon generation, degree of consolidation, and burial depth, among others (Morley and Guerin, 1996; Morley et al., 2017, 2018; Soto et al., 2021b).

As noted previously, weak mobile substrates, such as salt and mobile shale, have a profound influence on shortening styles in basins of all types globally (Davis and Engelder, 1985; Wood et al., 2004; Duerto and McClay, 2011; Morley et al., 2011; Duffy et al., 2018; Hudec et al., 2023; Dooley and Hudec, 2024). Not only do these mobile substrates form the main detachment horizons during shortening, facilitating significant vertical decoupling but they may also flow to fill fold cores, form diapirs, and create shallower level detachments by surface extrusion and reburial, as well as forming complex vertical and horizontal linkage zones between different detachments at different stratigraphic levels (Davis and Engelder, 1985; Harrison, 1995; Rowan et al., 2004; Morley et al., 2011, 2017, 2018; Harrison and Jackson, 2014; Restrepo-Pace, 2020; Rowan, 2020; Soto et al., 2021a, 2021b; Hudec et al., 2023; Jadoon et al., 2024).

Physical, analog, or sandbox modeling has long been used to investigate the geometries and kinematics of fold-thrust belts (see reviews and historical perspectives in, for example, Koyi [1997], Graveleau et al. [2012], and Butler et al. [2020] for more details). For a long time a ductile salt analog has been incorporated into contractional models run under normal gravity by the use of viscous silicone polymers (polydimethylsiloxanes [PDMS]; Weijermars, 1986), either as a pure detachment layer or a combination of detachments and allochthonous salt bodies (Cobbold et al., 1989; Cotton and Koyi, 2000; Bonini, 2003; Couzens-Schultz et al., 2003; Rowan and Vendeville, 2006; Dooley et al., 2007, 2009, 2015; Pichot and Nalpas, 2009; Konstantinovskaya and Malavieille, 2011; Buiter, 2012; Callot et al., 2012; Graveleau et al., 2012; Driehaus et al., 2014; Santolaria et al., 2015, 2021; Borderie et al., 2019; Darnault et al., 2016; Wang et al., 2022; Alania et al., 2023; Dooley and Hudec, 2024; among many others). However, to date, a similar analog material for mobile shales has not been universally used nor investigated. Cobbold and Castro (1999), Cobbold et al. (2001), and Mourgues and Cobbold (2003) use a complex air pressure system to simulate failure above a basal overpressured detachment but this required unrealistically steep dips to initiate failure, combined with a highly complex setup that was difficult to scale and maintain. Weaker granular materials, such as glass or ceramic microbeads have been used as a lower-strength basal decollement (Rossi and Storti, 2003; see the reviews of Buiter [2012] and Graveleau et al. [2012], for more details), as well as in more complex layered stratigraphies to successfully promote the development of multiple detachments during crustal-scale shortening experiments (Konstantinovskaya and Malavieille, 2011). Other workers have simply used traditional PDMS as a pseudo proxy for mobile shales in the laboratory in their studies on delta systems (McClay et al, 1998, 2000, 2003; Wu et al., 2015), mixed-salt shale systems in the Gulf of Mexico (Hudec et al., 2019), or more complex mixtures of sands and polymers as a proxy for impure rock salts or for organic-rich shales in fold-thrust belts (Callot et al., 2012; Pla et al., 2019), or as impure salts in gravity-driven systems (Cartwright et al., 2012; Ferrer et al., 2017). However, we know from the preceding discussions that salt and shale are mechanically very different. The PDMS used as a mobile-shale proxy in many of these studies is always ductile and weak, just as rock salt is, and behaves as a near-Newtonian fluid when pure (Weijermars, 1986; Weijermars et al., 1993). In contrast, shales possess a yield strength that must be exceeded before going mobile (Soto et al., 2021a), and thus any strict laboratory analog for mobile shales must also exhibit this behavior.

Yield stress fluids have a significant elastic component in addition to the viscous component of their rheology (Reber et al., 2020). One group of nonlinear yield stress fluids called Carbopol has become popular in some segments of the analog modeling community in the past 10–15 years (Schrank et al., 2008; Balmforth and Rust, 2009; Reber et al., 2015; Di Federico et al., 2017; Birren and Reber, 2019). The shear thinning properties of Carbopol together with the ability to independently adjust yield strength and bulk viscosity make it a useful analog material to model a wide variety of applications (Reber et al., 2020). At load levels below the yield strength, Carbopol behaves as a perfect elastic body and at stresses above the yield strength Carbopol deforms in a viscoplastic manner (Reber et al., 2020), and thus is an attractive option for a mobile shale analog. As a consequence of the mechanical behavior of Carbopol, a significant amount of work has been carried out on the physical properties of this fluid group (Di Giuseppe et al., 2015), and this material has been used to model a variety of deformation processes such as shear zone localizations (Schrank et al., 2008), gravity-driven flows (Di Federico et al., 2017), and semibrittle processes (Reber et al., 2015; Birren and Reber, 2019; among others). However, as far as we are aware, this material has never been incorporated into a traditional sandbox-style model investigating regional tectonic processes.

Our goal in this paper is to present results from a preliminary series of physical models (Table 1) that incorporate Carbopol as discrete layers within an otherwise layer-cake stratigraphy of traditional brittle materials to answer the following questions:

1. How does Carbopol embedded within a multilayered stratigraphy deform once yield strength is exceeded, and can we visualize this?

2. How do the resulting structures differ from experiments using our traditional salt analog?

3. What types of shortening structures are seen below and above our weak layers?

4. How efficient is Carbopol as a detachment and how does the presence of multiple detachments serve to decouple deformation during thin-skinned shortening?

5. How do the deformation patterns in our shale analog compare with contractional deformation of a natural multilayer?

Before presenting our model results, we will first address the materials used with more detailed information on the specific makeup of Carbopol and our other model materials, the design and scaling of the models (Table 1), and the monitoring, and image processing techniques used in this study.

As mentioned previously, we use a yield stress fluid called Carbopol, specifically Carbopol Ultrez 21, as our shale analog. This material has been studied quite extensively in the past decade or so (see Di Giuseppe et al. [2015] for a deep dive into the detailed properties of this material and scaling parameters; and Reber et al. [2020] for a summary of this kind of material and its applications). Carbopol is an acrylate cross-polymer gel commonly used as a rheologic modifier in personal care products such as hair gels, and thus easy to source and safe to use (Birren and Reber, 2019; Reber et al., 2020). Carbopol is a suitable “hybrid” fluid for use as a mobile shale analog as it is made up of micro gel grains and a viscous interstitial fluid phase (Figure 3a; Lubrizol Advanced Materials Inc., 2012; Di Giuseppe et al., 2015; Reber et al., 2015; Birren and Reber, 2019).

We follow the mixing procedures outlined in Birren and Reber (2019) and modified from Di Giuseppe et al. (2015), where we mix the powdered Carbopol with deionized water slowly for a period of 48 h before its natural acidity is neutralized (neutral pH) by the addition of a base (sodium hydroxide in this case). Adding the neutralizer leads to an immediate cross linking of the polymer chains and stiffening of the mixture (Figure 3a and 3b). The yield strength of the mixture is dependent on the ratio between micro grains and the interstitial fluid, or simply the weight percentage of Carbopol powder added to the distilled water. The viscosity of the Carbopol is dependent on the pH of the mixture. While the final mixture may look like a simple fluid on a macro scale, it is a mixture of elasto-plastic grains on a micro scale (Figure 3; Reber et al., 2020). Like a granular system, the microgel grains in Carbopol can form force chains that can jam, resulting in the yield strength of the fluid. With increased deformation, the force chains can fail, and individual microgel grains can break leading to the formation of anastomosing shear zones form of broken grains and linking viscous interstitial fluids (Reber et al., 2015). The localization of deformation and formation of anastomosing deformation patterns is a result of the shear-thinning material property of Carbopol. In our mixtures, minute quantities of food colorings are added to the Carbopol prior to the neutralizing process (green in this case, Figure 3b). This was done to distinguish different Carbopol layers from each other in models with multiple layers (see Table 1).

All Carbopol-bearing experiments described in this study used a mixture of 0.25 wt%, prepared in the manner briefly described previously (Table 1; see Birren and Reber [2019], for more details). Birren and Reber (2019) document yield stress values of 25, 144, and 357 Pa for Carbopol mixtures of 1, 2, and 3 wt%, respectively. For our 0.25 wt% mixtures we measured an expectedly low yield stress value of approximately 10 Pa at the strain rate used in our experiments (Table 1). However, as pointed out by Wood (2010), Morley et al. (2017, 2018), and Soto et al. (2021a, 2021b), natural mobile shale sequences can be very weak, and there is strong evidence that they can be much weaker than salt (see Hudec et al. [2023], for more details), justifying our choice of a weak Carbopol mixture. The advantage of Carbopol with such a low yield strength, and a concomitant low viscosity (see below), is that it is relatively pourable and therefore simplifies the introduction of our shale analog into our model stratigraphies. The bulk viscosity of Carbopol mixtures can be modified by tuning the pH of the mixture independently of the weight percentage (Birren and Reber, 2019). Here we keep the viscosity constant, to reduce the number of variables, and have measured dynamic shear viscosity of approximately 90 Pa s, at a pH of 7. Density measurements by Di Giuseppe et al. (2015) indicate a value for our 0.25 wt% mixtures of between 997 and 1016 kg m⁻³ (Table 1). Finally, Di Giuseppe et al. (2015) note that Carbopol sensitivity to temperature is relatively small for typical temperature variations within a laboratory — less than 1% decrease in viscosity for each 1°C rise in temperature. In most laboratories with temperature control this would not be a problem but our models are hardened on completion by adding a mixture of boiling water and gelatin, thus leading to issues we will discuss below.

As is common for most physical models of salt tectonics, we used a viscous silicone polymer PDMS to simulate rock salt (Weijermars, 1986; Cobbold et al., 1989; Vendeville and Jackson, 1992; Rowan and Vendeville, 2006; Santolaria et al., 2021; Dooley et al., 2024). On geologic time scales, salt flows as a fluid, and silicone polymers can provide a suitable analog for rock salt in laboratory-based experiments (Weijermars et al., 1993). Silicone polymers cannot replicate all mechanical aspects of salt movement but in combination with brittle model materials, such as sand and other granular materials, many processes related to salt tectonics can be investigated, especially deformation of sediments adjacent to, and above, the salt body (see reviews in Jackson and Hudec [2017], and Reber et al. [2020], for further details). The PDMS used here has a density of 970 kg m⁻³, a dynamic shear viscosity of 5 × 10⁴ Pa s, lacks a yield strength, and has near-Newtonian viscous characteristics at typical strain rates used in laboratory experiments (Weijermars et al.,1993).

The layered brittle overburden comprised different colored mixtures of silica sands (bulk density of approximately 1700 kg m⁻³; grain size of 300–600 μm; internal friction coefficient, μ, of 0.55–0.65; McClay, 1990; Krantz, 1991; Schellart, 2000), and hollow ceramic microspheres (“glass beads” or cenospheres) having a bulk density of 650 kg m⁻³, average grain size 90–150 μm, and typical μ = 0.45 (Rossi and Storti, 2003; Dooley et al., 2009). These granular materials were mixed to varying ratios to control the density of the overburden strata (Dooley et al., 2009). The densities of our Carbopol mixtures and PDMS are very similar, and when layering our overburdens atop our salt and shale analogs a mixture of sands and cenospheres with a similar density is required to avoid triggering unintentional gravity-driven deformation of the mobile layers during the model building stage (Dooley et al., 2009). In addition, the speckled nature of these granular mixtures makes them ideal for strain monitoring using the techniques outlined subsequently (Adam et al., 2005; Reber et al., 2020).

Four models are presented in this paper and their parameters and settings are shown in Table 1 and Figure 4. All were run in a 100 cm long and 20 cm wide deformation rig with shortening driven by a moving endwall that was linked to a variably geared motor (Figure 4). In models 1–3 a 1.5 cm-thick basal sequence was laid down in each model before the mobile layer was emplaced (Figure 4). In models 1 and 3 this was accomplished by pouring of the Carbopol mixture onto the surface between two spacer blocks and because of its yield strength, leveling using a long metal blade. In model 2 the weight of PDMS required to form a 1.2 cm thick layer was calculated, and this was placed on the model surface between the spacer blocks and allowed to settle over several days. In model 4 the same procedure as models 1–2 was followed but a second layer of Carbopol was emplaced within the brittle sequence at 4 cm height and laterally offset from the hinterland edge of the lower Carbopol layer by 20 cm (Figure 4).

The reader will note that the model stratigraphy in all these models is highly simplified and layer-cake, with tabular layers of our weak materials, PDMS, and Carbopol. As they currently stand, these models simply represent a portion of a regional fold-thrust belt, rather than a full thrust belt with all the inherent stratigraphic, structural, and decollement irregularities found in natural systems as have been investigated by experimental studies such as Konstantinovskaya and Malavieille (2011), Santolaria et al. (2015), Pla et al. (2019), and Wang et al. (2022), among many others. A simple setup was required to investigate the basic behavior of Carbopol embedded within a simple granular multilayer in a 2D deformation rig. In addition, Carbopol is quite difficult to handle, and despite being pourable at the concentrations used in these experiments, it is a yield stress fluid and requires careful leveling to embed it in these sequences.

Models are dynamically scaled such that 1 cm in the model approximates to 1 km in nature, and thus our models represent upper crustal thicknesses of between 5 and 7 km, and strike lengths of 50–70 km (Table 1 and Figure 4; see, for example, McClay [1990] for detailed discussions on scaling). Models were conducted with horizontal velocities of 1.7 × 10⁻³ cm/s, or 1 cm shortening in 10 min, similar to strain rates used in other modeling studies using Carbopol and PDMS (e.g., Reber et al., 2015; Roma et al., 2017). At this strain rate models are completed in 2–3 h of runtime. Upon completion of the experiment the model is buried under postkinematic sediments with a horizontal upper surface and wetted with a boiling gelatin solution. This is left overnight and cut into slabs on a bespoke slicing machine. As mentioned earlier, Di Giuseppe et al. (2015) note that Carbopol viscosity is only mildly sensitive to temperature in normal laboratory conditions — 1% decrease in viscosity with each 1°C rise in temperature. For normal operating conditions in a laboratory this will not present any noticeable issues. However, the addition of a boiling liquid within the pore spaces of the model has a much stronger effect and we have seen almost complete loss of our shale analog in models with very thin Carbopol layers (see below) as the fluid portions are leached out at these high temperatures. We have taken some mitigation measures for these models which include lowering the temperature at which the gelatin mixture is introduced into completed models but more work needs to be done to get around this limitation.

Models are recorded by computer-controlled high-resolution DSLR cameras through the glass sidewalls of the deformation rig. Photographs of these models were taken at 30 s intervals. Exported imagery was processed through digital image correlation (DIC) software (DaVis 10 by LaVision), and placing control points within the granular and ductile layers allowing us to track motions and strain in the models in the shale analog and our brittle sediments. The reader is referred to Adam et al. (2005) for more details on DIC monitoring techniques and processing. Note that there are always some artifacts at the edges of sandbox experiments due to friction along the sidewalls of the deformation box. These effects are somewhat mitigated by the application of cleaning solutions or water repellants on the glass (e.g., Rain-X™).

The results of model 1 (see setup in Figure 4a and Table 1) are illustrated by a series of cumulative maximum shear strain images (Figure 5), incremental maximum shear strain images (Figure 6), and the final side view combined with a section through the center of the completed model (Figure 7). After 6.5 cm shortening the model consists of a growing imbricate thrust stack below the shale analog and a forward-breaking sequence of forethrusts in the supra shale strata (Figure 5a). A minor amount of shear strain is seen emanating out from this hinterland system (Figure 5a). As shortening continued a low angle shear zone cuts through the Carbopol layer and ramps up through the overburden some 4 cm from the foreland edge of the Carbopol (Figure 5b). This shear zone links back to the leading edge of the sub-Carbopol imbricate thrust stack (Figure 5a). With continued shortening (12.5 cm, Figure 5c) deformation intensifies along the foreland-vergent system in our shale analog and the overburden, with a minor backthrust forming a conjugate to the steep foreland ramp in the system. After 20 cm shortening a hangingwall shortcut fault system develops in the main foreland-vergent thrust that links back through the Carbopol in an intense zone of shearing to the sub-Carbopol imbricate thrust stack (Figure 5d). Note that the shale analog above this thrust stack is highly deformed by an anastomosing shear zone that forms a broad backthrust to accommodate the uplift generated by this hinterland anticlinal stack, significantly increasing the taper angle (>12°, Figure 5c and 5d).

Incremental shear images allow us to track the progressive evolution of active deformations (Figure 6). After 6.5 cm shortening intense shearing is seen on one of the foreland-vergent thrusts in the hinterland, and shear strain is seen extending out along the base of the Carbopol layer from a zone of intense intra-Carbopol shear (Figure 6a). After 9.5 cm shortening activity waned on the hinterland fault systems and the low-angle shear zone cuts across and through the Carbopol linking to a steeper ramp in the overlying sediments (Figure 6b). With increasing shortening, intense shearing is seen both along the main foreland-vergent thrust system and a resurgence in activity along a complex backthrust system within the shale analog (Figure 6c and 6d). The frontal forethrust is cut by a minor hangingwall shortcut fault with activity along the original steeper fault segment in its footwall ceasing, generating a zig-zag structure linking down through the Carbopol (Figure 6d). Note that cumulative and incremental strain plots in Figures 5 and 6 show little to no strain in the Carbopol within the footwall of the main foreland-vergent thrust system (Figures 5d and 6d).

The final side view and a cross section through model 1 illustrate the structural features of this experiment (Figure 7a and 7b). These consist of: (1) A foreland breaking thrust sequence in the supra-Carbopol roof with evidence that our shale analog has been carried up and lubricated portions of the lower fault planes — some of the deeper portions of these thrust faults are very sharp, broadening up into wider shear zones that typify faults in sandbox models; (2) A major imbricate thrust stack developed below the Carbopol, and; (3) Development of a major backthrust system within the Carbopol and overlying strata, which served to increase the taper of the system to 17°. In the side view the minor hangingwall shortcut seen in the strain data is visible but this structure does not appear to be continuous across the entire model (Figure 7a and 7b). In the central cross section a relatively large proportion of the Carbopol was preserved (Figure 7b) but a comparison with the final side view (Figure 7a) reveals that some material has been lost during setting of the model. For example, the thickness of the Carbopol immediately in front of the imbricate thrust stack is less that that seen in the final side view, and the sagging of the roof near the foreland edge of this layer indicates some loss of material (Figure 7a and 7b).

Incremental strain imagery (Figure 8) and a cross-section (Figure 9) document the evolution and final geometry of model 2 (Table 1; Figure 4b). Although this experiment has less shortening than model 1 (see Table 1), Figures 6a–6c and 8a–8c are directly comparable and illustrate very different structural styles and amount of hinterland uplift. Model 2 consists of a series of forward-breaking thrusts with activity typically distributed across three to five individual structures during the entire model runtime (Figure 8). The strongest incremental shear is seen toward the foreland edge of the PDMS, and shear is seen at the sand-polymer interface although some of this could be noise due to smearing along the model sidewalls. No discernible shear is visible in the polymer layer toward its foreland edge (Figure 8d). Artifacts introduced by transport-parallel smearing are indicated on Figure 8d.

A cross section from the center of model 2 gives a much clearer picture of the structural style than the final sideview which is marred by significant smearing along the sidewall (Figures 8 and 9). In the sub-PDMS strata a short and rather complex thrust stack defines deformation, with some of the structures strongly impinging into the PDMS (Figure 9). Immediately above the PDMS there is a series of disharmonic structures that pass up into primarily foreland-vergent thrusts that are seen in the strain data (Figures 8 and 9). The main thrust is located right above the foreland edge of the PDMS layer, while minor backthrusts are seen above the hinterland edge of PDMS (Figure 9) and in the strain data (Figure 8). The taper in this model is lower, at 10.5° but this model had less shortening (Table 1). However, the taper in model 1 measured at 12.5 cm shortening was greater than that seen in model 2 at 15 cm shortening (compare Figures 5c and 9).

The evolution of model 3 with a thin Carbopol layer (Table 1; Figure 4c) is captured by cumulative and incremental maximum shear strain maps shown in Figures 10 and 11, respectively. A forward-breaking sequence of foreland-vergent thrusts is seen after 7 cm shortening with a zone of shear extending out along the Carbopol layer from the hinterland zone of intense shear in the cumulative shear strain data (S1 in Figure 10a). Incremental shear strain maps with lower clipped strain values clearly illustrate a zone of intense shear within the Carbopol linking to a ramp in the system cutting up through the overburden (Figure 11a). As shortening increased to 10.5 cm a clear asymmetric pop-up or box-like fold is developed at the leading edge of the system (S2 in Figures 10b and 11b). Incremental shear also shows the development of a hangingwall shortcut fault segment in this structure (Figure 11b). After 14 cm shortening shear along the Carbopol layer extends all the way to its foreland edge and ramps up to form a staircase system (S3 in Figures 10c and 11c). The hangingwall shortcut fault seen in Figure 11b linked to focus strain along this part of the system (Figures 10c and 11c). After 16 cm shortening deformation is now continuous and reaches up to the surface of the model along a complex fault system that uses the Carbopol layer as its main slip surface (Figures 10d and 11d). Incremental shear strain shows a waning of activity in the hinterland, and a steep backthrust system adjacent to the moving endwall (Figure 11).

The side view at the end of the model runtime and a cross section through the center of model 3 clearly illustrates the three main structural zones seen in shear strain imagery (S1–S3, Figure 12a and 12b). In the hinterland S1 consists of a series of closely spaced foreland-vergent thrusts in the supra-Carbopol strata underlain by an imbricate thrust stack (Figure 12). Moving outboard, S2 forms an almost box-fold geometry at shallow levels (Figure 12), and is separated from S1 by a tight syncline at shallow levels that is infilled by slumped material. At deeper levels, just above the Carbopol layer, the main foreland-vergent thrust of S2 is very sharp before diverging upward into a series of thrust strands (Figure 12b), that represent a part of the hangingwall shortcut fault seen in side views (Figures 10 and 11). At the foreland edge of the fold-thrust belt, structure S3 consists of a relatively symmetric pop-up structure located above the edge of the Carbopol and the taper measured is 11° (Figure 12). In this section most of the Carbopol was leached out during the setting process, although remnants remain, likely in areas where it had been thickened during thrusting (Figure 12b). This leaching process leads to weld development of strata once separated by our mobile shale analog (Figure 12a and 12b).

Cumulative and incremental shear strain imagery capture the side-view evolution of model 4 with two vertically offset but partially overlapping Carbopol layers (Figures 13 and 14; see Table 1 and Figure 4d for model setup details). The first group of structures to form are foreland-vergent thrusts in the supra-Carbopol overburden that are foreland breaking and define structural zone S1 (Figures 13a and 14a). Beneath the Carbopol, short thrust segments form and flatten into the lower detachment before ramping up through the overburden (C1 in Figure 13a). As shortening increases deformation propagates toward the foreland and shear is seen along the entire lower detachment layer (Figure 13b), and can be just seen ramping up in the incremental shear data (Figure 14b). The imbricate stack below the Carbopol steadily rises by basal accretion and incorporates heavily sheared Carbopol into the structure (C1 and C2 in Figure 13b). After 16.8 cm shortening a new main fault system begins to dominate linking up from the lower Carbopol detachment to the upper detachment along a ramp, as activity wanes in the hinterland (Figures 13c and 14c). A low amplitude anticline (structure S2) forms above the ramp and another fault related uplift (structure S3) at the foreland tip of this staircase fault system. More of the lower detachment is also incorporated into the deep imbricate stack system (C1–C3 in Figure 13c). At the end of the model run, major cumulative shear defines a staircase primary thrust detachment (Figure 13d). Incremental maximum shear strain maps shows that strain is localized along this structure with only minor strain seen in the hinterland (Figure 14d). Minor strains are seen in the fault bend fold of structure S2, possibly relating to flexural slip strain as this structure moves up the ramp and onto the flat (Figure 14d).

An oblique 3D view of the top surface, the final side view, and a section through the center of model 4 are shown in Figure 15, again consisting of three main structural zones, S1–S3. The lidar-generated 3D view of the model surface clearly shows these three zones from the oversteepened hinterland to undeformed flat-lying foreland (Figure 15a). The side view shows stacked lower Carbopol in the deep imbricate stack, foreland-vergent thrusts above this, and the major thrust detachment that links both Carbopol layers (Figure 15b). Above the ramp linking the two detachment layers, a tighter and more asymmetric fault-bend fold is separated by the upper detachment from a more open and symmetric fault-bend fold at shallow levels (S2 in Figure 15b). In the central cross section, we see that most of the Carbopol has been leached out during the setting process leaving just some remnants along the main fault system and imbricate stack but the structural geometries are preserved (Figure 15c and 15d). Pseudo welds now characterize the contacts between brittle layers that were once separated by the Carbopol (Figure 15c and 15d). The deeper fault-bend fold and overlying open anticline of S2 are now juxtaposed along these pseudo welds (Figure 15c and 15d). The thrusted green and white layers that lay below the upper detachment are now in faulted contact against the flat-lying strata in the footwall of the frontal thrust system, whereas the stratigraphy that lay above this upper detachment is continuous and thrust along what was the upper detachment and up the frontal ramp of S3 (Figure 15b–15d). A low taper of 9° is seen when linking the hinterland crest of S1 to the foreland ramp anticline of S3 but this increases to 14° when discounting S3 (Figure 15c).

In the Introduction section of this manuscript, we outlined five key questions we sought to answer in this study. In this section we revisit those questions briefly in the light of the results of the four new experiments presented previously.

Reber et al. (2015) describe the deformation process of Carbopol, after yield is exceeded, as a mixture of ductile flow, which is localized initially in the interstitial fluids (Figure 16), along with gradual fracturing of the microgel grains which also alters the ratio between the brittle and viscous phases along the shear zone. Soto et al. (2021a) state “At critical state, flow in shales generates an anastomosing network of highly sheared material, surrounding lenses in which previous structures are preserved. Thus, brittle and ductile structures may exist in mobile shales.” These two descriptions are complementary, and strongly suggest Carbopol is a good candidate for a mobile-shale analog by its behavior during yield. The degree to which brittle or ductile behavior dominates during yield in Carbopol is dependent on the percentage of the microgel grains in the Carbopol mixture (Reber et al., 2015), again similar to mobile shales whereby the degree of overpressured state, or the percentage of fluids present in the shale body (see grain- or fluid-supported flows, sensuWood, 2010), strongly controls the yield stress required for a shale to become mobile (e.g., Soto et al., 2021b).

In all of our physical models presented herein not only do we see shear zones developed in the brittle section, as expected for this type of DIC analyses of specular granular layers (Jagger and McClay, 2018), but we also see clear shear zones developed along and through our shale analog much as that described above, and by Reber et al. (2015) and illustrated in Figure 16 (see also Figures 5, 6, 10, 11, 13, and 14). Figures 17 and 18 shows more detailed views of model 1 and associated strain data sets derived by DIC analyses of the side views of this experiment. Grid deformations show strong and complex deformation within the Carbopol in the hinterland, which is heavily loaded by the moving endwall, with the gradual development of a foreland propagating shear zone within the Carbopol (Z1 and Z2 in Figures 17a, 17c, 17d and 18a).

Cumulative shear strain data also document the same pattern, with heavy deformation in the thickened and steepened Carbopol in the hinterland and a broad but discrete shear zone gradually propagating out from the hinterland high-strain zone that eventually ramps up into the overburden (Z1 and Z2 in Figures 17b, 17d, 17f and 18b). In all the shear strain data there is clear linkage between the tips of the sub-Carbopol structures and the base of the supra-Carbopol structures, through a broad region of shearing within our shale analog (Figures 17 and 18). Incremental shear data highlight which structures are active in model 1 at high shortening strains and here we see strain localized on the main foreland-vergent ramp-flat-ramp thrust and more minor strain along a passive roof backthrust (Figure 18c). It is important to note here that deformation is pervasive throughout the Carbopol layer, and not simply just slip at either the base or top of the detachment (Figures 17 and 18). Note that even in the complexly deformed hinterland “steep zone” of Carbopol the ductile shear zones are anastomosing and surround lenses of low strain (e.g., Figures 17b, 17d, 17f and 18b), that likely contain unbroken or fractured brittle microgel grains as demonstrated by Reber et al. (2015) (see Figure 16). Similarly, toward the foreland, there are considerable volumes of our shale analog with little to no strain adjacent to high strain zones (Figures 17f and 18). These are zones whereby yield strength has not been exceeded and thus the Carbopol remains immobile, like that seen to the left of the slip surface in Figure 16.

All the preceding points to intense deformation of our Carbopol mixtures, with complex shear zones forming where the yield strength has been exceeded, and the more ductile phase wrapping around lenses of intact microgel grains (Figure 18; Reber et al., 2015). In addition, portions of the shale analog undergo little or no strain where yield strength has not been exceeded, or where strain localization along shear zones bypassed portions of the Carbopol leaving them relatively undeformed (Figure 18b and 18c). As pointed out previously but bears reiterating, deformation of the Carbopol is not simply at its base or top but is, in places, pervasive through the entire mobile layer, and in other places the Carbopol remains almost undeformed (Figures 17 and 18). Strain analyses of model 1 reveal clear linkage between the tips of sub-Carbopol thrusts and the base of thrusts in the supra-Carbopol sequence, thus forming kinematically linked shortening structures (Figures 17 and 18), something that would be nearly impossible to visualize without the unique features of Carbopol and the use of DIC techniques.

Model 2 used the traditional salt analog, PDMS, in an experiment that had an almost identical setup to model 1 (Table 1). Strain data from this model demonstrated that the overburden synchronously deformed over a broad region, as bulk strain characterized deformation in the PDMS with increasing shortening (Figure 8). As noted earlier, no discernible shear zones can be observed within the PDMS, apart from noise generated by smearing along the sidewalls (Figure 8), a very different story to that seen in model 1 (Figure 17). This bulk deformation results in simultaneous activity along numerous forward-breaking thrusts in the overburden of model 2, rather than the discrete shear zones that cut through the Carbopol linking the brittle sub and supra-Carbopol packages seen in model 1 (Figure 17).

Figure 19a and 19b illustrates a comparison between sections from comparable models with Carbopol and PDMS (models 1 and 2, respectively), and clearly demonstrates the differences in structural styles between these two experiments with their very different mobile-layer behaviors. The lack of discrete shear zones within the PDMS layer is also consistent with passive marker deformation seen in the basal salt layer of previous physical model studies that displayed Poiseuille or asymmetric Poiseuille flow profiles in cross-sectional views (Bonini, 2003; Dooley et al., 2009, 2015; Santolaria et al., 2021), or in the more mobile units of a layered evaporite sequences (Cartwright et al., 2012; Weijermars et al., 2014).

Shortening in the supra-Carbopol sediments is accommodated along well-developed foreland-breaking thrusts (S1–S3, Figure 19a), with clear linkage between the tips of sub-Carbopol thrusts and the base of supra-Carbopol thrusts along a broad shear zone within the mobile layer (Figures 17, 18, and 19a). In contrast, bulk deformation of the salt analog is clearly illustrated in Figure 19b, with deformation encompassing the entire length of the model salt layer as gross inflation, combined with numerous small-displacement thrusts and folds in the overburden, and a lack of linkage of sub and supra-PDMS structures.

For completeness, two sections from a simple thrust wedge model with differing basal granular detachments are shown in Figure 19c and 19d. These sections illustrate a relatively simple forward-breaking thrust sequence, and a higher taper angle formed above a traditional sand detachment (Figure 19c), versus the microsphere detachment (Figure 19d). In both sections, there is structural continuity up through the entire overburden (Figure 19c and 19d), as opposed to partitioned and disharmonic deformation seen in our models with interbedded ductile detachments (Figure 19a and 19b).

In our models, imbricate thrust stacks characterize the deformation style beneath the shale-analog layer (Figures 7, 12, and 19a). These structures grow by basal and frontal accretion (see the model results in Bonini, 2001; Konstantinovskaya and Malavieille, 2011; Buiter, 2012; Graveleau et al., 2012; von Hagke and Malz, 2018), gradually jacking-up the hinterland portion of the model fold-thrust belt (Figure 19a). Broad and anastomosing shear zones link the tips of these imbricates across the shale analog to closely spaced foreland-breaking thrusts in the overburden, before a major low-angle shear zone propagates out into the Carbopol, eventually ramping up into the cover sediments, localizing deformation on a major ramp-flat-ramp thrust structure (Figures 5, 10, 17, and 18). In contrast, model 2 with a PDMS layer (Figures 8, 9 and 19b), formed a thrust stack in the sub-PDMS strata that was more fold-dominated and not as vertically stacked as those run with the shale analog (Figure 19a). Although the structures in this thrust stack heavily indent the PDMS, there is a lack of connectivity between the lower thrust stack and supra-PDMS units, with disharmonic structures seen just above the PDMS layer, and a lack of maturity of the fault systems in the supra-PDMS sequence (Figure 19b).

Our model 4, with dual model shale layers, is more complex, with shear enveloping the entire lower detachment outboard of the hinterland imbricate thrust stack, before ramping up at the tip of this lower layer and linking with a portion of the upper shale (Figures 13–15). This generated a spectacular, and geologically realistic, fault-bend fold below the upper detachment, with sharp stratal cutoffs along the fault plane, and decoupled from the open anticline-syncline pair in the shallow supra-Carbopol strata (Figure 15). Removal of the shale analog during the setting procedure resulted in artificial weld development that flank these structural cutoffs (Figure 15c and 15d).

The visualization of shear strain within our Carbopol models, especially incremental shear strains, allows us to track activity within the developing fold-thrust belts (Figures 6, 11, 14, and 18). Models 3 and 4 clearly demonstrate major shear developed along the most forward-breaking thrusts combined with concomitant but waning, activity on one or more hinterland thrusts (Figures 11 and 14), similar to results from discrete element modeling work by Hardy et al. (1998, 2009) and to results documented in Graveleau et al. (2012). Model 1 also shows this behavior (Figure 6b) but also demonstrates clear variations in passive-roof backthrust activity within the upper portions of the Carbopol layer, as this structure lights up in high strain colors to accommodate episodic vertical growth of the sub-Carbopol imbricate thrust stack due to basal accretion as the wedge advances (Figures 6 and 18).

An additional feature of our Carbopol models is the presence of very sharp fault surfaces in the lower portions of the supra-Carbopol sequences (Figures 7, 12, 15c, and 19b), as opposed to the broader shear zones that typify thrust faults in granular materials (e.g., compare faulting in S2 and S3 in Figure 12, with the domain between S1 and S2-S3 in Figure 19b). As indicated on Figure 7a this is attributed to lubrication, and local seepage, provided by the shale analog on which the lower suprashale strata sat prior to shortening. Similar “sharp faults” were documented by Dooley and Hudec (2020) whereby PDMS infiltrated pore spaces, allowing for efficient reactivation of these structures during inversion (see Figure 15 of Dooley and Hudec, 2020). Sharper faults are also seen in the simple thrust wedge model with a microbead basal detachment in Figure 19d, due to the finer grained nature of these materials (typically less than 100–150 μm; Reber et al., 2020).

Based on our experimental results, and discussion previously, our mobile-shale analog serves as a very efficient detachment even when the detachment layer or layers are very thin relative to the overburden thickness (Figures 10 and 13). As the yield strength of our shale analog is gradually exceeded, these layers light up brightly as strain propagates rapidly toward the foreland along these thin detachments (Figures 13 and 14). With a single thin detachment layer, strain ramped up at the foreland edge of the Carbopol layer into the overburden to form an asymmetric popup (Figure 12). With dual thin detachments that are offset in horizontal and vertical space things are more complex (Figures 13 and 14). Here, the shear strain ramps up at the tip of the lower mobile layer through the brittle overburden and flattens onto a portion of the upper mobile layer (Figures 13,–15). A distinct, and geologically realistic, fault-bend fold forms along this ramp-flat system below the upper detachment while above it the structures are largely decoupled but harmonic and form an open and symmetric anticline-syncline pair (Figure 15). Incremental strain maps from this model also demonstrate strain patterns that hint at flexural-slip deformation within the upper Carbopol layer and the strata directly above this layer (Figure 14d).

Although only one model with more than one shale-analog layer is presented here, we believe that our model 4 (Figures 13–15) represents a very promising start in investigating deformation styles in a multilayer containing a suitable mobile shale analog under contraction. Detachments realistically ramp up at the lateral tips of mobile layers and flatten onto mobile layers higher in the stratigraphy, with realistic structural cutoffs of hangingwall strata, facilitating strong decoupling in the overburdens (Figure 15c and 15d), and expanding on previous studies with multiple detachment levels consisting of PDMS and/or granular materials (Couzens-Schultz et al., 2003; Konstantinovskaya and Malavieille, 2011; Graveleau et al., 2012; Santolaria et al., 2015; Pla et al., 2019). Future modeling will continue to investigate deformation of multilayers under contraction, to further elucidate structural styles associated with these types of sequences, e.g., shale-cored detachment folds, zigzag thrusts, and shale diapir nucleation and growth (Soto et al., 2021a; Soto and Hudec, 2023; Soto et al., 2024), and to assess the role of flexural slip deformation in these multilayers.

An outcrop-scale example from the Betic Cordillera (Western Mediterranean) nicely illustrates the nature of deformation that occurs in a multilayer deformed by shortening (Figure 20). The example comes from an outcrop in central Betics, with an Upper Jurassic succession (Kimmeridgian-Tithonian; Díaz de Neira et al., 1992) consisting of a competent basal limestone sequence, overlain by a weaker middle marly unit, passing up into a heterogeneous unit consisting of variable-thickness layers of marls and marly limestones. This outcrop is situated in the upright limb of a major westsouthwest–eastnortheast trending anticline with southerly vergence (Díaz de Neira et al., 1992). The general geologic setting of this Alpine orogen and the structure of its external domains are documented in Flinch and Soto (2022).

The general structure of the outcrop consists of a succession of folds of variable amplitude and associated thrusts, with ramps developed in the more competent limestone layers and long flats along the less competent marly levels (Figure 20b). The lower part of the sequence is cored by the thicker, more competent, limestone layers and the structural style is southwest-vergent asymmetric folds, similar to that found in the sub-Carbopol sequences in models 1, 3, and 4 (Figure 19a). Above this sequence, the middle unit of less competent marls shows thickening toward the fold cores, whose geometry and position is decoupled with respect to the lower structures. This type of geometry is consistent with the structural style developed by the Carbopol layer in our models (see the thickening and style of deformation of the Carbopol layer in Figures 7 and 12b). In the overlying succession of marls and marly limestones deformation consists of asymmetric folds, disharmonic with respect to the structures in the intermediate marl layer, and with several cases of detachment folds and others thrust-related folds (Figure 20b).

Structural styles of the field example complement many of the findings of our new experimental models with an intermediate Carbopol layer (models 1 and 3 in Figures 7 and 12). For example: (1) the structural decoupling exerted by a weaker interlayer separating the shortening deformations below and above the mobile layer (Figure 19a); (2) the development of an antiformal stack beneath the Carbopol layer (Figure 19a); and (3) the creation of ramp-and-flat thrusts with associated folds in the sequences above the Carbopol layer (Figures 7 and 12).

On a smaller scale, the deformation observed within the weak marls in the core of folds broadly resemble the deformation within the Caropol layers in our models (Figure 18). Here we find evidence of penetrative deformation that typically characterizes extreme deformation in mobile shales (Soto et al., 2021a, 2021b), with the simultaneous development of foliations, tensional veins filled with calcite, and a complete destruction of the original sedimentary fabric (see the network defined by the fracture cleavage and the abundant cross-cutting tensional veins, as shown in Figure 20c and 20d). These observations are also consistent with the pattern of deformation and anastomosing shear zones found within the Carbopol layer in our models (Figures 17 and 18).

The similarity between observations from this natural example to those found in our Carbopol models confirm that the type of deformation seen inside mobile shales (see the deformation scheme of Soto et al., 2021a, their Figure 6b), may help to explain aspects of these types of complex structures seen in seismic profiles in mobile-shale basins (van Rensbergen and Morley, 2003; Soto et al., 2021b; Hudec et al., 2023; Morley et al., 2023; Soto and Hudec, 2023). For larger scale examples of mobile shale deformation, the reader is referred to our companion paper that presents results from a 3D seismic study in the northern Gulf of Mexico (Soto et al., 2024), with comparison to our preliminary modeling results presented herein.

Recent papers by Butler et al. (2020) and Rowan et al. (2022) provide thoughtful reviews of the inherent problems that are present in all physical models, particularly the lack of lithologic variability within the sediments that typically constitute the supradetachment overburdens. The addition of Carbopol as pre and syn-shortening layers to simulate natural overpressured shales that we know are present in many, if not all, basins may allow the generation of more realistic shortening geometries, including multiple stacked detachments commonly seen in shortened basins (Davison, 2021; Rowan et al., 2022; Jadoon et al., 2024). Dooley and Hudec (2024), in their study of the Sureste Basin (southern Gulf of Mexico), noted that the geometries of their model salt structures are quite realistic when compared with structures seen in seismic data from this basin but the weak deformation geometries within the thick, and strong, sedimentary packages between the model diapirs do not match natural geometries. The addition of Carbopol into similarly styled models would have the potential to rectify this imbalance in structural styles.

Our preliminary modeling results and comparison with an outcrop example of a shortened multilayer suggest that Carbopol, being a yield stress fluid, is a good analog for the behavior of mobile shales during shortening. Our shale analog only deforms in a visco-plastic fashion when its yield strength is exceeded. This allows for portions of the material to remain undeformed and immobile, while other portions deform, forming broad anastomosing shear zones where the shale analog is thick and more discrete shear zones where the mobile layer is thin. The hybrid nature of the Carbopol mixtures makes the visualization of strains within these layers possible by DIC analyses, providing a powerful quantitative and graphical tool for strain analysis and allowing for the full story of linkage between detachment layer and brittle overburden to be revealed.

Much more work needs to be done with Carbopol, including further testing in models run under contraction but also gravity-driven deformation such as that seen in delta systems along with extensional deformation. One long-term goal is to introduce this material into complex 3D physical models investigating, for example, the role of preexisting shale diapirs and their influence on contractional deformation. In addition, we envision studies investigating the possibility of producing mud volcanoes and their piercement mechanisms by using extremely weak Carbopol mixtures under certain conditions, as well as further investigation of mixed salt-shale systems such as those seen in parts of the Gulf of Mexico and other locations.

We feel confident that our current models and future experiments will provide important graphical structural templates as an aid for seismic interpretation in basins containing mobile shales or to those containing mobile shale and salt detachments. To conclude, we believe Carbopol is a powerful mobile-shale analog and that this preliminary study opens new and exciting future experimental modeling directions by introducing it into a variety of traditional 2D and 3D sandbox models.

T. Dooley thanks N. Ivicic, B. Williamson, R. Lucero, and E. Wheeler for assistance in his modeling laboratories. We thank D. Jerolmack for early discussions on the behavior of mobile shales and what type of material is a suitable analog for laboratory use. We gratefully acknowledge J. Duarte, K. McClay, C. Burberry, and P. Prince for their constructive and enthusiastic comments that helped to greatly improve the manuscript. C. Morley and M. Tingay are thanked for editorial handling and comments on the manuscript. This study was funded by the Applied Geodynamics Laboratory Industrial Associates program. Additional funding was provided by the Jackson School of Geosciences. Publication was authorized by the director of the Bureau of Economic Geology, Jackson School of Geosciences, The University of Texas at Austin.

Data associated with this research are available and can be obtained by contacting the corresponding author.

Biographies and photographs of the authors are not available.