Itzhak Hamdani

Einat Aharonov

Jean-Arthur Olive

Stanislav Parez

Zohar Gvirtzman

Initiating salt tectonics by tilting: Viscous coupling between a tilted salt layer and overlying brittle sediment Authors

Key points 1) We simulate and predict the formation of normal faults due to tectonic tilting of the shelf/slope region of a salt basin 2) The viscosities and thicknesses of the salt and the sediment overburden layers control fault initiation 3) The predicted faulting process agrees with observations of location and temporal evolution of normal faulting in the Levant basin margin Abstract Salt basins often exhibit a shelf/slope region of extension and a deeper domain of contraction. The upslope normal faults in such salt tectonics systems are often associated with the pinchout edge of the buried salt layer. Although the spatial correlation between the normal faults and the salt pinchout was previously observed, the mechanism was not fully explained. The Levant basin, which is a young and mildly deformed salt basin with a relatively thin overburden, provides an opportunity to analyze a simple salt tectonics system driven by basin margin tilt, and the formation of the normal faults. This work presents analytical and numerical modeling of the coupled viscous salt and overlying visco-plastic sediment layer. Results suggest the viscosities of both the sediment overburden and salt, as well as their thicknesses, control the deformation of the coupled layers. The visco-plastic deformation mechanism explains quantitatively the position of faulting and observations of the temporal evolution of brittle deformation in the Levant basin margin. It predicts that the largest stress in the overburden arises above the salt edge, driving normal faulting at that location. Our model also places quantitative constraints on the effective viscosity of the overburden (>2*10 20 Pa s), which is consistent with experimentally-determined creep laws.

Our results contribute to the understanding of halokinematics in salt basins during the early stages of deformation and will allow better assessment of geological hazards related to salt related deformation.

Introduction

The very weak rheology of buried salt layers makes them susceptible to rapid viscous flow when subjected to pressure gradients. This in turn commonly triggers deformation of any overlying rocks that are mechanically coupled to the salt. Such thin-skinned deformation, limited to the salt and its overburden, is broadly referred to as salt tectonics. Salt basins are widely studied using analogue experiments (e.g., Brun & Fort, 2012;[START_REF] Brun | Salt rollers: Structure and kinematics from analogue modelling[END_REF][START_REF] Dooley | Initiation and growth of salt-based thrust belts on passive margins: Results from physical models[END_REF][START_REF] Fort | Salt tectonics on the Angolan margin, synsedimentary deformation processes[END_REF]Ge, Rosenau, et al., 2019;[START_REF] Jackson | Retrospective salt tectonics[END_REF]Vendeville & Jackson, 1992a[START_REF] Vendeville | The fall of diapirs during thin-skinned extension[END_REF][START_REF] Warsitzka | Analogue experiments of salt flow and pillow growth due to basement faulting and differential loading[END_REF]and referenced therein) and numerical modeling (e.g., [START_REF] Aiban | An investigation of salt tectonic structural styles in the Scotian Basin, offshore Atlantic Canada: 1. Comparison of observations with geometrically simple numerical models[END_REF][START_REF] Allken | 3D numerical modelling of graben interaction and linkage: a case study of the Canyonlands grabens[END_REF][START_REF] Baumann | 3D numerical modelling of salt tectonics[END_REF][START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF][START_REF] Gemmer | Dynamic modelling of passive margin salt tectonics: Effects of water loading, sediment properties and sedimentation patterns[END_REF][START_REF] Gradmann | Coupled fluid flow and sediment deformation in margin-scale salt-tectonic systems: 2. Layered sediment models and application to the northwestern Gulf of Mexico[END_REF]and referenced therein).

The majority of salt tectonic studies focus on hydrocarbon-rich and highly deformed basins (e.g., the Gulf of Mexico, the Atlantic margins), where salt is buried under thick overburden, and remobilized from its original stratigraphic location (e.g., [START_REF] Rowan | Structural styles and evolution of allochthonous salt, central Louisiana outer shelf and upper slope[END_REF] through diapirs, forming complex structures and complicated stratigraphy (e.g., [START_REF] Dooley | Initiation and growth of salt-based thrust belts on passive margins: Results from physical models[END_REF][START_REF] Gradmann | Factors controlling the evolution of the Perdido Fold Belt, northwestern Gulf of Mexico, determined from numerical models[END_REF][START_REF] Letouzey | Evolution of Salt-Related Structures in Compressional Settings[END_REF][START_REF] Talbot | Extrusions of Hormuz salt in Iran[END_REF]. Young salt basins with relatively thin overburden, in contrast, provide a simpler setting for analysis: the original stratigraphic location of the salt is known, its initial geometry can be restored, the spatial and temporal evolution of deformation can be inferred, and specific deformation mechanisms can be assessed. Previous studies of these salt basins, many of them conducted along the margin of the eastern Mediterranean Sea, have shown an up-slope (landward) domain of extension, an intermediate zone of translation, and a down-slope (seaward) domain of contraction [START_REF] Allen | Gravity-driven deformation of a youthful saline giant: the interplay between gliding and spreading in the Messinian basins of the Eastern Mediterranean[END_REF][START_REF] Cartwright | Initiation of gravitational collapse of an evaporite basin margin: The Messinian saline giant, Levant Basin, eastern Mediterranean[END_REF][START_REF] Gradmann | Salt tectonics off northern Israel[END_REF]Netzeband et al., 2006). This geometry was inferred to have formed by basin-ward motion of the salt and its overburden, and to be driven by gravity gliding induced by tilting and/or spreading induced by differential-loading of sediments (see Section 3) [START_REF] Allen | Gravity-driven deformation of a youthful saline giant: the interplay between gliding and spreading in the Messinian basins of the Eastern Mediterranean[END_REF]Evans et al., 2020;Evans & Jackson, 2020;[START_REF] Oppo | Leaky salt: pipe trails record the history of cross-evaporite fluid escape in the northern Levant Basin, Eastern Mediterranean[END_REF] One of the noticeable features in the extensional domain of salt-basin margins are normal faults associated with the salt pinch-out, as documented in the Gulf of Mexico [START_REF] Rowan | Salt-related fault families and fault welds in the northern Gulf of Mexico[END_REF], offshore Israel [START_REF] Gradmann | Salt tectonics off northern Israel[END_REF][START_REF] Katz | Spatial and temporal relation of submarine landslides and faults along the Israeli continental slope, eastern Mediterranean[END_REF], offshore Egypt [START_REF] Allen | Gravity-driven deformation of a youthful saline giant: the interplay between gliding and spreading in the Messinian basins of the Eastern Mediterranean[END_REF] and more. It has been suggested that basinwards salt flow initiates these faults [START_REF] Allen | Gravity-driven deformation of a youthful saline giant: the interplay between gliding and spreading in the Messinian basins of the Eastern Mediterranean[END_REF]Bertoni & Cartwright, 2006, 2015;[START_REF] Cartwright | Direct calibration of salt sheet kinematics during gravity-driven deformation[END_REF][START_REF] Gradmann | Salt tectonics off northern Israel[END_REF][START_REF] Kirkham | 3D kinematics of a thick salt layer during gravity-driven deformation[END_REF]. These studies focus mainly on salt flow and have treated the faulting in the overburden as a trivial result of the flow, that causes a weak overburden to collapse [START_REF] Cartwright | Initiation of gravitational collapse of an evaporite basin margin: The Messinian saline giant, Levant Basin, eastern Mediterranean[END_REF]). Yet, the mechanics that couple salt flow and faulting in this setting has not been studied thoroughly.

The present work studies how regional-scale tectonic tilting initiates normal faults above the landward edge of a seaward-dipping salt layer. We will clarify several basic questions: 1) the conditions for normal faulting at a tilted basin margin (effects of tilt angle, rheology, sediment and salt layer thicknesses); 2) position and timing of fault initiation and progress, and relation to the salt pinch-out; 3) rate of faulting and connection to the salt flow rate.

