Characterization of Late Cretaceous to Miocene source rocks in the Eastern Mediterranean Sea: An integrated numerical approach of stratigraphic forward modelling and petroleum system modelling

Stratigraphic forward modelling was used to simulate the deposition of Upper Cretaceous, Eocene and Oligo‐Miocene source rocks in the Eastern Mediterranean Sea and, thus, obtain a process‐based 3D prediction of the quantity and quality distribution of organic matter (OM) in the respective intervals.


| INTRODUCTION
The eastern Mediterranean Sea (Figure 1) evolved to one of the world's most important hydrocarbon (HC) frontier exploration areas, after the discoveries of several large gas fields during the last decade. Next to the Nile Delta Cone, the eastern Mediterranean Sea features the Levant and the Herodotus basins, where the recent discoveries are either located in Miocene turbidites (such as the fields of Tamar and Leviathan offshore Israel, and the Aphrodite field offshore Cyprus) or in Mesozoic to Cenozoic shallow water carbonates (Zohr field offshore Egypt or Calypso and Glaucus offshore Cyprus). Both types are reported to be mainly sourced by gas of biogenic origin, although it is speculated that thermogenic petroleum systems exist in the region as well (e.g. Al-Balushi, Neumaier, Fraser, & Jackson, 2016;Bou Daher et al., 2016;Barabasch et al., 2019. In such new offshore frontier exploration areas, scientific as well as industrial challenges relate to high uncertainties and risks, which are associated with the general lack of well and core data. Seismic reflection data, which today are available for large parts of the eastern Mediterranean Sea, can certainly help to identify potential reservoir structures, yet the characteristics of reservoir rocks and sealing units are hardly predictable based on seismic data only. Another crucial element of potential petroleum systems is the distribution, quality and thickness of source rocks. Although there are possibilities to estimate the total organic carbon (TOC) content of a potential source rock based on special seismic reflection characteristics (Eastwell, Hodgson, & Rodriguez, 2018), a precise lateral and vertical prediction of the quantity and quality of sedimentary organic matter (OM) requires measured geochemical data and geological modelling. Therefore, additional information needs to be acquired from adjacent on-and offshore areas. Potential Upper Cretaceous and Upper Palaeocene source rock intervals have been proven along the eastern margin of the Levant Basin (e.g. Almogi-Labin, Bein, & Sass, 1993;Bou Daher, Nader, Müller, & Littke, 2015;Bou Daher, Nader, Strauss, & Littke, 2014) and in the Upper Palaeocene (Bou Daher et al., 2016). At the southern margin, source rocks are present in Jurassic to Miocene intervals in northern Egypt (Shaaban, Lutz, Littke, Bueker, & Odisho, 2006). At the western margin of the basin, Upper Cretaceous and Upper Eocene source rock units occur within the sedimentary cover of the Eratosthenes Seamount (ESM) south of Cyprus (Grohmann et al., 2018). Onshore eastern Cyprus, thermally immature OM-rich intervals have been identified within an Upper Miocene succession (Grohmann, Romero-Sarmiento, Nader, Baudin, & Littke, 2019).
Simulation results from Bou Daher et al. (2016) demonstrate that present-day oil generation only takes place at great depth of more than 4,300 m in the centre of the Levant Basin, within Palaeogene and older sedimentary strata. Therefore, in terms of thermogenic petroleum systems, the basin-ward extent in quantity and quality of OM is the crucial uncertainty, which determines whether these sediments can be considered as actual HC-producing source rocks. In addition, the lateral distribution of OM throughout the basin has to be known in order to achieve meaningful predictions of biogenic HC generation, which largely depends on a narrow temperature window and specific sedimentation rates (Schneider, Dubille, & Montadert, 2016).
One way to address these uncertainties can be provided by stratigraphic forward modelling (e.g. Granjeon & Joseph, 1999), which allows the simulation of sedimentary processes and, hence, to generate a process-based model of spatial and temporal distribution of the sedimentary infill of a basin. Such models have recently been provided for the Miocene (Hawie et al., 2015) and for the whole Mesozoic to Cenozoic geological history of the Levant Basin (Barabasch et al., 2019). However, both models focus on the general sedimentary infill of the basin and do not consider the distribution of OM. Further, both models include only the central part of the basin and its eastern margin, but not the western and southwestern margins, featuring the ESM and the Nile Delta Cone area respectively.
For an accurate petroleum system evaluation in the Eastern Mediterranean, missing information on source rocks distribution and their properties is an essential knowledge gap. Therefore, the scope of this contribution is to present one large-scale numerical model including both the eastern and western margin of the Levant Basin, the latter covering the greater area of the ESM. In addition to the distribution of clastic and carbonate sediments, the production and sedimentation of OM is simulated for Upper Cretaceous to Miocene periods. The resulted models provide a quantitative estimation of the lateral and vertical distribution of source rock quality. They are also incorporated into a classic petroleum system model in order to identify zones of thermal maturation and determine where thermogenic petroleum might occur. of stratigraphic forward modelling and petroleum system modelling on a unique numerical platform allowed a comprehensive understanding of the source rocks distribution as well as their varying capability to generate hydrocarbons. Hence, the approach presented in this contribution helps in de-risking exploration and could be applied on similar frontier hydrocarbon provinces.