After investigating these questions by analytical and numerical modelling, we compare our results to a relatively simple geological case study in the southernmost Mediterranean Sea (Levant Basin). The advantage of working on the Mediterranean salt giant is its relatively young age (~5.5 Ma) and shallow burial (a few hundred meters in deep basins) that allows structural reconstructions [START_REF] Hübscher | Global look at salt giants[END_REF]. In the eastern Mediterranean and particularly in the Levant Basin, salt is still in its original stratigraphic location and its internal stratification is preserved (Bertoni & Cartwright, 2006; Z. [START_REF] Gvirtzman | Intense salt deformation in the Levant Basin in the middle of the Messinian Salinity Crisis[END_REF][START_REF] Hübscher | Global look at salt giants[END_REF]. In the southern part of the Levant, offshore Israel, no deep-seated faulting is currently active, and the observed faults are solely related to the underlying salt [START_REF] Garfunkel | Constrains on the origin and history of the Eastern Mediterranean basin[END_REF][START_REF] Gvirtzman | Reactivation of the Levant passive margin during the late Tertiary and formation of the Jaffa Basin offshore central Israel[END_REF][START_REF] Gvirtzman | Inland jump of the Arabian northwest plate boundary from the Levant continental margin to the Dead Sea Transform[END_REF]. These circumstances make the Levant Basin an excellent natural laboratory for studying the early stage of salt tectonics [START_REF] Hübscher | Global look at salt giants[END_REF]. Our work will relate the formation of faults to coupled flow of salt and overlying sediments. This in turn constrains rheological properties of these layers and sheds light on the general geological and mechanical evolution of young salt basins with relatively thin overburden.

Faulting above a tilted salt layer: the Levant Basin example

The Levant Basin (Figure 1) and its passive margin developed during several episodes of continental rifting between the Late Permian and the Early Jurassic [START_REF] Bein | A Mesozoic fossil edge of the Arabian Plate along the Levant coastline and its bearing on the evolution of the Eastern Mediterranean[END_REF][START_REF] Freund | The Triassic-Jurassic structure of Israel and its relation to the origin of the Eastern Mediterranean[END_REF][START_REF] Garfunkel | Constrains on the origin and history of the Eastern Mediterranean basin[END_REF][START_REF] Garfunkel | Permian-early Mesozoic tectonism and continental margin formation in Israel and its implications for the history of the Eastern Mediterranean[END_REF]. Subsequently, post-rift conditions prevailed for ∼100 Myr [START_REF] Bar | Accelerated subsidence and sedimentation in the levant basin during the late tertiary and concurrent uplift of the arabian platform: Tectonic versus counteracting sedimentary loading effects[END_REF][START_REF] Garfunkel | Permian-early Mesozoic tectonism and continental margin formation in Israel and its implications for the history of the Eastern Mediterranean[END_REF][START_REF] Gvirtzman | The transformation of southern Israel from a swell to a basin: Stratigraphic and geodynamic implications for intracontinental tectonics[END_REF] until the beginning of the Santonian (∼86 Ma), when the "Syrian Arc" fold belt started forming alongside Tethys closure [START_REF] Eyal | Tectonic analysis of the Dead Sea Rift Region since the Late-Cretaceous based on mesostructures[END_REF][START_REF] Freund | The Triassic-Jurassic structure of Israel and its relation to the origin of the Eastern Mediterranean[END_REF][START_REF] Krenkel | Der Syrische Bogen[END_REF][START_REF] Reches | the Structure of a Monocline in the Syrian Arc System, Middle East-Surface and Subsurface Analysis[END_REF][START_REF] Sagy | m.y. of folding migration: New perspective on the Syrian arc from Levant Basin analysis[END_REF][START_REF] De Sitter | Structural development of the Arabian Shield in Palestine[END_REF]. In the Late Eocene (∼35 Ma), coeval with the final closure of the Tethys and breakoff of Arabia from Africa, the Levant continental margin was reactivated (Z. [START_REF] Gvirtzman | Reactivation of the Levant passive margin during the late Tertiary and formation of the Jaffa Basin offshore central Israel[END_REF][START_REF] Gvirtzman | Inland jump of the Arabian northwest plate boundary from the Levant continental margin to the Dead Sea Transform[END_REF]; northern Arabia and Eastern Africa were uplifted [START_REF] Adams | Dating the terminal Tethyan event[END_REF][START_REF] Bar | Accelerated subsidence and sedimentation in the levant basin during the late tertiary and concurrent uplift of the arabian platform: Tectonic versus counteracting sedimentary loading effects[END_REF][START_REF] Bar | The uplift history of the Arabian Plateau as inferred from geomorphologic analysis of its northwestern edge[END_REF]Buchbinder & Zilberman, 1997a; Z. [START_REF] Gvirtzman | Retreating late tertiary shorelines in Israel: Implications for the exposure of north Arabia and levant during Neotethys closure[END_REF] and sediment supply to the Levant basin increased significantly [START_REF] Steinberg | Origin and nature of the rapid late Tertiary filling of the Levant Basin[END_REF]. These late Cenozoic processes changed the paleogeography of the Middle East region and altered the nature of sedimentation in the Levant Basin from hemipelagic carbonates to mainly siliciclastics (Z. [START_REF] Gvirtzman | Reactivation of the Levant passive margin during the late Tertiary and formation of the Jaffa Basin offshore central Israel[END_REF]. About 6 Myr ago, during Messinian times, when the Atlantic-Mediterranean gateways narrowed, the entire Mediterranean Sea became a hypersaline basin and a km-thick salt layer was deposited on its floor (e.g., [START_REF] Hsü | The origin of the Mediterranean evaporites[END_REF]Ryan & Hsü, 1973). This dramatic event is termed the Messinian salinity crisis. When normal marine conditions resumed in the beginning of the Pliocene (5.3 Ma, Roveri et al., 2014 and references therein), younger sediments started burying the salt layer at the bottom of the Mediterranean.

In the Levant Basin, the Messinian salt reaches a thickness of ~1.8 km in the deep basin, and abruptly thins near the continental margin. Below the present day continental slope the salt wedges out to a thickness of a few tens of meters or below seismic resolution (Bertoni & Cartwright, 2006;[START_REF] Gradmann | Salt tectonics off northern Israel[END_REF]Z. Gvirtzman et al., 2013[START_REF] Gvirtzman | Intra-Messinian truncation surface in the Levant Basin explained by subaqueous dissolution[END_REF][START_REF] Katz | Spatial and temporal relation of submarine landslides and faults along the Israeli continental slope, eastern Mediterranean[END_REF][START_REF] Mart | Some Depositional Patterns at Continental Margin of Southeastern Mediterranean Sea[END_REF]Netzeband et al., 2006;[START_REF] Ryan | The nature and distribution of Messinian erosional surfaces -Indicators of a several-kilometer-deep Mediterranean in the Miocene[END_REF]. The Plio-Quaternary sediments offshore Israel mainly consist of Nile-derived clay-rich siliciclastic sediments that were transported by alongshore currents from the Nile delta to Sinai and Israel (Almagor & Garfunkel, 1979;[START_REF] Ben-Gai | Post-Messinian evolution of the Southeastern Levant Basin based on two-dimensional stratigraphic simulation[END_REF][START_REF] Buchbinder | Sequence stratigraphy of Miocene-Pliocene carbonate-siliciclastic shelf deposits in the eastern Mediterranean margin (Israel): effects of eustasy and tectonics[END_REF][START_REF] Garfunkel | Geology and structure of the continental margin off northern Israel and the adjacent part of the Levantine Basin[END_REF][START_REF] Gvirtzman | The Late Tertiary of the coastal plain and continental shelf of Israel and its bearing on the history of the Eastern Mediterranean. Initial Reports of the Deep Sea Drilling Project[END_REF]). This transport system built the continental shelf and the present-day continental slope. The slope (Figure 1) is dissected by normal faults [START_REF] Almagor | Salt-controlled slumping on the Mediterranean slope of central Israel[END_REF]Almagor & Garfunkel, 1979;[START_REF] Garfunkel | Large-scale submarine rotational slumps and growth faults in the Eastern Mediterranean[END_REF][START_REF] Garfunkel | Geology and structure of the continental margin off northern Israel and the adjacent part of the Levantine Basin[END_REF][START_REF] Mart | Quaternary tectonic patterns along the continental margin of the southeastern Mediterranean[END_REF], which develop above the salt wedge (Figure 2), displacing its top (red reflector) but not its base (yellow reflector). The normal faults form above, and mostly westward, of the salt pinch-out, which is where the salt thickness abruptly reduces from hundreds to tens of meters. i.e., the point of merging between the bottom and top of the salt layer, marked in Figure 2 by an arrow. Figure 1 depicts the salt pinch-out in map view as a purple line, as inferred from seismic data (following Bertoni and Cartwright (2006) and [START_REF] Katz | Spatial and temporal relation of submarine landslides and faults along the Israeli continental slope, eastern Mediterranean[END_REF]).

In addition, Figure 1 shows that indeed, as shown in the representative cross section of Figure 2, all along the Levant basin a belt of normal faults forms above or slightly westward of the salt pinch out line. Exceptions are near the Palmahim and Dor disturbances, which involve mega slides of the continental slope [START_REF] Almagor | Salt-controlled slumping on the Mediterranean slope of central Israel[END_REF][START_REF] Garfunkel | Large-scale submarine rotational slumps and growth faults in the Eastern Mediterranean[END_REF]Garfunkel et al., 1979). We shall not concentrate on the disturbances and will instead focus on the formation of the faults in the more typical margin configuration shown in Figure 2. [START_REF] Elfassi | Chronology of post-Messinian faulting along the Levant continental margin and its implications for salt tectonics[END_REF] showed that faulting was initiated about 2.7 Myr after the end of salt deposition, in the Gelasian (2.6-1.8 Ma). This timing suggests that salt tectonics was initiated by tilting and not by differential loading of sediments, because it coincides with coastal uplift [START_REF] Bar | The uplift history of the Arabian Plateau as inferred from geomorphologic analysis of its northwestern edge[END_REF] and predates the propagation of the continental shelf to its present location near the salt wedge [START_REF] Elfassi | Chronology of post-Messinian faulting along the Levant continental margin and its implications for salt tectonics[END_REF]. Similar conclusions were derived by [START_REF] Cartwright | Initiation of gravitational collapse of an evaporite basin margin: The Messinian saline giant, Levant Basin, eastern Mediterranean[END_REF] who reported a tilt angle of ~1° in central Israel and increasing northwards and recently, by Evans and Jackson (2020) who worked north of this area.

Rheological descriptions of salt and sediment flow

There are two end-member salt flow scenarios, in response to driving forces (see for example, Brun & Fort, 2011;Hudec & Jackson, 2007;[START_REF] Schultz-Ela | Excursus on gravity gliding and gravity spreading[END_REF]): (a) Gravity spreading due to a pressure gradient produced by differential loading. In this case salt flows from high to low pressure zones (thick and thin overburden, respectively). (b) Gravity gliding of a tilted salt layer. In this case, salt flows due to an elevation gradient from high to low regions. In nature, both driving forces usually coexist. These naturally occurring forces produce viscous salt flow in the subsurface which may lead to brittle deformation of the overburden. The first end-member was examined in detail by [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF]. Here we focus on the second end-member, gravity gliding due to elevation differences (tilting) while assuming a constant thickness of the overburden.

In order to analytically and numerically explore how salt flow creates faulting in overlaying sediments, we first need to decide upon the correct rheological description of both layers. This rheology will play a critical role in the problem of fault initiation in a brittle layer, driven by salt flow underneath.

Under typical geologic conditions, most evaporites deform visco-elastically [START_REF] Weijermars | Rheological and tectonic modeling of salt provinces[END_REF], where the elastic component of their deformation is only manifested on short time scales and is commonly ignored when considering geological deformation at slow strain rates and for geological time scales. In the Levant basin, the nearly 2-km-thick salt sequence contains 10-20 distinct clay layers, 1-5 m thick each [START_REF] Feng | Lithology and characteristics of the Messinian evaporite sequence of the deep Levant Basin, Eastern Mediterranean[END_REF]Z. Gvirtzman et al., 2013), that may change the rheology somewhat. But for simplicity, we will model the salt layer as composed of pure halite, which flows like a purely viscous, fluid-like material (Hudec & Jackson, 2007).

We further assume that salt is a Newtonian fluid with a constant viscosity, 𝜂. Depending on the deformation mechanism and the environment, salt viscosity dependence on stress, strain-rate, and temperature, may need to be considered (e.g., [START_REF] Li | Rheology of rock salt for salt tectonics modeling[END_REF][START_REF] Urai | Weakening of rock salt by water during long-term creep[END_REF][START_REF] Urai | Flow and Transport Properties of Salt Rocks[END_REF]. Here we neglect these dependencies: at the low strain-rates and stresses that arise in our modeled system (as shown in Figure 6 below) the viscosity of halite is predicted experimentally and theoretically to be nearly constant, depending on grain size and temperature (see Figure 2 in [START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF]. For such slow strain rates and low stresses, halite creep is predicted to be dominated by pressure solution [START_REF] Li | Rheology of rock salt for salt tectonics modeling[END_REF][START_REF] Spiers | Experimental determination of constitutive parameters governing creep of rocksalt by pressure solution[END_REF], with a viscosity that does not depend on stress (e.g., [START_REF] Li | Rheology of rock salt for salt tectonics modeling[END_REF].

The dependence on temperature can be shown to modulate viscosity by at most a factor of 2, comparing conditions at the top of the salt to those at the bottom of our salt layer (assuming Eq. 7 in van [START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF]). Thus, we argue that a Newtonian flow law is a reasonable 1 st -order approximation for salt deformation in our system conditions: such an approximation does not subtract from the physical validity of the study, and yet allows us to illuminate and analytically probe fundamental processes that occur during salt flow. Simple numerical tests (not presented) suggest that using a temperature-dependent viscosity slightly affects the channel flow profile, a variation that adds mostly a second order effect for faulting.

The value of the effective Newtonian salt viscosity has been assessed through experimental and modeling studies (e.g., [START_REF] Karcz | Stability of a sodium chloride indenter contact undergoing pressure solution[END_REF][START_REF] Karcz | Deformation by dissolution and plastic flow of a single crystal sodium chloride indenter: An experimental study under the confocal microscope[END_REF][START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF][START_REF] Mukherjee | Viscosity estimates of salt in the Hormuz and Namakdan salt diapirs, Persian Gulf[END_REF][START_REF] Nalpas | The structural evolution of the Messinian evaporites in the Levantine Basin[END_REF][START_REF] Weinberger | Mechanical modeling and InSAR measurements of Mount Sedom uplift, Dead Sea basin: Implications for effective viscosity of rock salt[END_REF], to range between 10 -10 𝑃𝑎 ⋅ 𝑠 for varying grain sizes and temperatures. Salt viscosity under 1-2 km burial load is often assumed to be ~10 𝑃𝑎 ⋅ 𝑠 (e.g., [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF][START_REF] Gradmann | Numerical modelling study of mechanisms of midbasin salt canopy evolution and their potential applications to the Northwestern Gulf of Mexico[END_REF].

We next consider the stress state and rheology of the sediment layer overlying the salt: Previous salt tectonics studies assessed stress in the overburden in various ways. Brun and Fort (2012) estimated the stress within the sediment using mechanics similar to that of a landslide (Klar et al., 2011 and references therein), with a detachment surface marking the base of the salt. Another approach is to consider the shear tractions that the flowing viscous salt exerts on the bottom of the brittle overburden. This was done by [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF] and [START_REF] Gradmann | Factors controlling the evolution of the Perdido Fold Belt, northwestern Gulf of Mexico, determined from numerical models[END_REF] to assess overburden stability. Here we further investigate this latter approach.

The overlying sediment layer is the marly Yafo formation [START_REF] Ben-Avraham | The structure and tectonic setting of the levant continental margin, Eastern Mediterranean[END_REF][START_REF] Ben-Gai | Post-Messinian evolution of the Southeastern Levant Basin based on two-dimensional stratigraphic simulation[END_REF][START_REF] Gvirtzman | The Late Tertiary of the coastal plain and continental shelf of Israel and its bearing on the history of the Eastern Mediterranean. Initial Reports of the Deep Sea Drilling Project[END_REF]) that was penetrated by many wells in the coastal plain, where it is composed of ~16% carbonates, 55% clays and ~29% quartz (G. [START_REF] Gvirtzman | The Late Tertiary of the coastal plain and continental shelf of Israel and its bearing on the history of the Eastern Mediterranean. Initial Reports of the Deep Sea Drilling Project[END_REF]. The lateral thinning of the sand layers towards the continental margin (G. [START_REF] Gvirtzman | The Late Tertiary of the coastal plain and continental shelf of Israel and its bearing on the history of the Eastern Mediterranean. Initial Reports of the Deep Sea Drilling Project[END_REF] indicates that in the area studied here the clay portion of the Yafo formation most probably exceeds 70%. Such sediments are likely to behave in a visco-plastic manner that combines viscous deformation and fracturing. Following many previous experimental and modeling works (e.g., [START_REF] Aharonov | Interaction between pressure solution and clays in stylolite development: insights from modeling[END_REF][START_REF] Baker | Pressure solution and hydrothermal recrystallization of carbonate sediments -An experimental study[END_REF][START_REF] Dewers | Influences of clay minerals on sandstone cementation and pressure solution[END_REF][START_REF] Geng | Time and Temperature Dependent Creep in Tournemire Shale[END_REF][START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF][START_REF] Hellmann | Experimental pressure solution compaction of chalk in aqueous solutions Part 1 . Deformation behavior and chemistry. Water-Rock Interactions, Ore Deposits, and Environmental Geochemistry: A Tribute to David A[END_REF][START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF][START_REF] Safaricz | Pressure solution in chalk[END_REF][START_REF] Zhang | Compaction creep of wet granular calcite by pressure solution at 28°C to 150°C[END_REF], we assume the viscous part is dominated by pressure solution creep. For simulating such pressure solution we assume (following the logic explained above for halite) a Newtonian viscous flow (e.g., [START_REF] Gratier | How pressure solution creep and fracturing processes interact in the upper crust to make it behave in both a brittle and viscous manner[END_REF][START_REF] Lehner | A model for intergranular pressure solution in open systems[END_REF][START_REF] Shimizu | Kinetics of pressure solution creep in quartz: theoretical considerations[END_REF], with a constant sediment viscosity that is larger than the viscosity of salt (this is required in order to initiate faults, see section 8.4). Different geological conditions predict sediment viscosity with a wide range of values, between ~10 -10 𝑃𝑎 ⋅ 𝑠 [START_REF] Connolly | Temperature-dependent viscoelastic compaction and compartmentalization in sedimentary basins[END_REF][START_REF] Geng | Time and Temperature Dependent Creep in Tournemire Shale[END_REF][START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF][START_REF] Schneider | Mechanical and chemical compaction model for sedimentary basin simulators[END_REF], depending on grain size, strain rate, porosity, pore size, temperature, sediment and water composition and degree of dissolution localization. If deformation localizes, strain rates increase, and so do shear stresses. Once the sediments reach a threshold stress, they will fault. The criteria for sediment faulting includes cohesion and Coulomb friction: cohesion was experimentally measured for marls to be quite low, at C~50kPa and friction angle is measured within a range of 22°-35° [START_REF] Aiban | An investigation of salt tectonic structural styles in the Scotian Basin, offshore Atlantic Canada: 1. Comparison of observations with geometrically simple numerical models[END_REF][START_REF] Anagnostopoulos | Geotechnical properties of the Corinth Canal marls[END_REF][START_REF] Shaqour | Geotechnical and mineralogical characteristics of marl deposits in Jordan[END_REF].

In order to understand the mechanisms by which normal faults form in a brittle-plastic layer coupled to an underlying flowing viscous layer (e.g., [START_REF] Nagel | Channel flow and the development of parallel-dipping normal faults[END_REF], we first develop an analytical model that predicts the conditions (e.g. tilt angle, viscosity values, flow patterns, overburden thickness) that favor the development of normal faults. We then further investigate the effects of these parameters in a more realistic setup using numerical modeling.

Analytical model of faulting above a tilted viscous (salt) layer

A simplified description of the physical setting we explore is a two-layer model on a slope of length 𝐿 and angle 𝛼: a lower viscous layer (of thickness ℎ ) that starts flowing when it is tilted, and an overlying viscoplastic-brittle layer, of thickness ℎ , that is mechanically coupled to it (Figure 3). The two layers are covered by water.

Predicting the onset of faulting

In order to predict normal fault initiation in the brittle layer, we first assess the coupled viscous flow field of the salt and the overburden. Then, we calculate the stresses this flow exerts on the brittle sediment layer and compare them to the failure stress required to initiate a fault.

where 𝑃 and 𝑃 are lithostatic and pore fluid pressure, respectively, 𝛿 is the Kronecker delta, 𝜂 is viscosity of sediments and 𝜀 is the strain rate tensor. Initiation of faulting is assumed to occur for small deformations, where the geometry remains nearly unperturbed, so that the most and the least compressive stresses are 𝜎 and 𝜎 respectively, i.e., stresses perpendicular and parallel to the slope. Therefore, the differential stress, ∆𝜎 , driving faulting is ∆𝜎 = 𝜎′ -𝜎′ . Furthermore, we assume incompressible deformation of the sediment layer: 𝜀 = -𝜀 , so that 𝜎′ = -𝜎′ .

As a result

∆𝜎 = 4𝜂 𝜀 = 4𝜂 𝑑𝑢 (𝜒, 𝜁) 𝑑𝜒 = 4𝜂 𝑛𝜋 𝐿 𝑈 (𝜁) 𝑐𝑜𝑠 𝑐𝑜𝑠 𝑛𝜋𝜒 𝐿 (3)
Eq. ( 3) predicts that the maximum strain rate will be attained at the edges of the model (𝜒 = 0 and 𝜒 = 𝐿). These locations (above salt edges) are thus expected to experience the greatest differential stresses and be the most prone to brittle failure. We recognize that assuming 𝑢 = 0 at the edges of the model forces the strain rate to increase there, but we will justify this assumption by comparison to the full numerical model (Section 6).

The brittle behavior of the sediment is modeled by the Drucker-Prager failure criterion, 𝜎 , (the super script F stands for "failure") assuming hydrostatic pore pressure:

𝜎 = (𝑃 -𝑃 )𝑠𝑖𝑛𝜙 + 𝐶 ⋅ 𝑐𝑜𝑠𝜙 (4)
Where 𝜙 and 𝐶 are the friction angle and cohesion of the sediment, respectively (See Table 1 for a summary of notations). The Drucker-Prager failure criterion should be compared to the square root of the second invariant of the deviatoric stress in sediments, 𝜎 , which is related to the differential stress as:

𝜎 = ∆𝜎 2√2 (5)
This expression is readily obtained from the definition of the second invariant of the deviatoric stress tensor expressed in principal coordinates, under the assumption that the out-of-plane stress is the average of the least and most compressive principal stresses.

The brittle sediment layer will fail if:

𝜎 ≥ 𝜎 (6)
Let 𝑅 be a non-dimensional "failure number", such that:

𝑅 = 𝜎 𝜎 (7)
Therefore, failure occurs as 𝑅 approaches 1.

Analytical results

To identify the primary controls on whether failure occurs, we need to examine 𝜎 and 𝜎 at the bottom of the sediment layer, where the largest deformation takes place due to underlaying salt flow. 𝜎 (Eq. ( 4)) is directly related to the depth below the seafloor (𝑑 = ℎ /𝑐𝑜𝑠(𝛼), since the layer is tilted), so the sediment thickness, ℎ , is expected to be a key control on 𝑅.

The second main control is the ratio of viscosities of the sediment and the salt, , which can be shown to appear in the differential stress term (∆𝜎 ), in 𝜎 (Eq. ( 5)). This is shown from an order-of-magnitude assessment of our model equations: since sediments deform viscously, ∆𝜎 should scale with the strain-rate 𝜀 and the viscosity in the sediments. 𝜀 can be approximated by a velocity difference over a characteristic downslope distance 𝑙 < 𝐿, where 𝑙 is the distance over which flow develops from zero at the salt edge to a fully developed Poiseuille flow. Since the sediment upslope are pinned and do not move, ∆𝜎 may be approximated as:

∆𝜎 ( , ℎ ) = 4𝜂 𝜀 ( , ℎ ) ≈ 8𝜂 ( ( , ) ) , ( 8
)
where we have estimated the stress at the interface between salt and sediment, and at a small down-slope distance, , from the salt edge. The velocity across the boundary between salt and sediment is continuous, i.e., 𝑢 (𝑙, ℎ ) = 𝑢 (𝑙, ℎ ), so, Eq. ( 8) may be written as:

∆𝜎 ( 𝑙 2 , ℎ ) ≈ 8𝜂 (𝑢 (𝑙, ℎ ) -0) 𝑙 (9)
The velocity at the top of the salt scales as a viscous flow down a slope [START_REF] Turcotte | Geodynamics[END_REF]. This assumption simulates flow with a free rigid top (the sediment layer)

that is rafted freely downslope by the fluid layer:

𝑢 (𝑙, ℎ ) ≈ 𝜌 𝑔𝑠𝑖𝑛(𝛼)ℎ 2𝜂 (10)
where 𝜌 𝑔𝑠𝑖𝑛(𝛼) is the driving force for the salt flow downslope. With Eq. ( 10), Eq. ( 9) becomes:

∆𝜎 𝑙 2 , ℎ ≈ 4 𝜂 𝜂 𝜌 𝑔𝑠𝑖𝑛(𝛼)ℎ 𝑙 (11)
Using Eq. ( 4) with 𝑑 = ℎ /𝑐𝑜𝑠(𝛼), together with Eqs. ( 7) and ( 11), the non-dimensional failure number is expected to scale as:

𝑅 ≈ 𝜌 ℎ √2(𝜌 -𝜌 )𝑠𝑖𝑛(𝛷)ℎ 𝑙 𝜂 𝜂 𝑠𝑖𝑛(2𝛼) (12)
Where we have neglected cohesion, since its value is small relative to the shear strength.

Our estimate suggests that although the sediment strain rate depends on the viscosity of the salt, the ratio between sediment and salt viscosities controls ∆𝜎 , and therefore controls 𝑅 (see Appendix B). In the following description of the model results, we therefore emphasize the effect of viscosity ratio, along with sediment thickness and tilt angle. Since 𝑅 is a non-dimensional parameter, in fact the length-scale ratios * , also control the system approach to failure. However, we did not further explore the effect of the ℎ , or 𝑙 in this work, and only explore the effect of ℎ which varies with time as sediments accumulate.

While Eq. ( 12) captures the effect of various salt and sediment properties on failure onset, we provide a more accurate prediction based on the full analytical solution of momentum balance of salt and sediments given in Appendix A. This analytical solution models the full coupling between salt and sediment layers, whilst Eq. ( 12) only assumes that flow is similar to a free fluid layer down a slope. It should be noted that in reality, the failure number (R) cannot exceed 1, because once failure occurs, differential stresses become capped at the failure envelope, as will be shown below in the numerical results. The full analytical calculation predicts that 𝑅 will increase with viscosity ratio (warmer colors at the right-hand side of the diagram) and decrease with overburden thickness (warmer colors towards the bottom of the diagram). The dependence of 𝑅 on viscosity ratio is further demonstrated in Figure B1.

The analytically predicted effect of the viscosity ratio on sediment failure is straightforward:

increasing sediment viscosity increases the differential stress for a given strain rate, while decreasing salt viscosity promotes faster strain rates. For a low viscosity ratio, 𝑅 diminishes slowly with increasing ℎ (barely crossing sub-vertical isolines on left-hand side of Figure 4), while with greater viscosity ratios, 𝑅 diminishes rapidly with increasing ℎ (crossing many sub-horizontal isolines on the left side of Figure 4). The effect of overburden thickness is clear from Eq. ( 4); larger thickness of ℎ strengthen the overburden (𝑑 in Eq. ( 4)) and thus decreases 𝑅 . The 𝑅 = 1 curve in the failure diagram illustrates that each thickness has a critical viscosity ratio that allows faulting, and vice versa, each viscosity ratio has a critical overburden thickness above which faulting is prevented.

Figure 4 also illustrates the effect of tilting. Greater tilt angles increase the sediment strain rate resulting in greater viscous stress for the same failure stress and therefore greater 𝑅 values (Eq. ( 12)). The patterns of 𝑅 in Figure 4 are similar for both angles, with a shift of the 𝑅 = 1 curve that enables failure of a thicker overburden when tilt angle increases. In

summary, the analytical model predicts that fault initiation is promoted by a greater viscosity contrast between the two layers, a greater tilt angle, and a thinner sediment cover.

6. Numerical modeling

Numerical method

To explore a more realistic setting, we numerically study the geometry of a salt layer with a tip that is not pinned at the edge of the domain, plus a slope angle that decreases seaward from a steep margin to a flat basin (Figure 5). Such modelling also allows exploration of failure development beyond the initiation of the first fault. The simulations are carried out using SiStER (Simple Stokes solver with Exotic Rheologies), a 2-D finite-difference/particlein-cell code [START_REF] Gerya | Characteristics-based marker-in-cell method with conservative finite-differences schemes for modeling geological flows with strongly variable transport properties[END_REF], 2007;[START_REF] Harlow | Numerical Calculation of Time-Dependent Viscous Incompressible Flow of Fluid with Free Surface[END_REF] written in MATLAB and based on the methodology of [START_REF] Gerya | Introduction to numerical geodynamic modelling[END_REF]. The code was developed by [START_REF] Olive | The role of elasticity in simulating long-term tectonic extension[END_REF] and described in detail in their paper, while here only the main features of the code are summarized.

Conservation of mass and momentum in a 2-D continuum are solved assuming material incompressibility and steady-state flow (Stokes' equation):

𝜕𝑣 𝜕𝑥 = 0 (13) 𝜕𝜎 𝜕𝑥 - 𝜕𝑃 𝜕𝑥 + 𝜌𝑔 = 0 (14)
where 𝑣 and 𝜎 denote velocities and deviatoric stresses, respectively, and P is pressure (See Table 1 for a summary of notations). Note that we use 𝑣 to denote velocity in a horizontal/vertical coordinate system, (𝑥, 𝑦), as opposed to 𝑢 in a down-slope/slopenormal coordinate system, (𝜒, 𝜁). Repeated indices imply summation. Inertial forces are negligible relative to viscous forces due to the slow rates of deformation considered. Eqs.

(13) and ( 14) are discretized on a fully-staggered Eulerian (non-deforming) grid using a conservative finite difference scheme [START_REF] Duretz | Discretization errors and free surface stabilization in the finite difference and marker-in-cell method for applied geodynamics: A numerical study[END_REF][START_REF] Gerya | Introduction to numerical geodynamic modelling[END_REF]. In what follows we present results of 30 simulations 465 exploring the effect of the tilt angle (𝛼), thickness of sediment cover (ℎ ), and the 466 sediment/salt viscosity ratio (

). The list of simulations and their parameters are 467 presented in Table 2. 468 469 The failure stress for brittle behavior in sediments follows the Drucker-Prager criterion, in which brittle strength increases with pressure 𝑃 according to Eq. ( 4).

Brittle strain, 𝜀 , accumulates where the failure criterion is met. To promote strain localization, cohesion and friction are decreased linearly (from their initial values 𝐶 , 𝜇 )

with increasing plastic strain 𝜀 [START_REF] Olive | The role of elasticity in simulating long-term tectonic extension[END_REF] until they reach a minimal value (𝐶 , 𝜇 ) when a critical strain, 𝜀 , is reached, corresponding to a critical amount of fault slip, ℎ (e.g., [START_REF] Lavier | Factors controlling normal fault offset in an ideal brittle layer[END_REF].

Density of sediments, 𝜌 , is depth dependent:

𝜌 (𝑑) = 𝜌 + 𝑤 ⋅ 𝑑 (15)
Where 𝜌 is the density at the seafloor, and w=0.2 𝑘𝑔/𝑚 . This relation approximates density logs of the Yafo Formation (Plio-Quaternary) in wells (Z. Gvirtzman, unpublished).

Numerical Results

Fault initiation

An example of a simulation that initiates faults by localizing strain within a few grid cells (S14, = 10 , Table 2), versus a simulation that only develops diffuse strain without faults (S11, = 10, Table 2) is shown in Figure 6. Both simulations were run with a tilt angle of 1° and overburden thickness of 600 m. In both simulations, extension (positive 𝜎 ) occurs upslope and contraction (negative 𝜎 ) occurs downslope in the overburden (Note that 𝜎 is deviatoric, not differential stress). When the viscosity ratio is larger (S14, left panels in Figure 6), the strain rate is much larger in both the salt and the sediment (yellow colors in green circles correspond to a transitional case (width within ±50 m of 𝑤 ). The solid line represents the 𝑅 = 1 curve from our analytical model (Eq. ( 7)) for reference (see Figure 4).

We find that the transitional cases approximately coincide with the 𝑅 = 1 curve indicating a good agreement between the analytical and the numerical results. The two independent methods imply that greater tilt angles, lower viscosity ratios, and thinner overburden favor faulting, and vice versa. Further comparison between numerical and analytical results is shown in Figure 8, comparing the "failure number", 𝑅 , from the analytical model with that of the numerical model, 𝑅

, estimated by dividing the square root of the second invariant of the deviatoric stress by the failure stress given in Eq. ( 4). The numerical values of the deviatoric stress and the failure stress used in 𝑅 refer to a location at the center of the shear zone developed above the salt edge, where the strain rate is maximal -see the point of maximum strain rate in Figure 6d. To directly compare the numerical model to the analytical prediction, we turn off brittle failure in the simulations such that 𝑅 can grow beyond 1 and can be compared to 𝑅 in the full range (Figure 8a). This comparison shows a general agreement, which improves for greater viscosity ratios. This demonstrates that our simplified analytical model captures the stresses, which develop in the sediment as a result of salt flow, qualitatively and quantitatively. When failure is allowed, 𝑅 are bound from above by approximately 1 (Figure 8b).

In Figure 8c we show the width of the shear zone, normalized by the grid size, as a function of 𝑅 . The width clearly decreases with 𝑅 until reaching the limit of 1 imposed by the numerical resolution. For the critical value 𝑅 = 1, the width is about 3 grid cells, justifying our choice of the localization criterion introduced above.

Evolution of faulting with time

Unlike the analytical model, the numerical model allows us to track the evolution of brittle strain past the onset of normal faulting. Figure 9 shows the evolution of faulting over a few Myr. In the presence of salt flow, faulting initiates rapidly, in a few kyr (after brittle failure is enabled), with a single normal fault above the salt edge (Figure 9b). Approximately 5000 years later, reverse faulting develops at the bottom of the slope (Figure 9c). Then, the locus of both new normal and reverse faults migrates basinward, in the flow direction (Figure 9d).

For the simulation presented in Figure 9 (S10), extension is accommodated by fewer normal faults compared to the number of reverse faults formed during contraction; thus, on average, the strain rate on a normal fault is higher than that on a reverse fault. Evolution of plastic strain in Figure 9 is complemented by the evolution of the velocity field,

shown in Figure 10. Velocity profiles are also plotted landward of the salt edge (x=2km) and in the sediments below the salt but are practically zero at these locations. It is noticeable that velocity increases downslope at all times, moving away from the "pinned" sediment wall at the salt edge. The velocity difference between the pinned salt edge and the relatively large downdip velocity produces the strain that causes extensional 𝜎′ , as predicted in Section 5. Prior to faulting (F i g u r e 10a) ,Poiseuille flow is observed in the salt "channel" .

After failure )Figure 10b,c (the flow accelerates and evolves towards a mixed Couette -Poiseuille flow )i.e. boundary-driven viscous flow with a free surface at the top) in the salt, and a relatively constant velocity within the overburden that translates downslope as one block. This means that after faulting sediments behave like a "raft" on top of the salt, having detached from their up-slope section through slip on the fault. After a few hundred kyr, the system gradually decelerates, returning to a Poiseuille flow (pressure gradient-driven) like pattern (Figure 10d). At this stage, Poiseuille flow is not limited to the tilted domain as in Figure 10a, but is also visible in horizontal portions of the layers (Figure 10e). Overall, brittle and viscous deformation affect the bathymetry of the basin. Grabens are formed in the extension domain above the salt edge, and small mounds are formed in the contractional domain (Figure 10e). Another way to examine the fault activity is through the velocity of the overburden. Figure 11 shows the horizontal component of velocity, 𝑣 , for S10, at mid depth of the overburden, and three different locations downslope: just above the salt edge (blue curve), and 2 and 5 km down the slope away from the salt edge (red and yellow curves, respectively). Because of the low slope angle, the horizontal velocity component is a good approximation of the downslope velocity.

All three curves exhibit velocity decay after an initial rise. Close to the salt edge (blue and red curves in Figure 11), the velocity shows more variability due to the proximity to faults.

The velocity also decays faster and sooner upslope (red and blue curves) than downslope (yellow curve). An exponential decay fitted to this data gives a decay timescale of 70 kyr for the velocity upslope, and 1.4 Myr for the downslope velocity. It is notable that different salt viscosities will result in different timescales, but here we use a salt viscosity that is generally accepted (10 𝑃𝑎 𝑠 (e.g., [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF]) and the resulting timescale matches observations from nature, as will be discussed in Section 8. Finally, we assess the temporal and spatial development of faults in Figure 12. Strain rate at the mid-depth of the overburden is plotted as a function of the horizontal coordinate and time since instantaneous tilting. For clarity, points are plotted at every timestep but only where brittle strain is larger than 𝜀 (i.e., our proxy for strain localization). Such a plot was constructed for simulations S10, S15, S25 and S30 (Table 2) to examine the significance of overburden thickness and tilt angle.

In all four simulations, the first fault is formed (at small t) at x=5 km, exactly above the salt layer tip (Figure 6d). This fault forms at the transition from no-flow to flow zones. Almost immediately after normal faulting initiates, contraction and reverse faulting occurs at the bottom of the slope (x=25 km), where flow slows down (Figure 10). As time progresses, both normal and reverse faulting propagate basin-ward. After some time, at t~100 kyr for α=2°, new reverse faults also form further upslope (Figure 12b,d). In all simulations, the strain rate is highest when a fault is created and then decays, as is also shown in Figure 11.

Note that low tilt angles (Figure 12a,c) result in less faults relative to high tilt angles (Figure 12b,d). Similarly, thin overburden (Figure 12a,b) leads to more faults relative to thick overburden (Figure 12c,d).

Discussion

This work interprets two aspects of salt tectonics in a tilted salt basin margin with relatively thin overburden: 1) the onset of faulting and 2) the subsequent evolution of the basin and particularly its margin. Analytical and numerical models of viscous salt covered by viscoplastic sediments were used in this study and provided a physical explanation for fault initiation. Viscous salt flow, driven by tectonic tilting, may generate large enough strain rates for stresses to meet the failure criterion in overlying sediments, assuming the displacement rate is continuous across the interface. Numerical simulations support the analytical prediction that a larger viscosity contrast, greater tilt angle and thinner sediment cover promote faulting.

The results show that 1) normal faults initiate in sediments right above the salt edge, in agreement with observations in the Levant margin (Figure 1 and Figure 2) and in other regions [START_REF] Rowan | Salt-related fault families and fault welds in the northern Gulf of Mexico[END_REF]. 2) faults initiate rapidly (~1 kyr) in the presence of a sufficient driving force, here produced by tilted layers, followed by a gradual strain rate decay on them (Figure 11 and Figure 12), in agreement with observations in the Levant basin [START_REF] Elfassi | Chronology of post-Messinian faulting along the Levant continental margin and its implications for salt tectonics[END_REF]. 3) contraction in the deep basin is a result of the tilting and is coupled to the extensional faulting at the margin. The origin of these results is discussed at length below.

8.1. Why do normal faults form at the tilted margins of salt basins?

Normal fault formation at the edge of salt basins is commonly explained as a result of slope instability on steep continental margins. [START_REF] Garfunkel | Large-scale submarine rotational slumps and growth faults in the Eastern Mediterranean[END_REF] argued that the unstable Pliocene-Quaternary section along the Levant margin slides basin-ward over Messinian evaporites that act as a lubricant. [START_REF] Garfunkel | Large-scale submarine rotational slumps and growth faults in the Eastern Mediterranean[END_REF] also suggests that over-pressured shale beds confined between the impervious evaporites promote this process. Similarly, [START_REF] Almagor | Salt-controlled slumping on the Mediterranean slope of central Israel[END_REF] explained that the lubricating evaporitic layers allowed detachment formation followed by faulting and downslope rotation of large overlying fault-blocks. Brun and Fort (2012) similarly attribute faulting in tilted margins to slope instability and estimated the stress within the sediment using landslide mechanics calculations, assuming a detachment surface at the base of the salt (equations 7,8 therein). In these models, the salt has a passive role in that it merely facilitates the sliding of the overburden. A different perspective, attributing an active role to the salt, is that a flowing viscous salt layer exerts shear tractions on the bottom of the overburden and actively modulates its stress state [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF][START_REF] Gradmann | Factors controlling the evolution of the Perdido Fold Belt, northwestern Gulf of Mexico, determined from numerical models[END_REF].

Here we tested the mechanism for how salt causes faulting in the overlying brittle layer, by assigning a viscoplastic rheology to the overburden, along with a viscous rheology to the salt. A visco-plastic rheology for sediments was previously used to model sedimentary basin deformation, and has been shown to capture well the compaction-depth curve of sedimentary basins and other sedimentary features (e.g., [START_REF] Audet | A mathematical model for compaction in sedimentary basins[END_REF][START_REF] Connolly | Temperature-dependent viscoelastic compaction and compartmentalization in sedimentary basins[END_REF][START_REF] Fowler | Pressure solution and viscous compaction in sedimentary basins[END_REF][START_REF] Shi | Pore pressure generation in sedimentary basins: Overloading versus aquathermal[END_REF]. We chose to not add elasticity to our model, since adding elasticity does not affect the stresses that develop after a long time in flowing sediment or salt. This is because a visco-elasto-plastic material, i.e., a Maxwell material, accumulates elastic stress only until reaching the Maxwell characteristic time, 𝜏 , determined by the ratio of viscosity to elastic moduli:

𝜏 = 𝜂 𝐺 (16)
After time 𝜏 , viscous strain-rate and stress dominate, and stresses can be considered constant with time, following Eq. ( 3) (as long as strain-rate is constant and the geometry does not change significantly through the creation of new faults). If the resulting viscous stress within the sediments exceeds the failure strength, faults are produced. In sediments, using 𝐺 = 30 𝐺𝑃𝑎 and 𝜂 = 10 𝑃𝑎 ⋅ 𝑠 (maximal viscosity in our simulations), 𝜏 is calculated to be <10 kyr. Since this rise-time of stress is very short relative to our simulation time scales, stress and strain in our system (following a brief adjustment period) are primarily determined from the strain rate (viscous) and not from the accumulated elastic strain. Thus, a flow-derived stress is the cause of faulting in the brittle layer above the salt tip. The fluidmechanics reason for faulting to occur at the salt tip location is simply explained in Section 5, and briefly summarized next:

Our simulation results emphasize the coupling between the salt and the overburden and their mutual influence. The influence of the salt on the overburden is straightforward: due to velocity continuity between salt and sediment, the flowing salt forces a downslope extensive stress on the bottom of the sediment, which is proportional to the viscosity of the sediment and the strain rate of the salt. The salts' downslope velocity is in turn proportional to the pressure gradient that drives the salt, divided by the salt viscosity. Thus stress (and thus faulting) within the sediments depends on the sediment to salt viscosity ratio, as well as on length scale ratios in the system and the slope tilt angle. An effect that is not included in the above simple argument, yet is accounted for in the analytical and numerical models, is that the flow of the overburden influences the salt. A very viscous overburden that hardly moves (nearly zero flow on the salt-overburden interface), forces a channel-like flow pattern in the salt (Figure 10a). In contrast, a low-viscosity overburden that flows relatively easily (a nearly free salt-overburden interface) allows a close to Couette flow pattern to develop in the salt.

Similarly, a Couette-like flow pattern can be formed even under high viscosity ratios, after faulting when the overburden starts sliding downslope almost freely and forms a nearly free top salt surface (Figure 10b).

In reality, the salt has a wedge shape (Figure 2), but for simplicity we did not use this shape in the model. Our analysis shows that the large stress that leads to sediment faulting near the salt edge does not primarily arise from a stress concentration due to the sharp tip of the salt edge geometrical shape (all the simulations share the same geometrical shape and exhibit stress concentration around the pointed tip). Instead, this large stress arises mostly from the fact that the largest viscous strain occurs near the salt edge, due to the lateral viscosity contrast. The physics causing large strain at the sediment-salt transition (at the salt edge) is robust and will form faults given the right conditions, with only 2 nd order effects of the geometry of the salt's edge. The fact that the faulting potential, 𝑅, depends on salt thickness (Eq. ( 12)), may suggest that thinning of the salt to a wedge shape will decrease failure, or cause it to initiate somewhat downslope, in thicker parts of the salt. However, it is not possible to carry this conclusion, since Eq. ( 12) was derived assuming fully developed 1D flow, and is thus not directly applicable for a wedge type geometry. Yet, using a straight vertical edge, like we have done in this work, is suspected to amplify faulting, and indeed we have done some preliminary tests with a wedge-shaped salt, that resulted in decreased (but still present) faulting. Also, in agreement with the logic outlined above, some seismic profiles from the Levant margin show isolated thin salt pockets upslope of the primary fault, yet further examination regarding the effects of the salt edge geometry is needed in the future.

Finally, we note that the viscosity of the salt may exert some feedback on the geometry and spacing of the normal faults that cut through the sediment. This was suggested by [START_REF] Nagel | Channel flow and the development of parallel-dipping normal faults[END_REF] who studied the most energetically efficient fault geometries that can develop over a viscous layer subjected to far-field extension.

Timing of fault initiation

An important conclusion arising from our simulations is that if faults occur, they initiate rapidly in the presence of gravitational flow. This result directly stems from the viscoplastic rheology: As discussed in subsection 8.1, the charact eristic time f or stress growth is 𝜏 . After time 𝜏 stresses reach their maximum value, and may be approximated as purely viscous.

Only an elastic component of the rheology causes the overburden stresses to accumulate with time, i.e., to depend on strain. This raises the question of how long it would take to initiate faults in a viscoelastic medium. But it is easy to show that the time in which viscoelastic behavior is dominant is short in terms of our simulated system:

If the stress reached at (or before) time 𝜏 is higher than the failure stress, a fault will be created; otherwise, faults will not form. Thus, if our model finds that faulting occurs by time 𝜏 , this timing of faulting may be off by no more than an order of 𝜏 . The viscosity values used for the simulations here, together with 𝐺 = 30 𝐺𝑃𝑎, give 𝜏 <10 kyr, which is relatively short compared to the period of fault activity (Gelasian, 0.8 Myr). Although fault initiation occurs very rapidly after tilting, further faulting occurs at later times, as seen in Figure 12.

The displacement rates measured by [START_REF] Elfassi | Chronology of post-Messinian faulting along the Levant continental margin and its implications for salt tectonics[END_REF] show that faulting peaked in the Gelasian (2.6-1.8 Ma) and decayed since. That peak coincides with the onset of marginal tilting [START_REF] Bar | The uplift history of the Arabian Plateau as inferred from geomorphologic analysis of its northwestern edge[END_REF]. Our model indeed predicts that deformation will peak immediately after tilting. Models show that the driving force wanes as time progresses. Salt topography flattens and the velocity of the sediments and strain rate on all faults decays (Figure 11, Figure 12). This result is also found in analogue experiments (Ge, Warsitzka, et al., 2019) and is further expected when sediments accumulate and the overburden thickens and strengthens (Figure 4). In the Levant margin, the tilt angles are larger in the North [START_REF] Cartwright | Initiation of gravitational collapse of an evaporite basin margin: The Messinian saline giant, Levant Basin, eastern Mediterranean[END_REF], which means that faulting there is expected to intensify (e.g., Figure 4a vs. 4b and Figure 7a vs. 7b) (Ben Zeev & Gvirtzman, 2020).

Locus of fault initiation

The location of the faults is not random; According to Eq. ( 3), the largest stress is expected above the salt edge, where large flow gradients are generated. This is exactly where the first fault is formed in the numerical model (𝑥 = 5 km in the simulations, Figure 9-12). Moreover, the simulations show that the locus of subsequent faulting later migrates by 5-10 km basinward, depending on the tilt angle (Figure 12). This modeling result is consistent with observations from the Levant continental margin, where normal faults are confined to a 10-20 km belt whose eastern boundary approximately coincides with the salt pinch-out line (Figure 1, Figure 2; Bertoni and Cartwright, 2006;[START_REF] Gvirtzman | Intra-Messinian truncation surface in the Levant Basin explained by subaqueous dissolution[END_REF][START_REF] Katz | Spatial and temporal relation of submarine landslides and faults along the Israeli continental slope, eastern Mediterranean[END_REF][START_REF] Zucker | Salt tectonics in the Eastern Mediterranean Sea: Where a giant delta meets a salt giant[END_REF]. The correlation between normal faults and salt edges is well known and is often used as an indicator for identifying salt wedges [START_REF] Rowan | Salt-related fault families and fault welds in the northern Gulf of Mexico[END_REF]. Here we provide a physical explanation and a quantitative model for this phenomenon. The flow in the analytical and numerical solutions in the present study is driven by the elevation gradient at the salt top, but the Stokes equation can be adapted to also include other driving forces and similarly predict the location of onset of faulting in such settings (e.g., [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF] showed that faulting initiates at the shelf edge for a differential loading situation).

In addition to extensional faults at the tilted margin, downslope flow in the salt and sediments causes contraction within the flat basin, simply due to mass conservation. The analytical results (Section 5) predict that contractional strain will be maxima l at the lowest edge of the flow. Thus, thrusts are expected to form where downhill flow is forced to stop, which in our model is where the slope flattens (Figure 10a, Figure 9c). Indeed, in the numerical model contraction occurs where flow stops, i.e., where the salt flattens. In reality, there could be 3D interactions with neighboring salt tectonic systems affecting salt flow. In the Levant basin for example, the location of the down-dip compressional belt along the continental margin is not clear from seismic observations (Ben Zeev & Gvirtzman, 2020) and may have been overprinted by the younger deformation of the circum-Nile belt (Z. [START_REF] Sagy | The enigma of the Jonah high in the middle of the Levant basin and its significance to the history of rifting[END_REF][START_REF] Zucker | Salt tectonics in the Eastern Mediterranean Sea: Where a giant delta meets a salt giant[END_REF], which is not considered in our 2D modeling. The difficulty in identifying the compressional belt is not surprising in light of the modeling results presented above, showing that the compressional belt is wider and milder than the extensional belt (Figure 9, Figure 12). Yet, an indication for this mild compression is illustrated by Ben Zeev and Gvirtzman (2020), who showed that by removing the younger and dominant circum-Nile folds, mild folds with a ~5 km wavelength still remain close to the extension zone (Figure 18 there-in). These mild folds probably express the down-dip shortening that occurred when the updip extension had begun.

In a basin-wide perspective, the extension increases from south to north along the continental slope [START_REF] Cartwright | Initiation of gravitational collapse of an evaporite basin margin: The Messinian saline giant, Levant Basin, eastern Mediterranean[END_REF], which in the terms of this paper implies higher R values in the north. [START_REF] Cartwright | Initiation of gravitational collapse of an evaporite basin margin: The Messinian saline giant, Levant Basin, eastern Mediterranean[END_REF] indeed report higher tilt angles on the north that comply with higher R values. Moreover, south of the Palmahim disturbance (Figure 1), faults are absent. Basin-wide works (e.g., [START_REF] Sagy | The enigma of the Jonah high in the middle of the Levant basin and its significance to the history of rifting[END_REF][START_REF] Steinberg | Origin and nature of the rapid late Tertiary filling of the Levant Basin[END_REF][START_REF] Ben Zeev | When Two Salt Tectonics Systems Meet: Gliding Downslope the Levant Margin and Salt Out-Squeezing From Under the Nile Delta[END_REF][START_REF] Zucker | Salt tectonics in the Eastern Mediterranean Sea: Where a giant delta meets a salt giant[END_REF] show that the salt in this area is covered by a relatively thick sediment layer that complies with lower R values.

Viscosity of the overburden

As discussed above, the viscosity of the sediments plays a major role in controlling the stresses that develop in the sediment layer and the onset of faulting. The value of the viscosity stems from the dominant creep process, which at the expected low temperature and slow strain rates (see fig 6), is likely pressure solution [START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF][START_REF] Rutter | A Discussion on natural strain and geological structure -The kinetics of rock deformation by pressure solution[END_REF].

The effective viscosity of the sediments undergoing pressure-solution is highly variable, and depends on many factors such as sediment and water composition (e.g., [START_REF] Macente | 4D porosity evolution during pressure-solution of NaCl in the presence of phyllosilicates[END_REF][START_REF] Renard | Enhanced pressure solution creep rates induced by clay particles: Experimental evidence in salt aggregates[END_REF][START_REF] Tada | Pressure solution during diagenesis[END_REF], grain size and temperature [START_REF] Rutter | A Discussion on natural strain and geological structure -The kinetics of rock deformation by pressure solution[END_REF], pore size [START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF] , porosity [START_REF] Tada | Pressure solution during diagenesis[END_REF], mechanical damage [START_REF] Ben-Itzhak | Precipitation of CaCO3 in pressure solution experiments: The importance of damage and stress[END_REF][START_REF] Croizé | Experimental calcite dissolution under stress: Evolution of grain contact microstructure during pressure solution creep[END_REF][START_REF] Gratier | The Role of Pressure Solution Creep in the Ductility of the Earth's Upper Crust[END_REF] and degree of pressure solution localization (i.e., whether it is diffuse or localized into cleavage planes or stylolites) [START_REF] Angheluta | Compaction of porous rock by dissolution on discrete stylolites: A one-dimensional model[END_REF][START_REF] Toussaint | Stylolites: A review[END_REF].

A lower bound on the sediment viscosity in the Levant margin can be inferred by applying our analytical model to the dimensions of the salt tectonics system in the Levant. When faulting initiated (~2.6 Ma), the thickness of the overburden was ~400m and salt thickness is estimated to have been 500-1000 m (Figure 2). The margin is tilted by ~1° in central Israel (at the whereabouts of the seismic section in Figure 2). Given these numbers, in order to initiate faults under a 400m overburden (to be on the right-hand-side of Figure 4a and Figure 13), the overburden/salt viscosity ratio needs to be at least 50 (red cross in Figure 13). Assuming salt viscosity of 10 𝑃𝑎 ⋅ 𝑠 (e.g., [START_REF] Gemmer | Salt tectonics driven by differential sediment loading: Stability analysis and finite-element experiments[END_REF], the overburden's minimal effective viscosity is 5 * 10 𝑃𝑎 ⋅ 𝑠. Since faulting is still active today, under 1000m thick overburden (e.g., Figure 2), Figure 13 implies that the sediments viscosity is at least 2 * 10 𝑃𝑎 ⋅ 𝑠. This value is well within the viscosity values inferred for sediments, which range between ~10 -10 𝑃𝑎 ⋅ 𝑠 [START_REF] Connolly | Temperature-dependent viscoelastic compaction and compartmentalization in sedimentary basins[END_REF][START_REF] Geng | Time and Temperature Dependent Creep in Tournemire Shale[END_REF][START_REF] Van Keken | The effective viscosity of rocksalt: implementation of steady-state creep laws in numerical models of salt diapirism[END_REF][START_REF] Schneider | Mechanical and chemical compaction model for sedimentary basin simulators[END_REF]. However, with a viscosity ratio of ~10 , in a system with 𝛼 = 1° and ℎ = 1000𝑚 (solid line in Figure 13), if the overburden continues to thicken beyond the present day 1000 m, faults will cease to grow. A greater tilt angle (dashed line in Figure 13) will allow even thicker overburden to fault, while thinner salt (dotted line in Figure 13) will allow faulting only in thinner overburden.

Conclusions

The main conclusions of this work are:

1) Following tilting of the margin, stress increases in the sediment overburden due to coupled viscous flow in the salt and sediment layers (as opposed to through the build-up of elastic deformation).

2) The primary controls on overburden failure are the thickness of the overburden and the viscosity ratio between the sediments and the salt. 3) Another control on faulting is the extent of the gravitational driving force, which in our case is determined by the thickness of salt and the tilt angle. 4) The largest differential stress in the overburden is predicted to arise above the salt edge. This drives normal fault initiation at that location. Faulting then progresses with new faults forming sequentially basinwards. 5) A contractional belt forms in the deep basin, due to the downslope flow. The first thrust fault forms far from the margin. New thrusts form sequentially basin-wards. 6) The mechanism suggested here quantitatively explains the location of normal faulting and observations of the temporal evolution of brittle deformation in the Levant basin margin. 7) The overburden sediment viscosity is estimated from modeling to be > 2 * 10 𝑃𝑎 • 𝑠 . Given the salt thickness in the levant, Brittle faulting is favored for overburden thicknesses < 1 km.

FigureB1.

Figure 1 :

 1 Figure 1: Salt tectonic extensional faults (black) and contractional folds (red) in the (a) Levant Basin offshore Israel (data from Z. Gvirtzman et al. (2015)). White line represents the location of the seismic line shown in Figure 2. The violet line approximately represents the tip of the salt wedge, where salt thickness abruptly thins below ~150 m (following Bertoni and Cartwright (2006) and Katz et al. (2016)). Note that this line approximately bounds the fault belt from the east except for the abnormal area of the Dor disturbance. (b) extensional faults (black) and contractional folds (red) in the larger view of the eastern Mediterranean (Map coordinates are in UTM 36N). This larger view shows that the deep basin folds are a part of a circum-Nile fold belt (following Zucker et al., (2020)). ESM -Eratosthenes Seamount.

Figure

  Figure 2: (a) Un Figure 1. Red a Pliocene, Q -Q Average seismi i.e., two-way ti

Figure

  Figure 3: Geom (exaggerated in lines. Red circle The flow d boundary shear stre (details in 𝐿, so that: coordinate The analyt general fo

Figure 4

 4 Figure 4 presents the full analytical prediction for 𝑅, using ℎ = 1𝑘𝑚, 𝐿 = 40𝑘𝑚, for (a) 𝛼 = 1° and (b) 𝛼 = 2°. This prediction, 𝑅 , is calculated from Eq. (7) at the salt-overburden interface at the left boundary of the model (location is marked by a red circle in Figure 3) with stress that derives from the analytical solution for the flow (Appendix A). Specifically, ∆𝜎 was calculated from Eq. (3) (with 𝜒 = 0; 𝜁 = ℎ ). 𝑈 was calculated by Eqs. (A23)-(A26) from Appendix A, with 𝜌 = 𝜌 𝑑𝑦, (i.e., the mean density of sediments in depth 𝑑), and 𝜎 is given by Eq. (4) (with 𝜙 = 28°; 𝐶 = 50𝑘𝑃𝑎; 𝑑 = ℎ /𝑐𝑜𝑠(𝛼)). The resulting 𝑅 is presented as a function of sediment overburden thickness vs. sediment/salt viscosity ratio. 𝑅 > 1 predicts faulting while 𝑅 < 1 predicts no faulting. The dots in the figure represent the parameter values for which we ran numerical simulations, as detailed in Section 6.

Figure 4 :

 4 Figure 4: Results of the full analytical solution for failure number, 𝑅 , above the upper edge of the salt, for varying sediment layer thickness and viscosity ratios, plotted as (𝑅 𝑎𝑛 ) for 2 slope tilt angles (a) 𝛼 = 1° and (b) 𝛼 = 2°. Parameter values are given in the main text. The bold isoline is 𝑅 = 1. Black circles mark the cases numerically simulated in this study. The outcome of our numerical runs, i.e., cases that produced localized deformation (faults) vs. cases that did not, are shown in Figure 7.

  Figure 5: sche respectively.

  Figure 6a, This localiz that is form Conversely salt is rath case, even

Figure

  Figure 6: Fault faults (rightstress at times rectangles in a) c) and g) by wh

Figure 7 :

 7 Figure 7: Results of the numerical simulations plotted on two failure diagrams as functions of overburden thickness and viscosity ratio for (a) 𝛼 = 1° and (b) 𝛼 = 2°. Thirty circles represent the thirty simulations of Table 2, color coded according to a characteristic width of the shear zone developed above the salt wedge. 𝑤 here is 150 m. The black line represents the 𝑅 = 1 line from Figure 4. Note that yellow circles (narrow shear zone), representing fault formation, fall to the right of the R=1 curve, whereas blue circles (diffuse shear zone; no faults) fall left of the R=1 curve.

Figure 8 :

 8 Figure 8: Comparison of 𝑅 values predicted by the analytical model, 𝑅 , and (a) the numerical model, 𝑅 , without failure (brittle failure turned off), (b) the numerical model with failure turned on, and (c) characteristic width of the shear zone, above the salt edge, normalized by the grid cell size. Colors refer to viscosity contrast, shape to overburden thickness, and fill to tilt angle. The black vertical line marks R=1. The dashed blue line in (a) and (b) is a 1:1 line.

Figure 9 :

 9 Figure 9: Plastic strain evolution for simulation S10. Normal faults initiate in the steep part of the slope (above the salt edge, not shown). Reverse faults develop where the slope flattens. Both normal and reverse fault systems grow basin-ward with time.

Figure 10 :

 10 Figure 10: Velocity field evolution for simulation S10. Time frames correspond to those in Figure 9. a) just before faulting with Poiseuille-like flow pattern in the salt. b) and c) shortly after fault initiation with Couette-like flow pattern in the salt. d) and e) show gradual decay of the deformation approaching a nearly flat top salt surface. Scale is indicated by a white arrow at the top right corner of the figure.

Figure 11

 11 Figure 11: x-component of velocity of the overburden at three different locations along the slope. Above the salt edge (blue, x=5 km), 2 km (red, x=7km) and 5 km (yellow, x=10km) downslope of the salt edge.

Figure 12 :

 12 Figure 12: Magnitude of the second invariant of the strain rate (𝜀 ) plotted on a time-distance chart for simulations S30 (top left), S25 (top right), S15 (bottom left), and S10 (bottom right). Only large strain rates are plotted to illustrate when and where faults form and how they propagate with time and space. Note that the salt edge is at x=5km. The region of 5<x< 25 km is the slope, where strain is extensional, forming normal faults. The region of x>25kms is the flat region beyond the slope, experiencing contractional strain and reverse faults. Note that thinner sediment layers and larger tilt angles promote denser faulting in both space and time. Each trace shows a single fault that persists in space and time. The shorter traces seen in c&d show faults that initiated at later times. These newer faults form in a basin-ward sequence and in-between older faults. Ext. and Com. indicate the position of extensional and compressional faults, respectively.

Figure 13 :

 13 Figure 13: The Levant region in the faulting phase-diagram. The faulting-no faulting boundary (𝑅 = 1) is plotted for 3 sets of model parameters: Dashed line -𝛼 = 2°; ℎ = 1000𝑚 (same as Figure 4b and Figure 7b), solid line -𝛼 = 1°; ℎ = 1000𝑚 (same as Figure 4a and Figure 7a), dotted line -𝛼 = 1°; ℎ = 500𝑚. The blue lines represent estimates for the bounds of overburden

Figure 1 .

 1 Figure 1.

  Figure 2.

  Figure 3.

  Figure 4.

  Figure 5.

  Figure 12.

Table 1 :

 1 Summary of parameter notations and values.

	Notation Definition		Value
	𝜌	Seawater density	1000 𝑘𝑔 𝑚
	𝜌	Density of sediments at sea floor	1877 𝑘𝑔 𝑚
	𝜌		Sediment density	See Eq. (15)
	𝜌		Salt density	2100 𝑘𝑔 𝑚
		𝑣	Velocity Field with components
			(𝑣 , 𝑣 )
		𝑃	Pressure (mean stress) Field	-(𝜎 + 𝜎 )/2
	𝜎 , 𝜎	Stress and deviatoric stress tensors
	𝜀	Strain rate tensor	1 2	𝜕𝑢 𝜕𝑥	+	𝜕𝑢 𝜕𝑥
	𝜀	Second invariant of strain rate tensor
	𝜂	Viscosity of seawater	𝜂		* 10
	𝜂		Viscosity of sediments	10 -10 𝑃𝑎 ⋅ 𝑠
	𝜂		Viscosity of salt	10 -10 𝑃𝑎 ⋅ 𝑠
		𝑑	Depth from seafloor
	𝜇	, 𝜇	coefficient of sediments maximum and minimum friction	0.53, 0.265
	𝜙	Angle of internal friction of	𝑡𝑎𝑛 𝑡𝑎𝑛 (𝜙) = 𝜇
			sediments
	𝐺	Shear modulus of sediments		30 𝐺𝑃𝑎
	𝐶	, 𝐶	sediments Maximum and minimum cohesion of	50000, 10 𝑃𝑎
		𝑔	Gravitational acceleration	9.8 𝑚 𝑠
	𝜀	plastic strain
	𝜀		Critical plastic strain		0.01
	ℎ		Critical offset

Table 1 .

 1 The model consists of a Newtonian viscous salt layer (yellow) sandwiched between viscoplastic (i.e., viscous and brittle) sediments (green). For simplicity, the salt edge is modeled with a straight edge rather than a salt wedge (see Section 8.1). The topmost layer (blue) is low-viscosity water, which acts as a "sticky layer" to enforce a free interface between sediment and water[START_REF] Crameri | A comparison of numerical surface topography calculations in geodynamic modelling: An evaluation of the "sticky air" method[END_REF] . Boundary conditions are no-slip at the bottom edge of the model and no shear stress (free slip) on the left, right and top edges.

Table 2 :

 2 List of numerical simulations and their parameter values. See Table 1 for the values of the constant parameters in all

	470 471	simulations.
		Simulatio n S1 S2 S3 S4 S5 S6 S7 S8 S9 S10 S11 S12 S13 S14 S15 S16 S17 S18 S19 S20 S21 S22 S23 S24 S25 S26 S27 S28	𝜶 (°) 𝒉 𝒔𝒆𝒅 (𝒎) 1 200 1 200 1 200 1 200 1 200 1 400 1 400 1 400 1 400 1 400 1 600 1 600 1 600 1 600 1 600 2 200 2 200 2 200 2 200 2 200 2 400 2 400 2 400 2 400 2 400 2 600 2 600 2 600	𝜼 𝒔𝒂𝒍𝒕 (𝑷𝒂 𝒔) 10 18 10 18 5 ⋅ 10 𝜼 𝒔𝒆𝒅 (𝑷𝒂 𝒔) 10 19 10 18 10 20 10 16 10 19 10 18 10 22 10 18 10 19 10 18 5 ⋅ 10 10 18 10 20 10 16 10 19 10 18 10 22 10 18 10 19 10 18 5 ⋅ 10 10 18 10 20 10 16 10 19 10 18 10 22 10 18 10 19 10 18 5 ⋅ 10 10 18 10 20 10 16 10 19 10 18 10 22 10 18 10 19 10 18 5 ⋅ 10 10 18 10 20 10 16 10 19 10 18 10 22 10 18 10 19 10 18 5 ⋅ 10 10 18 10 20	Grid size (nodes) 10 43x1051 𝜼 𝒔𝒆𝒅 𝜼 𝒔𝒂𝒍𝒕 50 43x1051 10 43x1051 10 43x1051 10 43x1051 10 47x1051 50 47x1051 10 47x1051 10 47x1051 10 47x1051 10 51x1051 50 51x1051 10 51x1051 10 51x1051 10 51x1051 10 49x1051 50 49x1051 10 49x1051 10 49x1051 10 49x1051 10 53x1051 50 53x1051 10 53x1051 10 53x1051 10 53x1051 10 57x1051 50 57x1051 10 57x1051

Acknowledgments

EA thanks ISF grant #910/17. SP acknowledges the support of Grant No. 19-21114Y from the Czech Science Foundation (GA CR). This research was partially funded by the Israel National Committee for Earthquake Preparedness. It is further based upon work of the European Cooperation in Science and Technology (COST) Action "Uncovering the Mediterranean salt giant" (MEDSALT). The authors would like to thank HIS Markit for providing Kingdom academic licenses for seismic interpretation. We thank Roger Buck and Garrett Ito for insightful discussions. We thank S. Gradmann and an anonymous reviewer for helpful comments and suggestions that much improved this manuscript.

The SiStER code used for this work is available in an online archive (Hamdani et al., 2020).

. 3D seismic analysis of circular evaporite dissolution structures, Eastern Mediterranean. Journal of the Geological Society, 162(6), 909-926. https://doi.org/10.1144/0016-764904-126 Bertoni, C., & Cartwright, J. A. (2006). Controls on the basinwide architecture of late Miocene (Messinian) evaporites on the Levant margin (Eastern Mediterranean). https://doi.org/10.1016/j.sedgeo.2006.03.019 Bertoni, C., & Cartwright, J. A. (2015). Messinian evaporites and fluid flow. Marine and Petroleum Geology, 66, 165-176. https://doi.org/10.1016/j.marpetgeo.2015.02.003 Brun, J.-P., & Fort, X. (2011). Salt tectonics at passive margins: Geology versus models. Marine and Petroleum Geology. https://doi.org/10.1016/j.marpetgeo.2011.03.004 Brun, J.-P., & Fort, X. (2012). Salt tectonics at passive margins: Geology versus models -Reply. Marine and Petroleum Geology, 37(1), 195-208. https://doi.org/10.1016/j.marpetgeo.2012.04.008 Planetary Science Letters, 379, 108-119. https://doi.org/10.1016/j.epsl.2013.07.018 Gvirtzman, Z., Reshef, M., Buch-Leviatan, O., Groves-Gidney, G., Karcz, Z., Makovsky, Y., & Ben-Avraham, Z. (2015). Bathymetry of the Levant basin: interaction of salt-tectonics and surficial mass movements.

Appendix A -Analytical derivation

The tilted margin is simplified to a stratified flow problem. Salt and sediments are layers of thickness ℎ and ℎ , respectively. The layer overlay a slope of length 𝐿 and angle 𝛼 and the slope is under water (see Figure 3).

We chose a coordinate system, in which χ-coordinate runs along the slope and ζ-coordinate measures the distance from the bottom of the salt layer, as shown in Figure 3.

This setup leads to a viscous flow in the χ direction with velocity field, 𝑈, 𝑈(𝜒, 𝜁) = {𝑢 (𝜒, 𝜁) 0 ≤ 𝜁 ≤ ℎ 𝑢 (𝜒, 𝜁) ℎ ≤ 𝜁 ≤ ℎ + ℎ (A1)

The velocity field satisfies the Stokes equation 𝜂𝛻 𝑈 -𝛻𝑃 + 𝜌𝑔𝑠𝑖𝑛𝛼 = 0, (A2) No tractions are assumed on the interface between the salt and the sediments.

Since the slope is submerged to a horizontal water level, weight of water above a point on the slope varies along the slope. Hence, a pressure gradient, 𝛻𝑃, in the χ direction is created

And we get two differential equations: one for salt (0 ≤ 𝜁 ≤ ℎ ) and the other for sediments (ℎ

Where 𝑘 ∈ 𝑠𝑎𝑙𝑡 = 1, 𝑠𝑒𝑑 = 2. To get the velocity field 𝑈(𝜒, 𝜁), solving both equations of 𝑢 is required. Although the two equations are independent, we require the solutions to match at the interface.

Boundary conditions:

We impose no-flow boundary conditions in the 𝜒 direction, 𝑢 (0, 𝜁) = 𝑢 (𝐿, 𝜁) = 0.

Expanding the flow velocity into a Fourier series, this boundary condition implies the following form of the solution 𝑢 (𝜒, 𝜁) = ∑ 𝑈 (𝜁)𝑠𝑖𝑛( ).

(A5) Similarly, the right-hand side of Eq. (A4), 𝛾 , is written as

For simplicity, 𝛾 is treated as a constant, neglecting the variation of density with depth.

Adding a linear variation of density would still be analytically feasible, but since the relative contribution of the linear term to the total density is within 6% for the present setup, we prefer to neglect it and simplify the expressions below.

Using the inverse Fourier transformation, the coefficients of the series 𝛤 may be expressed as

Thus, all even coefficients are zero. Using the expansions (A5) and (A6), Eq. (A4) implies the following identity for the flow velocity coefficients

The solution to the corresponding homogenous equation, 𝑈 ,

) and a particular solution, 𝑈 , satisfying (A9) is

(A12 ) Therefore, a general solution to (A9) has the following form 𝑈 (𝜁) = 𝐴 𝑒 + 𝐵 𝑒 + ( ) ((-1) -1), (A13 ) and the flow velocity becomes:

(A14 )

The coefficients 𝐴 and 𝐵 are to be determined from the boundary conditions in 𝜁 direction.

Assuming no slip at 𝜁 = 0, the velocity of salt at the bottom vanishes ) leading to

Since the top of the sediments is a free surface, shear stress and thereby vanish at

The resulting values (i.e. (A23), (A24), (A25) and (A26)) were put back in Eq. (A14).

The analytical solution was compared to a numerical simulation with 𝐿 = 50 𝑘𝑚, ℎ = 1000 𝑚, ℎ = 1500 𝑚, and 𝛼 = 3°.

The series solution (A14) converges quickly with n, as shown in Figure A1. Addition of more than five modes does not change the solution significantly (Figure A1), so only the first five modes are used in what follows. Clearly, the analytical solution somewhat overestimates the numerical flow velocity, particularly in the sediments. The difference arises from:

1) No vertical velocity in the analytical solution, while it is considerable in the numerical model. 2) No shear traction at the interface between the sediments and water is assumed in the analytical model. The discrepancy between the analytical and numerical flow velocity in the sediments suggests that the traction is significant along the slope.

Despite the deviations, the analytical prediction and the numerical model are qualitatively similar, while the values differ by a factor of up to 4.

Appendix B -Viscosity ratio

𝑅 depends on the ratio of viscosities of the sediments and the salt, while the magnitudes of the viscosities may change. To verify that, we have calculated 𝑅 for a wide range of viscosities of the salt and the sediments. Figure B1 shows 𝑅 as a function of the viscosities of the salt and the sediments. Contour lines of constant 𝑅 (solid lines) increase along the slope of 1, lining up with contour lines of constant viscosity ratio. Hence, 𝑅 linearly increases with with a coefficient independent of the viscosities.