| GEOLOGICAL SETTING
The geodynamic evolution of the Eastern Mediterranean region started with the Pangea breakup in the Late Palaeozoic with the disintegration of Gondwana and the opening of the Neotethys Ocean in the north of the Afro-Arabian plate ( Figure 2; Robertson, Parlak, & Ustaömer, 2012   and Gardosh and Druckman (2006) were used to construct the eastern part of the static model. Red-dashed polygon shows the area with new isopach maps of the Eratosthenes Seamount, provided by Papadimitriou et al. (2018). Profiles A-A' and B-B' of the respective petroleum system model are shown in Figure 9d,e. The burial history of locations E1 and E2 (pink dots) is illustrated in Figure 13. (b) Depth horizons providing the base structure for the static model rifting phase is proposed to have started in the Triassic followed by a main rift-drift event during the Early Jurassic (Ghalayini, Daniel, Homberg, Nader, & Comstock, 2014;Segev & Rybakov, 2010) leading eventually to the formation of the Herodotus Basin in the west and the Levant Basin in the east (Figure 1). Whereas the Herodotus Basin is believed to be floored with oceanic crust (Granot, 2016), the basement of the Levant Basin is supposed to consists of thin, stretched continental crust (Inati et al., 2016(Inati et al., , 2018. Both basins are separated by a thick fragment of continental crust, the Eratosthenes Continental Block (ECB). On its top the ESM evolved, which represents a huge carbonate platform (Emeis, Robertson, & Richter, 1996). Time-equivalent shallow water carbonate deposits have also been documented along the Levant Margin .
In the Late Cretaceous, the overall extensional setting changed to a compressional one when the Afro-Arabian plate moved northwards and a north-verging subduction zone evolved along the southern margin of Eurasia Montadert, Nicolaides, Semb, & Lie, 2014). Slab pull-related forces along the southern Eurasian Margin caused the formation of new oceanic crust and the obduction of several ophiolitic complexes ( Figure 2), such as the Troodos Ophiolite, which today represents the basement of the island of Cyprus (Pearce & Robinson, 2010). However, at that time, Cyprus was in a deep marine environment. A general sea level rise in the Late Cretaceous together with post-rift subsidence caused the drowning of the ESM as well as the whole Levant region leading to a more pelagic carbonate sedimentation throughout the area. Further, a large-scale upwelling system affecting the whole Levant region is proposed to be responsible for the deposition of Upper Cretaceous source rocks along the Levant Margin (Almogi-Labin et al., 1993;Bou Daher et al., 2014 and the ESM (Grohmann et al., 2018). Besides, upwelling-like processes above the ESM are believed to have taken place in the Eocene, resulting in the deposition of source rock intervals around the seamount (Grohmann et al., 2018).
In the Oligocene, the breakup of Afro-Arabia with an accelerated northward drift of the Arabian Plate initiated the ultimate isolation of the Mediterranean Basin and caused a strong uplift of the Levant Margin, thereby changing the overall carbonate dominated towards siliciclastic-dominated depositional environments (e.g. Gardosh, Druckman, Buchbinder, & Calvo, 2008;Nader, 2014). Huge amounts of fine and coarse-grained clastic material were shed into the basin-forming stacked turbidite features reaching far into the Levant Basin (Hawie et al., 2015). Mainly Oligocene but also Miocene source rocks are believed to be mainly responsible for the biogenic gas which was discovered during the last decade in the Eastern Mediterranean Sea. Since the Miocene, the ECB is in collision with the Troodos Ophiolite, which caused a strong local uplift and the rise of Cyprus above sea level (Galindo-Zaldivar, Nieto, Robertson, & Woodside, 2001;Robertson, 1998a).

| WORKFLOW
In the present study, DionisosFlow ® and TemisFlow ® software packages (Beicip Franlab) were used for stratigraphic forward and petroleum system modelling respectively. For the 3D, present-day structural representation of the model area, 3D depth maps are used ( Figure 1b). For the eastern Levant Basin and the onshore Levant Margin these maps are derived from seismic interpretation by Gardosh and Druckman (2006) and  and new 3D maps of the greater ESM area are provided by Papadimitriou et al. (2018). The geometrical construction of the northern part of the Nile Cone Delta interpolated from a schematic cross-section (Tassy, Crouzy, Gorini, & Rubino, 2015).

| Stratigraphic forward modelling
Stratigraphic forward modelling is based on the 'source to sink' principle, where the spatial and temporal transport and deposition of sediments can be simulated per distinct time intervals on an initial palaeo bathymetry ( Figure 3). The available net accommodation space for sediment deposition is provided by the combination of sea floor subsidence and the variations in palaeo sea level (as provided by Haq, Hardenbol, & Vail, 1988). Subsidence maps are calculated based on the initial and final bathymetry of a simulated seismic package as well as its thickness. For the present model, a lateral resolution of 5 × 5 km was chosen, whereas the vertical resolution represents time steps of 0.5 Ma. In the model, different sediments are defined by distinct properties, such as grain size, density and their transport behaviour expressed by an empirical diffusion equation, which vary between continental and marine environments, as well as between gravitational or fluvial driven transport: where Q s is the sediment flux (km 2 /a); k s i and k w i are the diffusion coefficients for different processes (km 2 /a); i refers to the respective sediment; h is the elevation (m); Q w is the local water discharge in each cell (dimensionless); S is the slope gradient in each cell; and n and m are the constants, between 1 and 2 (Tucker & Slingerland, 1994). Six different sediment types and sediment components were simulated, which are grouped into three classes: (1) carbonates, (2) clastics and (3) OM. Whereas grain size and density are kept constant, the diffusion coefficient is individually adjusted for each simulation. The respective sediment characteristics are shown in Table 1. Diffusion coefficients throughout the literature range from <0.01 km 2 /ka (Gawthorpe, Hardy, & Ritchie, 2003) to >10 7 km 2 /ka (Burgess, Lammers, van Oosterhout, & Granjeon, 2006). These diffusion coefficients are chosen to represent different transport behaviours of each sediment class relative to each other, rather than representing real physical parameters. Clay for example is easily transported than sand grains and thus has a higher diffusion coefficient. For numerical simplification, each seismic package was simulated individually. Compaction of these packages was simulated based on present burial depths without considering any density differences of the overburden lithological column and thus has to be considered as only a rough approximation.

| Simulating carbonate production
The relative production of carbonate takes place in situ as a function of water depth (Figure 4), whereas the absolute production can vary throughout the geological history and is used as calibration parameter. However, absolute production rates (Table 2) were chosen in a way to account for changing climatic (temperature) conditions during the simulation, so that there is, for example, a significant increase in the production rate towards the Oligo-and mainly Miocene, when the climate in the Mediterranean Region was significantly warmer (Gatt & Glyas, 2012). For the present model, two different types of carbonate were simulated: (a) Carbonate mud, representing mainly pelagic carbonate sedimentation with the highest production rate between 100 and 200 m water depth and an exponential decrease down to 3000 m; and (b) carbonate grains representing shallow water, platform carbonates with the highest production rate between 20 and 40 m water depth and no production below 180 m.

| Simulating siliciclastic sediment input
Siliciclastic sediments enter the model via sediment input sources ( Figure 3b-d), which are placed at the model's border. In the model, only riverine sources are considered, whereas aeolian input of fine clastic material is neglected. Nine sediment sources were defined to simulate the clastic input during the Oligo-Miocene (Figure 3b-d): Sources 1, 2 and 3-sediments originating from the Arabian Plate in the east; source 4, 5, 6 and 7-a broad sediment input originating from an adjacent Nile river(s) in the southwest of the model. Note, that the model area does not include the onshore area in the south. Therefore, the Nilotic sediment sources were set relatively broad, to avoid a punctual, radial sediment distribution around these sources. Further, the source 8 was placed at the northwestern boarder of the model to simulate the sediment input derived from the western Anatolian area, and the source 9 was placed at the northeastern boarder to simulate the sediments originating from the Latakia region since the Upper Miocene. The sediment supply (in km 3 /Ma) is composed of a defined percentage of different sediments and determines the total amount of provided sediment, whereas the fluvial discharge (in m 3 /s) influences the transport distance. Both parameters were chosen according to global drainage relationships ( Figure 5) and expected sediment volumes inferred from isopach maps.

| Simulating production and preservation of organic matter
Two different types of OM were defined in the model: (a) Terrestrial OM, which is treated like any other clastic sediment and enters the model as a fraction of the total sediment supply via the input sources. Hydrogen index (HI) values of terrestrial OM are predefined and do not change throughout the simulation. (b) Marine OM is assumed to be derived mainly from phytoplankton, which is produced in situ according to specific primary production (PP) laws at the water surface and is defined as the mass of produced organic carbon over one square meter per year (gC m −2 a −1 ). In order to respect different possible scenarios of marine OM production, three different 'subclasses' for marine OM were created: (a) 'normal-marine-OM', where the PP is defined as a function of distance to the shoreline; (b) 'upwelling-marine-OM', where an increased PP is attributed to an upwelling zone defined around a specific bathymetry ( Figure 6); and (c) 'defined-marine-OM', where the PP zones are defined manually by productivity maps. The latter one can be used to test different scenarios where increased PP rates are expected in specific areas. The relative PP rate for each subclass can vary for different ages allowing a separated or combined use of all three scenarios. It has to be mentioned here, that upwelling systems are complex systems that are affected by many factors. Wind and ocean current might affect the location and the shape of upwelling areas, especially if these are located above isolated seafloor features, such as the ESM. However, for simplification, no wind-or ocean currents are simulated here, and thus the simulated outcomes represent only the rough area of upwelling and thus the zone of increased primary productivity.
Only a small fraction of the OM that is produced at the water surface is expected to reach deeper environments (Grabowski, Letelier, Laws, & Karl, 2019). Once the marine OM is generated at the water surface, it is transported by the carbon flux (e.g. Honjo, Manganini, Krishfield, & Francois, 2008) towards the sediment water interface (SWI). The exported fraction of PP reaching the SWI (CF) is calculated according to Martin, Knauer, Karl, and Broenkow (1987): with CF being the exported organic flux (gC/m 2 a), PP the primary productivity (gC/m 2 a), z the water depth (in m) and n assumed to be constant and equal to 0.86.
The degradation (decrease of HI) of the marine OM is simulated in relation to the dissolved oxygen content in the water column. To calculate the oxygen concentration, a mass balance equation between oxygen consumption and oxygen renewal is used after Mann and Zweigel (2008). A purely vertical (1D) relationship is assumed for the present model: with Ox z as the oxygen level at the sea floor (no dimensional between 0 = anoxic and 1 = oxic), A z (z) as the mixing rate coefficient, Ox REF the reference oxygen level and D z (z) the oxygen consumption linked to the OM degradation by organisms. In the present model the mixing rate is constant, thus the oxygen level is depending on bathymetry and PP only.
Another parameter controlling the preservation of marine OM in the first few meters in the sediment is the burial efficiency, which depends on the sedimentation rate and the oxygen concentration at the SWI (Betts & Holland, 1991;Tyson, 1995). The respective equation is: with C org as the available organic carbon (in wt.%), k the degradation coefficient (function of the burial depth into the sediment, k = k 0 e −z/z0 ) and S the sedimentation rate (in m/Ma).
By combining these different equations, an initial marine TOC content is calculated during the early diagenesis stage. Note that no thermal degradation of OM is taken into consideration at this point (this is done in the petroleum system model).
The simulated bulk TOC is then the sum of marine TOC plus the terrestrial TOC. In the following paragraphs, the term 'simulated TOC' will be used to describe the modelling results (expressed in %), whereas the term 'TOC' refers to the TOC content measured on actual rock samples (expressed in wt.%).
More detailed descriptions of the above-mentioned processes and their numerical implementation are described in Chauveau, Granjeon, and Huc (2013) and Bruneau, Chauveau, Baudin, and Moretti (2017).

| Model calibration
The calibration of the model was achieved in two steps. First, the thickness of each simulated interval was compared to the actual thickness obtained from the respective isopach maps.
In the second step, the sediment distribution was compared qualitatively to gross-depositional-environment maps from the literature (Barrier & Vrielynck, 2008;. For calibration, the above-mentioned parameters for carbonates and clastics were tested in a reasonable range until one best-fit model was achieved. On this best-fit model, the deposition of OM was simulated.

| Geological framework/static model
The export of the DionisosFlow ® results to TemisFlow ® is one crucial part in this integrated workflow. The DionisosFlow ® stratigraphic forward model is first used to provide a general distribution of the different simulated sediment classes throughout the model area, rather than simulating precise internal sedimentary architectures. For the TemisFlow ® petroleum systems modelling, however, specific facies maps are required, which are correlated with distinct physical properties rather than sediment distribution maps (TemisFlow ® output). Therefore, based on specific threshold values, the different sediment compositions in the DionisosFlow ® model are grouped to facies types, which can be read by TemisFlow ® (using a common facies library-i.e. same lithological characteristics). For the present model, 13 different facies were created (based on standard, predefined TemisFlow ® facies types). Their respective threshold values and physical properties are shown in Table 3.
The second choice that needs to be made addresses the vertical resolution of the TemisFlow ® model. In DionisosFlow ® , 0.5 Ma intervals were simulated. However, regarding the scarce availability of well data, this temporal resolution is too detailed and could not raise a claim of representing reality. Hence, the DionisosFlow ® model was investigated for sequences where significant changes in facies and/or OM content are observed. In TemisFlow ® , the main seismic packages are defined by 3D depth maps obtained from seismic interpretations (Gardosh & Druckman, 2006;Papadimitriou et al., 2018;Tassy et al., 2015). The DionisosFlow ® grouped sequences are then used to split the TemisFlow ® packages according to the proportion of each sequence relative to the total simulated thickness in DionisosFlow ® . The facies of each sequence were determined considering the dominant facies in the  included simulation steps, whereas the TOC and HI values were calculated as the average of each sequence. The selected sequences providing the base for the TemisFlow ® model are shown in Table 4. The Lower Cretaceous consist of four sequences, where the two oldest ones represent the clastic input from the southwest during the Early Cretaceous. The two younger ones represent changes in the palaeo coastline, and therefore a shift from pelagic to shallow water carbonates along the margin-basin transition. The Upper Cretaceous to Eocene interval was subdivided into nine sequences. The four oldest ones represent lateral and vertical changes in the extent of deposited source rock intervals. The following groups reflect lateral changes of pelagic and shallow water carbonates along the Levant Margin coastline. Further, one extra group was chosen to represent the Upper Eocene source rocks. For the Oligocene, Lower-and Upper Miocene, 6, 5 and 13 groups were defined respectively. The intercalation of more fine-grained and more coarse-grained intervals in the basin centre as well as along the margin was maintained to represent reservoir heterogeneities as well as intercalations of source and reservoir rocks.
The two youngest intervals (Messinian 6-5.3 Ma and Pliocene to Quaternary 5.3-0 Ma) were drawn manually, since no detailed information on the post-salt sequences is required for the purpose of this work. For the Messinian, pure salt was assigned to the whole basin. To the margin, sandstones were assigned, and on top of the ESM and the Latakia Ridge, limestone was favoured. For the Quaternary, a mix of 30% sand and 70% shale is assigned to the onshore Levant Margin and the Island of Cyprus, whereas the transition to the basin composes a mix of platform and pelagic carbonates. The ESM is composed of pelagic carbonate whereas the basin and the Nile Delta area are composed of sand, which is mixed with clay and mud towards the deeper basin.
Based on this static model, 3D simulations were performed considering compaction, pressure, thermal maturation and kerogen transformation. HC migration was not taken into account in this model. Classical backstripping techniques were used to reconstruct the basin structural evolution. No faults are included in the present model.

| Geochemical properties/source rocks
Since the exact onset of thermal cracking of kerogen and, thus, the beginning of HC generation depends on the type of OM (Tissot & Welte, 1984), the correct assignment of source rock-specific kinetic parameters is of great importance in petroleum system modelling. For the present model, non-compositional kinetic parameters were adopted from Bou Daher et al. (2016) for the Upper Cretaceous (mainly Type II/IIS kerogen) and Eocene/Palaeocene (Type II/III kerogen) source rock intervals. The latter one was also used for Oligocene-Miocene source rocks, as these source rocks show also a mixed Type II and Type III kerogen composition. Initial TOC and HI values were taken from the respective intervals in the stratigraphic forward model.

| Boundary conditions
In the TemisFlow ® software, the upper thermal boundary is defined by palaeo surface temperature maps. For the present model, respective temperatures were calculated according to Wygrala (1989) based on the palaeo latitude and surface temperatures of the Arabian plate and respective palaeo water depths (taken from the DionisosFlow ® model) for the offshore areas.
For the lower thermal boundary, the basement, the lower crust and the upper mantel are simulated with a heat source of 1333°C at the upper mantle (Allen & Allen, 2013). Thickness values for these intervals are adopted from the model by Bou Daher et al. (2016) and slightly modified according to Inati et al. (2018).

| Model calibration
Model calibration is mainly based on the near-shore well, and pseudo well data published in Bou Daher et al. (2016) for the Levant Margin. Note that these used pseudo wells are based on an independent study, where calibration with real onshore data was achieved and independent crustal-tectonic modelling was used to deduce the evolution of heat flow (Halstenberg, 2014). For the area of the ESM, borehole temperatures are taken from the initial ODP Leg 160 reports, precisely from well 967E (Emeis et al., 1996). Vitrinite reflectance values from the ODP wells (Grohmann et al., 2018) were also used for calibration in this part.

| Stratigraphic forward modelling
The stratigraphic forward modelling results are illustrated in Figures 7-10. They are representative for the carbonatedominated (Late Cretaceous to Eocene) and the clastic-dominated (Oligo-Miocene) depositional environments. Figure 11 shows the thickness error maps for the Late Cretaceous to Eocene, the Oligocene, the Lower and the Upper Miocene.  Figure 7a). In order to obtain any sedimentation on top of the ESM, an almost flat initial bathymetry had to be used, which is different compared to its present-day anticlinal shape. In addition, the palaeo diameter of the platform had to be larger than today. Both observations are in accordance with the tectonic history of the ESM, which underwent significant folding and faulting of its flanks since the Cretaceous and mainly during the Miocene collision of the ESM and the Cyprus Arc . It is Important to note that the ESM is characterized by relatively steep flanks, resulting in erosion and resedimentation of platform carbonates creating a ring of a mixed carbonate facies and mass-transport deposits around the seamount (Figure 7a). In terms of simulated thickness (Figure 11a,b), good results are achieved in the basin, and along the margin (100%-140% of actual thickness). In contrast, locally on top of the ESM and at the southern basin-margin transition, too large thickness values are obtained (140 to >300% of actual thickness based on seismic data interpretation). On top of the ESM these spots with a too large thickness are related to areas where the actual thickness of the isopach maps is very low («10 m), which explains these high percentage values. In general, these patchy thickness errors, especially on top of the ESM emphasize the importance of a well-defined initial palaeo bathymetry for such specific seafloor features, which can differ a lot compared to the present-day state.

| Upper Cretaceous to Upper Eocene
For OM simulation, three different types of marine OM production were used for the Upper Cretaceous to Eocene interval. From 89 to 66 Ma, the 'upwelling-marine-OM' was used in order to simulate the upwelling system that was active during the Upper Cretaceous (e.g. Almogi-Labin et al., 1993;Bou Daher et al., 2015). The PP was defined with the highest rate occurring at a bathymetry of 200 m water depth with DZ of ±50 m and a maximum productivity of 600 g C m −2 a −1 (Figure 6), which is in accordance with primary productivity rates that are observed in recent continental upwelling areas (Lalli & Parsons, 1997). In this way, the deposition of Upper Cretaceous source rocks could be simulated along the Levant Margin as well as on top of the ESM (Figure 7e). In the simulation, PP values reach the maximum of 600 g C m −2 a −1 close to the margin, whereas PP rates above the ESM are lower with 100-500 g C m −2 a −1 (Figure 7c). These production rates result in the formation of an oxygen minimum zone between 100 and 500 m water depth.
Along the Levant Margin, simulated TOC contents average around 1% to 2% in the northernmost part, around 3% to 4% in the central part and around 11% to 12% in the southern part respectively (Figure 7e). These values are in good agreement with the measured organic carbon contents of the immature source rocks in northern ( Figure 12a) and southern Lebanon (1.5 wt% TOC and 3.8 wt% TOC on average respectively; Bou Daher et al., 2014Daher et al., , 2015 and in Israel (6.8 wt% to 14 wt% TOC, Almogi-Labin et al., 1993). Simulated HI values (Figure 7f) average around 300-400 mg HC/g TOC in Lebanon and around 500 mg HC/g TOC in Israel. The highest simulated TOC contents and HI values are found at water depths between 150 and 350 m and are thus located within the centre of the oxygen minimum zone. Below the oxygen minimum zone (max. 600 m water depth) only verylow simulated TOC and HI values can be found. These results are comparable to Quaternary organic-rich sediments related to a recent upwelling system offshore Peru, where highest TOC contents (4.5 wt.% on average) and HI values (418 mg HC/g TOC) are observed at a water depth of about 460 m (Ten Haven, Littke, Rullkötter, Stein, & Welte, 1990). Along the Levant Margin, these simulated intervals span a thickness of 200-300 m, whereas the deeper, basin-ward extents (which show much lower simulated TOC contents and HI values) reach up to 1000 m, without considering compaction at present-day burial (Figure 10d,e). Simulated TOC contents decrease within the Upper Cretaceous from bottom to top (Figure 10e), but are still higher and wider distributed than in the Eocene (Figure 10d).
On top of the ESM, OM is preserved along the outer rim of the seamount, with simulated TOC contents mainly between 0.5% and 3%, which is in agreement with measured TOC contents for the Upper Cretaceous interval in ODP Leg 160 well | 859 967 (0.9 wt% TOC on average, Grohmann et al., 2018). At the location of well 967, simulated organic-rich intervals show similar, irregular spacing with comparable TOC contents, although the exact location does not fit as well as at the margin (Figure 12b), pointing out again the importance of the initial palaeo bathymetry and its influence on the depositional At all locations, a sharp decrease in simulated TOC and HI towards the deeper basin can be observed (Figure 7c,e,f). Figure 7d indicates that the oxygen concentration is almost zero directly below the main PP zone, and much higher in the bottom water of the deeper basin. In addition, Figure 7b and e shows that locations, where OM is accumulated, are characterized by low-sedimentation rates. In turn, for instance in the northern part of the basin/slope, some OM was deposited as well, but sedimentation rates are very high, so that the OM is diluted-resulting in low TOC contents. F I G U R E 9 (a) Long-(green) and short-term (blue) sea level curves after Haq et al. (1988) and the composite curve (red), which is used for the present model. High stand (HS) and low stand (LS) periods are indicated, which are used to define the relative clastic sediment composition during the Oligo-Miocene period. LS periods bring more coarse-grained material into the basin (b), while HS periods are defined by more fine grained deposits (c). (d and e) Shows the resulting intercalation of fine and coarse-grained intervals in the basin at present-day burial. For location of the profiles, see Figure 1 | 861 EAGE GROHMANN et Al.

| Oligocene to Miocene
The simulation of several sediment input sources (Arabia, Nile, Western Anatolia, Latakia, see 3.1.2) illustrates the change from carbonate towards clastic-dominated sedimentation along the Levant Margin and within the basin (Figure 8a). Small and patchy accumulations of a mix of pelagic and shallow water carbonates along the Levant Margin remain only in local spots, where no clastic sediments reach the basin (Figure 8a). The ESM and the area of Cyprus remain carbonate-dominated with mainly pelagic carbonate sediments, mixed with some shallow water carbonates on top of the ESM (Figure 8a).
With the stepwise establishment of the Levant Fracture system from south to north (e.g. Bartov, Steinitz, Eyal, & Eyal, 1980;Ghalayini et al., 2014), sediments coming from the east were blocked and did not further feed the basin. Hence, the sediment input sources 1-3 ( Figure 3) were successively deactivated from south to north during the Oligocene to Miocene: 34-22.5 Ma-sources 1, 2 and 3 active; 19-7 Ma-sources 1 and 2 active, 7-6 Ma-source 1 active. This approach could produce a shift of the depocentre from south to north from the Oligocene, over the Lower-, to the Upper Miocene, which is in accordance with the isopach maps of the respective intervals. Consequently, the sediment supply coming from the Nile had to be increased with each deactivated source in order to still fill the complete basin. Sediment supply as well as fluvial discharge were mainly adopted from Hawie et al. (2015) and Barabasch et al., (2019).
The sediment supply delivered by either fluvial systems or by the erosion of an uprising margin is not constant over geological periods and depends on several factors, such as climatic conditions (precipitation), tectonic activities or sea level fluctuations (e.g. Einsele, Chough, & Shiki, 1996;Mulder & Syvitski, 1995). Changes in these conditions can affect the total amount of transported sediment, as well as the composition of the sediment load. Therefore, such changes can affect the frequency of occurring mass flow events, such as turbidites. For a rough consideration of this relationship, total sediment supply and fluvial discharge were adjusted by +10% during sea level low stands, and by −10% during sea level high stands (Figure 9a). The sediment composition changes from 40% sand, 59.5% clay and 0.5% terrestrial OM during high F I G U R E 1 1 Absolut and relative thickness errors maps for the Upper Cretaceous to Eocene (a and b), the Oligocene (c and d), the Lower Miocene (e and f) and the Upper Miocene (g and h) intervals. For the relative thickness error '0' means 0% of the actual thickness (from isopach maps) was reproduced by the numerical simulation. This means that the model resulted in no sediment deposition at this location. '1' means perfect match of simulated and actual thickness values. '2' means simulated thickness value is two times larger than the actual thickness value at this location | 863 EAGE GROHMANN et Al. stands to 80% sand, 19.5% clay and 0.5% terrestrial OM during low stands respectively. In that way, stacked sand bodies that are intercalated with more clay-rich units could be produced (Figure 9b-e), which is in accordance with seismic interpretation of the offshore area (e.g. . In terms of simulated thickness, good results are achieved within the basin, with 90%-110% of the actual thickness (Figure 11c-h). Best results for the deeper basin part are obtained for the Upper Miocene interval (Figure 11g,h). For this interval, representing the period with the strongest uplift of Mount Lebanon (e.g. , it is important to mention that a relatively thick substratum is required, which can be eroded during the uplift and fill the accommodation space towards the northwest. Sediment supply only from the sediment input sources could not fill this gap. Nevertheless, one area with a large thickness deficit of 200-600 m occurs at the northern basin margin along the Latakia Ridge in the Oligocene as well as in the whole Miocene (Figure 11c,e,g), which represents 20%-80% of the actual thickness in the Oligocene and Upper Miocene respectively (Figure 11d,f,h). This is reasonable, since DionisosFlow ® does not consider any structural deformation, which, however, is to be expected along the ridge since the Oligocene (Hall, Calon, Aksu, & Meade, 2005).
The area of the ESM shows moderately good results. In the Oligocene, it shows a deficit of −50 m in the northwestern and a too large thickness of +50 m in the southeastern part ( Figure 11c), which are between 60% and 200% of the actual thickness (Figure 11d). In addition, there are some spots with high-relative additional thicknesses along the Levant Margin in the Oligocene (Figure 11d), which can be explained by a too high-carbonate build-up in these areas. Running the simulation with different carbonate production rates has shown that these build-ups do not significantly affect the overall sediment distribution within the Basin. They do, however, affect sediment structures within the surrounding clastic deposits, which would affect, for example, reservoir bodies and properties in these intervals. These reservoir properties are not the focus of this work and are, thus, not considered problematic in this workflow. In the Lower and Upper Miocene, the ESM shows higher thicknesses with 200 to partly 300% of the actual thicknesses (Figure 11f,h), which, however, represents an additional average thickness value of less than 100 m (Figure 11e,g).
Oligocene and Miocene source rock samples from the Nile Delta Cone show that OM in these samples is either composed of kerogen Type III or a mix of Type II and III (Shaaban et al., 2006). Hence, a combination of 'specific-marine-OM' as well as the input of terrestrial OM via the sediment input sources (0.5% of total sediment supply) was chosen. Since the production of marine OM is usually controlled by nutrients derived from clastic input, the PP of 'specific-marine-OM' was related to the sediment distribution from the input sources, with PP decreasing from 200 g C m −2 a −1 at the coast towards 50 g C m −2 a −1 in the offshore area (Figure 8c). The input of terrestrial OM as 0.5% of the total sediment supply could produce good results in comparison to the literature values. The simulated TOC contents are generally «1% throughout the basin (Figures 8 and  10a-c). Whereas the initial PP of marine OM is relatively high, clastic sedimentation rates within the basin are also high and dilute the preserved marine OM. Furthermore, the vertical distribution of terrestrial TOC is rather homogenous, showing only minor differences between more clay or more sand-rich intervals (Figure 10a-c), which is in accordance with observations within the Nile Delta, where thick (up to >1000 m) Oligocene intervals of sand-and siltstones with relatively homogenous TOC contents are observed (Rasoul & Khaled, 2019).

| Petroleum system model
The main focus of the integrated petroleum system model was to investigate if the simulated source rock intervals reach far enough into the basin to be affected by thermal maturation, and hence, if they might contribute to viable thermogenic petroleum systems in the Levant Basin. It has to be emphasized here, that the simulation of source rock maturation mainly depends on the maximum temperatures reached during the burial history of a source rock. Next to the spatial source rock distribution, the temperature simulation represents one of the major uncertainties in the presented model. Due to the limited amount of available calibration data for both the present day as well as palaeo temperatures, the results presented in the following paragraphs have to be updated as soon as new calibration data become available for the offshore area. Further has to be noted here, that only Easy%Ro  (Emeis et al., 1996) are used. Note that the grey temperature dots (B) represent most likely inaccurate temperature measurements (see Emeis et al., 1996). (b) Temperature calibration, (c) vitrinite reflectance calibration. Depth given in metres below sea level. In (B) and (C) depth values represent metres of burial, relative to the sediment surface  (Sweeney & Burnham, 1990) was used for simulating the vitrinite reflectance which has been the standard calculation scheme for almost three decades. Other recently developed algorithms, such as Basin%Ro (Nielsen, Clausen, & McGregor, 2017) or Easy%RoDL (Burnham, Peters, & Schenk, 2016), will lead to different maturity ranges.
Recently, all three algorithms were compared in the context of petroleum generation kinetics for coal/type III kerogen (Froidl, Zieger, Mahlstedt, & Littke, 2020). Figure 13 shows the present-day simulated temperature throughout the model area (Figure 13a) as well as the temperature (Figure 13b) and vitrinite (Figure 13c) calibration, whereas Figure 14 shows the burial history for two representative locations, one offshore northern Lebanon (Figure 14c,d) and one southeast of the ESM (Figure 14a,b). Temperatures at the top basement are around 120°C along the onshore Levant area. In the basin, the temperature is around 200°C in the southern part and reaches up to >300°C in the northern Basin, where the deepest burial of about 14 km is reached. At the ESM, instead, temperatures at the top basement are very low with values < 100°C. The thermal evolution and the resulting simulated vitrinite reflectance during the burial history are shown in Figure 12 for two respective sites in the basin (one at the eastern flank of the ESM and another one offshore northern Lebanon). values were compared to those obtained from ODP Leg 160 wells 966 and 967. It is believed that the ESM was uplifted by about 1 km during the Miocene (Robertson, 1998b). However, vitrinite reflectance measurements of the Upper Cretaceous source rocks that are present along the ESM (Grohmann et al., 2018) show that the maximum burial depth cannot have been much deeper than today. Hence, no uplift/erosional event was considered in this model. Overall, good fits were obtained for all calibration points in the model (Figure 13b,c). Figure 15 illustrates the lower (Figure 15a,c,e) and upper (Figure 15b,d,f) boundary of the Upper Cretaceous source rock interval at present-day burial. Along the Levant Margin, source rocks with high-simulated TOC contents (»1%) form a narrow, continuous band of about 50 km parallel to the margin. At the ESM, highest simulated TOC contents form a patchier, halo-like distribution around the summit of the seamount. At both locations the source rocks reach farther into the basin at the lower boundary of the interval compared to the upper boundary (Figure 15a,b).

| Upper Cretaceous source rocks
In terms of temperature, Figure 15c and d shows that the highest values are reached in the deepest parts of the basin, which is located offshore northern Lebanon. Temperatures reach 250°C at the lower boundary and around 210°C at the upper boundary. Along the margin, temperatures are much lower with 70°C on average whereas the ESM is even colder with <50°C. Comparing the resulting Easy R o values (Figure 15e,f) with the simulated TOC distribution (Figure 15a,b) shows that the high-quality parts (>1% TOC) of the source rocks are not affected by thermal maturation or are, if at all, early mature ( Figure 15e,f, location 2 and 3). Nevertheless, the basin ward, low-TOC (<1%) source rock facies reaches higher maturities. These zones are located ~50 km offshore northern Lebanon, where both the lower and upper boundary reach into the late dry gas window (Figure 15e,f, location 1), ~100 km east of the ESM, where early dry gas maturity is reached (Figure 15e,f, location 5 and 6), and ~150 km south of the ESM, where the outermost parts of the source rocks reach into late wet gas window (Figure 15e,f, location 4). Solid bitumen can be found in Upper Cretaceous intervals onshore northern Lebanon (Bou Daher et al., 2015), as well as in Upper Cretaceous and Eocene intervals of the ESM (Grohmann et al., 2018), supporting the model outcome of deeper located HC kitchen areas adjacent to these two locations. Figure 14 shows the burial history for two locations, one offshore northern Lebanon and one southeast of the ESM. In addition, Figure 15g,h shows the evolution of the kerogen transformation ratio through time for six-specific locations (see above). The high-simulated TOC source rocks (>1% TOC) along the margin, which are only slightly affected by thermal maturation, started to generate HC 20 Ma ago and have converted less than 20% of its initial kerogen up to now (Figure 15g,h, locations 2 and 3).
Offshore northern Lebanon, petroleum generation started at the lower source rock boundary in the late Palaeocene/ early Eocene (ca. 60 Ma) and during the middle Eocene at the upper boundary (Figure 15g,h, location 1). This is ca. 15-20 Ma earlier than proposed by Barabasch et al. (2019) or Bou Daher et al. (2016, who, however, did not distinguish between the lower and upper boundary of these source rocks. Source rocks falling into this zone of maturation have converted 50% of its initial kerogen between 20 and 15 Ma, and almost 100% at present day (Figure 15g,h, location 1).
Approximately between 50 and 100 km east of the ESM (Figure 15g,h, location 5) petroleum generation started at the lower source rock boundary in the early Eocene around 50 Ma ago and before 40 Ma at the upper boundary. They have converted 50% of its kerogen between 18 and 15 Ma ago and reach a transformation ratio of 1 and 0.8 at present day. South of the ESM (Figure 15g,h, location 4), temperatures are lower, resulting in a slightly later initiation of petroleum generation and a transformation ratio of 0.6 and 0.4 at the lower and upper boundary respectively.
Based on the distribution of TOC and the area of maturation, migration pathways in the eastern part of the basin are (without considering any faults) expected towards the direction of the Levant Margin, as well as towards the Latakia Ridge. With respect to potential reservoirs offshore Lebanon, anticlinal structures of the Syrian Arc might be considered as potential targets. These folds formed during late Turonian to early Eocene time (Eyal, 2011) and were, thus, already present before the main phase of hydrocarbon generation.
In the area around the ESM, migration should be either vertical or, depending on the availability of potential impermeable layers, towards the seamount. This interpretation fits well to the presence of solid bitumen, which was observed in the Upper Cretaceous intervals in the ODP Leg 160 wells (Grohmann et al., 2018).

| Eocene source rocks
For the Eocene interval ( Figure 16) temperatures along the margin and the ESM are almost the same as in the Upper Cretaceous interval. In the northern basin, temperatures reach values of 150-170°C. There is no significant difference between the lower and the upper boundary.
Along the margin, the main parts of the source rocks are immature (Figure 16a,c, location 1). At the northern margin-basin transition only rocks with very low-simulated TOC contents (< 0.3%), which are not considered as thermogenic source rocks (Hunt, 1967), reach into the zone of oil and wet gas generation (Figure 16a,c). For the source rocks along the ESM, most zones with simulated TOC contents >0.3% are immature. Only the southernand eastern-most extents reach into the onset of the oil-and early wet gas window respectively (Figure 16a,c, locations 2, 3 and 4). At these locations, petroleum generation started ca. 25 Ma ago. Up to today, they have converted less than 20% of their kerogen. As seen in Figures 14 and 16d, enhanced thermal maturation was reached after the Messinian salt was deposited, increasing the chance of HC trapping under the salt. Potential reservoir structures might be present as structural traps that are related to the formation of the Levant Fracture System during the late Miocene.

| Oligo-Miocene source rocks
Present-day temperatures in the Oligocene to Miocene intervals are clearly lower than in the underlying intervals. Temperatures in the basin are around 100°C (Figure 17d) reaching up to 150°C (Figure 17c) in the deepest, northern part in the lowest Oligocene source rock interval, around 100°C in the lowermost Miocene (Figure 18c), and 50-70°C in the uppermost Miocene source rock interval ( Figure 18d).
Generally, it has to be mentioned that sedimentation rates during the Oligo-Miocene were high, which can have significant effects on temperature equilibrium in sediments, especially if sedimentation is non-uniform or when overpressures are generated leading to pulsed fluid flow. Most pronounced transient effects occur in subduction zones, where sediments are brought to great depth 10-100 times faster than in fast subsiding basins. This leads to significant temperature disequilibrium (Lutz, Littke, Gerling, & Bönnemann, 2004). A detailed treatment of such effects was, however, beyond the scope of this study due to the limited calibration data.
Oligocene as well as Miocene source rocks occur almost throughout the entire basin, but show mainly very low-simulated TOC contents «1% (Figure 17a,b; Figure 18a,b). Oligocene source rocks lie within the oil window in the southern model area with their lowermost boundary reaching into the onset of the wet gas window in the northern part of the basin (Figure 17e,f). Miocene source rocks, instead, just entered the onset of the oil window at their lowermost boundary in the northern part of the basin (Figure 18e,f). Oligocene source rocks started petroleum generation around 15 Ma ago (Figure 17g,h). In the deepest part of the basin, only 10%-20% of the initial kerogen have been converted before the deposition of the Messinian salt ( Figure 17g,h). Thus, further HC generation, which might be still ongoing today, is likely with a good chance of trapping under the salt. However, as these intervals are also believed to be the source for biogenic gas generation (Schneider  (Figures 1 and 13)) at the lower and upper boundary respectively. The cumulative thickness of the source rock at the respective six locations is as follows: 1 = 1,300 m, 2 = 1,100 m, 3 = 1,100 m, 4 = 100 m [Corrections updated on 12 November 2020, after initial online publication: Figure 17 has been corrected in this version.] 870 | EAGE GROHMANN et Al. et al., 2016), it is uncertain how much of the initial TOC is left for thermogenic generation. Miocene source rocks (Figure 18g,h) instead have generated almost no thermogenic HC until today, and are, therefore, primarily interesting as biogenic source rocks.

| CONCLUSION
The present stratigraphic forward model provides for the first time a process-based 3D prediction of Upper Cretaceous, Eocene and Oligo-Miocene source rocks in the Eastern Mediterranean region (Levant Margin and ESM). Simulated TOC contents and HI values are in good agreement with published geochemical data of the respective, immature source rocks. Upper Cretaceous and Upper Eocene source rock models support the idea of upwelling-induced enhanced primary bioproductivity, consequent oxygen depletion and deposition and preservation of OM along the Levant Margin, as well as on top of the ESM. Source rocks along the margin form a narrow band of ~50 km parallel to the palaeo shelf break. Within this narrow band, simulated TOC contents and HI values are high with 1%-11%, and 300-500 mg HC/g TOC, respectively, and decrease rapidly towards the deeper basin. At the ESM, OM accumulations show a more patchy distribution and are restricted to the top of the seamount and its closer vicinity. Simulated TOC contents are between 0.5% and 3% and HI values between 150 and 250 mg HC/g TOC, and show a rapid decrease towards the deeper basin.
Inferred Oligo-Miocene source rocks show that marine OM is diluted by high-sedimentation rates in most locations of the basin. The highest contribution to the simulated TOC content, which occurs relatively homogenously throughout the entire basin but is mostly «1%, is made up by terrestrial OM which is brought into the basin via the sediment input sources. This leads to low quality, but thick kerogen Type III source rocks with a cumulative thickness of 1-2 km.
The incorporation of these simulated source rocks into a classical burial history model shows that some parts of these source rocks are affected by thermal maturation and, hence, suggests the presence of several thermogenic petroleum systems in the Levant Basin: • Upper Cretaceous source rocks started petroleum generation in the late Palaeocene/early Eocene with peak HC generation between 20 and 15 Ma ~50 km offshore northern Lebanon. Southeast of the ESM, such source rocks started petroleum generation in the early Eocene with peak generation between 18 and 15 Ma. At both locations these source rocks show a transformation ratio of almost 1. • Eocene source rocks started HC generation ca. 25 Ma ago between 50 and 100 km southeast of the ESM and reach into the oil to wet gas window at present day. However, until today they have converted less than 20% of their initial kerogen. • Although the Miocene source rocks are mostly immature, Oligocene source rocks lie within the oil window in the southern Levant Basin and reach into the onset of the wet gas window in the northern Levant Basin. However, only 10%-20% of their initial kerogen has been transformed yet.
Although the presented results are in good agreement with onshore observations and seismic interpretations of the offshore area, it has to be emphasized here that this is still a highly conceptual model, based on many assumptions. Such models of the Eastern Mediterranean region have to be updated according to new information, which will come up in the future. With respect to the presented model, among others, simulated source rock distribution and simulated temperatures represent major uncertainties. Certainly, new offshore well data would provide valuable information in this context. A highly recommended well location, in terms of gaining new insights on the offshore area, would be close-shore, offshore Lebanon, with a depth reaching at least down to the Upper Cretaceous. Presence, thickness, TOC and HI data for the source rock intervals would provide data of great relevance for understanding the petroleum system, improving the presented stratigraphic forward model. Maturity (vitrinite reflectance) data would allow a better palaeo temperature calibration. In addition, information on the younger strata could provide evidence on gas generation and also provide information on the quality of the Oligo-Miocene reservoir intervals.