Stratigraphy, Paleogeography and Depositional Setting of the K–Mg Salts in the Zechstein Group of Netherlands—Implications for the Development of Salt Caverns

: The 1 km thick evaporitic Permian Zechstein group in the Netherlands is subdivided into 5 halite rich evaporitic sequences including K–Mg salts (polyhalite, kieserite, sylvite, carnallite and bischoﬁte) for which the position in the Zechstein stratigraphy is still poorly constrained. Understanding the repartition of K–Mg salts is especially important for the development of salt caverns which require a salt as pure as possible in halite. By compiling well log and seismic data in the offshore and onshore domains of the Netherlands, regional cross-sections and isopach maps were performed in order to update the lithostratigraphy of the Zechstein group by including the K–Mg salts. Results enable (i) to propose paleogeographic maps representing the spatial repartition and the thickness variations of one to two K–Mg rich intervals in each evaporite cycle, (ii) to constrain the depositional setting of the different type of salts and the hydrological conditions which inﬂuenced the Zechstein stratigraphic architecture and (iii) to develop over the Netherlands risking maps assessing the risk of encountering K–Mg salts in salt pillows or salt diapirs eligible in term of depth and thickness for the development of salt caverns.


Introduction
The Netherlands subsurface displays large accumulation of salt belonging to the Zechstein group. The evaporites deposited during the Permian over a salt giant basin that was extending from the east of the Scottish coasts to the westernmost Belarus and southern Latvia in the east ( Figure 1A) [1][2][3]. The Zechstein group was up to 2 km thick and drove intensive halokinetic deformations since the Late Triassic [4,5]. The Zechstein evaporites include four main types of lithologies which are marking 5 to 7 basin-scale evaporite sequences: shale, carbonate, sulfate, halite and K-Mg salts (dominantly including polyhalite, kieserite, carnallite and bischofite) ( Figure 1B) [6][7][8]. Carbonate and sulphate deposits formed prograding platforms in the marginal part of the basins, whereas halite accumulated in basin centers during drawdown periods. K-Mg salts accumulated in time periods marked by drawdown maxima at salinity concentrations that were more than 70-90 times that of the original seawater [9]. The Zechstein group played a key role for many petroleum systems: it acts as the main seal of pre-salt reservoirs, intra-salt petroleum systems have been reported with the Zechstein carbonates acting as the main reservoirs and the halokinetic deformations of the salt shaped the structural trap of many post-salt reservoirs [10]. Magnesium and potash mineral resources are also exploited from the mining and leaching of K-Mg salts in the Zechstein [11][12][13].  [4] and reference therein) and localization of the Netherlands onshore and offshore domain (red line). (B): General lithostratigraphic chart of the Zechstein in Netherlands (after [14]).
For the storage of oil and natural gas, salt caverns have also been successfully created in diapirs of the Zechstein evaporites, such as in the Zuidwending salt dome, where five salt caverns about 550,000 m 3 each (radius ~25 m * height ~300 m) have been developed at 1000 m depth and are able to provide a working gas volume of ~45 million m 3 of natural gas per cavern [15,16]. In the coming decades, several scenarios may cause the development of more salt caverns in the Zechstein, for the storage of green hydrogen especially (hydrogen produced by electrolysis of water) [16,17]. Indeed, H2 is a promising energy resource, allowing us to significantly reduce the emission of greenhouse gas. The production of green H2 would rely on an eolian electricity resource, dependent of climate conditions, and thus prone to intermittency availability. Salt caverns would then offer the opportunity to (i) store the excess of produced H2 in periods marked by favorable wind conditions and (ii) use the stored H2 to compensate the deficit of produced H2 in periods marked by unfavorable wind conditions [17,18]. To date, all existing salt caverns were developed onshore. However, recent studies highlighted the opportunity of developing salt caverns in the offshore domain [19][20][21] in the Netherlands, especially where half of  [4] and reference therein) and localization of the Netherlands onshore and offshore domain (red line). (B): General lithostratigraphic chart of the Zechstein in Netherlands (after [14]).
For the storage of oil and natural gas, salt caverns have also been successfully created in diapirs of the Zechstein evaporites, such as in the Zuidwending salt dome, where five salt caverns about 550,000 m 3 each (radius~25 m * height~300 m) have been developed at 1000 m depth and are able to provide a working gas volume of~45 million m 3 of natural gas per cavern [15,16]. In the coming decades, several scenarios may cause the development of more salt caverns in the Zechstein, for the storage of green hydrogen especially (hydrogen produced by electrolysis of water) [16,17]. Indeed, H 2 is a promising energy resource, allowing us to significantly reduce the emission of greenhouse gas. The production of green H 2 would rely on an eolian electricity resource, dependent of climate conditions, and thus prone to intermittency availability. Salt caverns would then offer the opportunity to (i) store the excess of produced H 2 in periods marked by favorable wind conditions and (ii) use the stored H 2 to compensate the deficit of produced H 2 in periods marked by unfavorable wind conditions [17,18]. To date, all existing salt caverns were developed onshore. However, recent studies highlighted the opportunity of developing salt caverns in the offshore domain [19][20][21] in the Netherlands, especially where half of the storage potential is located offshore (>50 km away from the coast [18]). Offshore storage in salt caverns could have the benefits of (i) increasing the potential number of salt caverns without being embarrassed by surface installation already in place; (ii) reducing the environmental and human risks, in case of leakage, earthquake-induced events, H 2 contamination of underground drinking water or any operational issues; (iii) avoiding the long-term monitoring of the salt cavern once they are abandoned and (iv) increasing the social acceptability of underground H 2 storage [22].
Salt cavern development can be performed in bedded or diapiric evaporites, provided that some geological safety criteria are fulfilled. The most important ones are (i) an evaporite thickness of at least~200 m [17,18,23], (ii) a depth of the top salt not exceeding 1700-2000 m [24] and (iii) a composition as pure as possible in halite [25,26]. Modern seismic imaging tools are relatively efficient to constrain the thickness, shape and depth of salt deposits and salt structures. However, constraining the intrasalt lithologies requires a good knowledge of the targeted evaporite formation and of its stratigraphic architecture especially [27,28]. In the case of the Zechstein group, a very large number of publications constrained the primary stratigraphic architecture of the evaporites, their composition, their spatial repartition and their halokinetic deformations [29][30][31][32][33][34][35][36]. However, if marginal platforms have been the focus of many stratigraphic and sedimentologic studies [35,[37][38][39][40][41][42][43], the spatial distribution and thickness variations of the K-Mg salts within halite remain poorly constrained. This is illustrated by the fact that, in the published stratigraphic charts of the Zechstein group in the Netherlands, K-Mg salts and halite are always undifferentiated ( Figure 1B). However, predicting the occurrence of K-Mg salts in bedded or diapiric salt is critical when dealing with the development of salt cavities as this can influence the shape and stability of the cavern during leaching and storage operations.
Accordingly, this study uses bibliographical, well log and seismic data to propose a geological synthesis of the Zechstein salt basin in the Netherlands with a focus on the K-Mg salt deposits. The main goals of this data synthesis are (i) to improve the stratigraphy of the Zechstein salts by constraining the spatial architecture of the halite and K-Mg salts in the onshore and offshore domains of Netherlands, (ii) to review the depositional model of the Zechstein salts by discussing the influence of tectonic and hydrological conditions on the stratigraphic architecture and (iii) to provide maps allowing us to assess the risks of encountering hyper-soluble K-Mg salts in salt structures eligible for salt cavern development.

Tectono-Stratigraphic Setting
The Central European Basin developed from the End of the Variscan Orogeny with a post-orogenic collapse and a following rifting [3,44]. Crustal cooling drove a regional subsidence during the Permian coupled with minor and localized extensional faulting [45,46]. In the Dutch area, this Permo-Carboniferous tectonics activity subdivided the study area in structural lows and highs bounded by large Variscan fault zones which strongly influenced facies and thickness variations of the Permian to Mesozoic deposits [5,[47][48][49].
Early to Late Permian facies (Rotliegend Formation) deposited in lacustrine to fluvial and arid continental depositional environments [50]. The Zechstein salt basin developed in the Late Permian, when a global glacio-eustatic sea-level rise, coupled with active rifting in the North Sea, allowed the development of a seaway between the Permian Basin and the Artic Ocean. During 2.8 to 7 My [3,44,[51][52][53], marine influxes in the Central European Basin provided the brines necessary for the development of the Zechstein group and the thick accumulation of evaporites under arid conditions. The Zechstein evaporite factory was influenced by eustatic variations which subdivided the stratal architecture of the formation in five evaporites cycles commonly named Z1 to Z5 and further detailed below ( Figure 1B).
Post-salt sedimentation includes up to 4000 m of Triassic to Neogene marine and continental deposits. Halokinesis of the Zechstein initiated in the Early Triassic under extensional strains. An important Late Jurassic uplift event induced the erosion of Triassic to Zechstein deposits above structural high, especially at the Texel Ijsselmeet High [47] and in the northern salt province of the Netherlands offshore domain. Compressive strains occurred from the Late Cretaceous and induced the deformation and quizzing of many salt structures [4,33]. Finally, a Cenozoic inversion reactived the halokinetic deformations [54]. Nowadays, this multi-phased halokinesis results in a large number of salt pillows, diapirs and walls all over the Permian Basin [5,54] (Figure 2). Post-salt sedimentation includes up to 4000 m of Triassic to Neogene marine and continental deposits. Halokinesis of the Zechstein initiated in the Early Triassic under extensional strains. An important Late Jurassic uplift event induced the erosion of Triassic to Zechstein deposits above structural high, especially at the Texel Ijsselmeet High [47] and in the northern salt province of the Netherlands offshore domain. Compressive strains occurred from the Late Cretaceous and induced the deformation and quizzing of many salt structures [4,33]. Finally, a Cenozoic inversion reactived the halokinetic deformations [54]. Nowadays, this multi-phased halokinesis results in a large number of salt pillows, diapirs and walls all over the Permian Basin [5,54] (Figure 2).

Stratal Architecture of the Zechstein group in Netherland-State of the Art
The Zechstein evaporites display a succession of shales, carbonates, sulfates, halite and K-Mg salts organized in evaporite sequences which were at first order controlled by the amount of the seawater inflows in the basin, itself linked to third or fourth-order sequences of eustatic variations [55][56][57]. Three formation boundaries act as maximum flooding surfaces: base Z1 with the Coppershale Member; base Z3 with the Gray Salt Clay Member and base Z4 with the Red Salt Clay Member ( Figure 1B). The carbonates represent

Stratal Architecture of the Zechstein group in Netherland-State of the Art
The Zechstein evaporites display a succession of shales, carbonates, sulfates, halite and K-Mg salts organized in evaporite sequences which were at first order controlled by the amount of the seawater inflows in the basin, itself linked to third or fourth-order sequences of eustatic variations [55][56][57]. Three formation boundaries act as maximum flooding surfaces: base Z1 with the Coppershale Member; base Z3 with the Gray Salt Clay Member and base Z4 with the Red Salt Clay Member ( Figure 1B). The carbonates represent the high-stand system tracks. The sulfate and capping halite and K-Mg salt represent the lowstand system tracks. The depositional thickness of the Zechstein group ranges from less than 50 m in the southern Netherlands to more than 1200 m in the northern offshore area [8]. It is hereafter sum up the stratigraphy of the 5 Zechstein sequences reported in the Netherlands, based on [8,58] and references therein.

Z1-Werra Fm
The Z1 comprises the basal Coppershale (Z1E), the Z1 Carbonate (Z1B), Z1 Anhydrite (Z1A) and Z1 Salt (Z1H). The formation starts with a 0.5 m thick, finely laminated claystone (Coppershale or Kupferschiefer), capped by the Z1 Carbonate, a carbonate platform up to 200 m thick which developed onshore-ward and along the Mid North Sea High [30] and which was laterally grading to 8-10 m thick limestones and dolomite in sediment-starved basinal setting [35]. With increasing saline conditions due to the lowering of the sea level, the Z1 carbonates graded at the basinward side to an anhydrite platform which was locally up to 300 m thick. Offshore ward, in the main basin, the precipitation rate of sulphates was limited under bottom anoxic conditions by bacterial sulphate consumption and the anhydrite was~40 m thick. On the anhydrite platform, a series of fault bounded depressions allowed the local accumulation of up to 300 m of halite during the last evaporitic stages [4]. K-Mg salts in the Z1 were locally reported in the eastern Netherlands [59]. A sheet-like upper anhydrite unit, developed in shallow water to sub-aerial conditions, ends the Z1 and marks the marine transgression of the following Z2 evaporite cycle.

Z2-Stassfurt Fm
The Z2 comprises the Z2 Carbonate (Z2C), the Z2 Basal Anhydrite (Z2A), the Z2 Salt (Z2H) and the Z2 Roof Anhydrite (Z2T). It is represented in the southern onshore by anhydrite bearing claystones, against sandstone deposits in the western offshore [60]. The Z2 carbonate unit is a large platform, up 80 m thick and 70 km wide, having prograded northward on the Z1 platform [61] and grading to~10 m thick bituminous limestone/dolomite in the main basin. With increasing salinities, the Z2 anhydrite formed a platform up to 40 m thick laterally grading to a 3 m thick layer in basin center. More than 600 m of halite interlayered with K-Mg salts accumulated in the main basin, filling the depression beyond the former carbonate-anhydrite platform, above which salt deposits remained relatively thin (10 to 100 m thick). Based on intrasalt correlative markers, Geluk and Rohling [34] subdivided the Z2 salt in three members that they named the Lower, Middle and Upper Z2 salts. During the depositional stages of Lower and Middle salt member, water depth in the central basin was on the order of 140-200 m [62,63] and "deep-water" polyhalite and carnallite have been reported [8]. During the last stages of halite infill, the basin was flattened and sylvite and carnallite salts regionally developed at the top of the Z2 salt. Finally, the Z2 ends with an anhydrite layer (Z2 Roof Anhydrite), recording the progressive inflow of seawater.

Z3-Leine Fm
The Z3 comprises the Gray Salt Clay (Z3G), Z3 Carbonate (Z3C), Z3 Main Anhydrite (Z3A) and Z3 Salt (Z3H). The Gray Salt clay is 5 to 10 m thick. The Z3 carbonate platform developed in a lowered topography as compared to the previous Z2 carbonate [30]. Slope facies thickness reach 40 m whereas basinal organic-rich limestones are only few metersthick. The Z3 Main Anhydrite form a platform up to 100 m thick north of the Z3 Carbonate and about 50 m thick all over the basin center. The Z3 salt is estimated to have been up to 400 m thick. It includes a basal part made up of halite and an upper part involving K-Mg salt deposits with kieserite, carnallite and sylvite. It also includes bischofite deposits up to 10 m thick, but only reported in the north-eastern onshore and north-western offshore area [14]. In the onshore domain, part of the depocenters in the Z3 were controlled by salt tectonics activity and related salt movements of the Z2 salt [36].

Z4-Aller Fm
The Z4 comprises the Red Salt Clay (Z4R), the Z4 Pegmatite Anhydrite (Z4A) and the Z4 Salt (Z4H). The two basal units, about 2 m thick each, are regional whereas the Z4 salt, up to 150 m thick, is only found in depocenters which were possibly controlled by tectonic movements [36,64]. The Z4 includes K-Mg salt in its middle part, whereas the upper part is characterized by an alternation of halite and claystone marking the evolution of the salt basin into playa-type conditions [8]. Along the southern basin margin, the Z4 grades to sabkhas deposits and fluvial sandstones.

Z5-Ohre Fm
The Z5 comprises a several meters-thick basal claystone (Z5R), followed by up to 15 m thick halite deposits (Z5H). Z5 deposits are limited to very local occurrence in the northeast of the country and to the north-western offshore, outlining the last depocenters [65].
The Zechstein evaporites are finally disconformably capped by the Zechstein Upper Claystone Fm (ZEUC), up to 50 m thick and composed of red and gray anhydritic claystones and sandstones, deposited in a lacustrine to mudflat setting.

Materials and Methods
The following study relies on wirelines logs, drilling reports, seismic data and structural to isopach maps provided by the Geological Survey of the Netherlands and freely available on https://www.nlog.nl/ (11 April 2022). The 2110 wells which reached the Zechstein group were investigated. A set of wireline log including gamma ray (GR, especially useful to identify K-Mg salts because of the radioactivity of the 40 K), caliper, density, sonic (DT) and neutron porosity were interpreted on the Schlumberger's Petrel software for the identification of the intra-salt lithologies. These interpretations also considered the descriptions of cuttings produced during well drilling operations and reported in well reports.
The lithological interpretation coupled with the electrofacies patterns were used to perform well-by-well strata correlations. Considering the intensive halokinetic deformations which affected the Zechstein salts, none of the studied Zechstein intervals with thick halite deposits could be considered as being strictly undeformed. However, I have selected 518 wells for which the internal stratigraphy, over the entire Zechstein group, or in some strata intervals only, was estimated as being relevant for the study (Figure 2). This sorting considered: (i) the reference stratigraphic succession established by the bibliographical review in the Netherlands and extended from well to well, (ii) anomalous symmetric stratigraphic intervals suggesting the occurrence of folds [66], (iii) localized facies or thickness anomalies which could not be reported in surrounding nearby wells, (v) the isopach map and top salt depth map of the Zechstein group, allowing me to localize the main diapirs and welds [5], and (vi) 2D and 3D seismic lines passing through the wells, allowing me to assess the instrasalt halokinetic deformations.
The correlated beds were then used to establish isopach maps of the investigated stratigraphic intervals. The isopach maps do not involve the stratigraphic intervals that were identified as highly anomalous, thus limiting the influence of halokinetic deformations as much as possible. The spatial thickness variations of these maps were built following the isochore interpolation algorithm proposed by the Petrel software. This method is a control-point orientated algorithm that converges upon the solution iteratively, adding more resolution with each iteration. With this method, general trends are retained in areas with little thickness data, while detail is retained in areas where the data exist. Each map is presented with its constraining points, thus allowing me to assess the precision of the maps over the study area.

Results: Stratal Evolution and Facies Repartition the Zechstein Salt Units
Two representative cross-sections complemented by isopach maps in the different subunits composing the cycles are hereafter presented (Figures 3-6). As illustrated, it was possible to identify and correlate the Z1 to Z4 evaporite cycles of the Zechstein group over most of the study area. However, significative stratigraphic interpretations could not be performed in the half eastern part of the offshore domain, due to intensive halokinetic deformations coupled with a lack of representative wells in the evaporites. We used the isopach maps and cross-sections to detail the stratal evolution of the Zechstein group in the halite and K-Mg-rich intervals. As previously noted in the methodology section, the data compilation cannot clearly identify stratigraphic intervals in which nondestructive salt inflation or deflation occurred due to post-or syn-depositional flow. Accordingly, it is important to note that only the large-scale thickness and facies variations of the evaporites cycles were of interest. Small-scale localized anomalies in each section and map are beyond the scope of this study and were not investigated.

Z1 Salt Unit
The isopach map of the Z1 salt reproduces the previously published spatial and thickness repartition of this unit in the central onshore domain of the Netherlands [8] ( Figure 3; Figure 5B). In the main depocenters, the unit is relatively halite-pure with the occurrences of K-Mg salt layers reported in three wells only ( Figure 5B). These K-Mg salts form a strata interval named Z1K, which is up to~70 m thick in the middle part of the halite member (see well WSK 01 on Figure 3). Laterally, away from the thickest depocenters, the wireline interpretations highlight a high proportion of clay and anhydrite interbedded or mixed with the halite (e.g., well P5-2 on Figure 4).

Z2 Salt Unit
Based on polyhalite regional markers, the three salt depositional stages previously highlighted by Geluk and Röhling [34] in the Z2 salt could be correlated over the crosssections ( Figure 3; Figure 4). The Lower Z2 salt member, hereafter named Z2H1, is about 400 m thick in the eastern side of the Texel Ijsselmeet High, on the onshore domain, and 100 m thick in the northwestern offshore domain ( Figure 5B), resulting in a sigmoidal platform-like shape on the cross sections. In the lower proportion, the Z2H1 also seems more developed in the southern part of the offshore domain ( Figure 4). The cross sections also highlight that the Z2H1 is relatively pure halite southward, where it is thick, whereas it includes thin alternations of halite, polyhalite and, in a lower proportion, anhydrite in the basinal offshore domain.
The Middle Z2 salt member includes a basal halite-dominated member, hereafter named Z2H2, capped by a potash-rich interval named Z2K1. The Z2H2 is once again thickening southward, where it is~200 m thick and~100 m thick in the distal domain ( Figure 3; Figure 4). The cross-section correlations also suggest that the Z2H2 aggraded over the Z2H1 and the anhydrite platform (Z2A) in the southern domain. The Z2K1 is up to 260 m thick and only well developed in the central part of the offshore domain, along an E-W oriented trend ( Figure 5C). Wireline log responses suggest very local and thin occurrences of bischofite in the uppermost part of this subunit (Figure 4; Figure 5C).        Figure 6A for the localization of the strata intervals on the lithologic column).  Figure 6A for the localization of the strata intervals on the lithologic column).
The Upper Z2 salt member includes a third halite interval named Z2H3, capped by a second potash-rich interval named Z2K2. On the opposite of the Z2H1 and Z2H2, the Z2H3 is further developed distal-ward with thickness reaching~300 m in the offshore, against 70 m on the onshore domain ( Figure 4). Finally, the Z2K2 forms a widespread deposit, up to~90 m thick in the offshore domain and subdivided in 1 to 50 m thick more or less isolated depocenters in the southern part of the study area, especially above the platform domain ( Figure 6A).

Z3 Salt Unit
The Z3 salt unit includes three halite-dominated subunits, named Z3H1, Z3H2 and Z3H3, interlayered with two potash-rich subunits, named Z3K1 and Z3K2 ( Figure 5A). Z3 salts are lacking in the central offshore domain above the Clever Bank High (Figure 4). This area was significantly uplifted during the Late Jurassic [5] and the lack of Z3 salts there probably results from post-salt erosion rather than non-deposition on topographic high. The Z3H1 subunit tends to thicken offshore-ward, where it is~200 m thick and locally interlayered with very thin layers of anhydrite to potassic salts ( Figure 3). The Z3K1 is few meters to~50 m thick and dominantly developed in the northwestern offshore and in the northeastern onshore, around the Lauwerszee Trough and Groningen High area ( Figure 6B). In accordance with previous publications [7,31], well log data highlight the occurrence of several meters-thick bischofite-rich layers in the main depocenters of this interval ( Figure 3).
The Z3H2 halite subunit is relatively tabular and~25 m thick. The Z3K2 interval is well extended over the study area, with a thickness up to~90 m in the northwestern offshore and the northeastern part of the onshore ( Figure 6C). In the western part of the offshore domain, cutting data suggests the potassic salt to be relatively shale rich. This shale content gives sonic log response which can be confused with bischofite layers. However, the latter are suspected to be locally present (Figure 3; Figure 6). The capping Z3H3 subunit displays a large extension, is up to 90 m thick and is also characterized a high shale content, according to well log and cutting data.

Z4 Salt Unit
The Z4 salt unit is~60 m thick over the main basin domain. The basal halite-dominated member, named Z4H1, is~20 m thick and locally interlayered with two to three, few metersthick potassic salt layers, mostly in the offshore domain ( Figure 3; Figure 4). The central K-Mg-rich interval, hereafter named Z4K, is regionally correlative and can be up to 40 m thick in restricted depocenters in the southern part of the offshore and in the east of the onshore ( Figure 6C). The capping shale-rich halite member, named Z4H2, is~25 m thick when complete.

Controlling Factors of the Stratal Evolution in the Z2, Z3 and Z4 Salt Units
Based on the stratigraphic correlations performed in this study, it is possible to propose an updated lithostratigraphic chart of the Zechstein group in the Netherlands including the repartition of the main K-Mg salts (Figure 7). Moreover, it is hereafter proposed and discussed new facies maps representing the depositional setting of each K-Mg salt deposits (Figures 8-11).        The presented geometries of halite and K-Mg salt bodies might have been controlled by (i) sedimentary processes linked to hydroclimatic conditions or eustatic variations (e.g., [67,68]), (ii) syn-to post-salt halokinetic deformations (e.g., [30,69]) and/or (iii) basement tectonic activity (e.g., [70]). The alignment of the Z1H unit with NW-SE faults ( Figure 2) and seismic observations already demonstrated the tectonic control on this halite depocenter ( Figure 8A) [4]. We hereafter focus the discussion on the stratigraphic architecture of the Z2, Z3 and Z4 salt units in order to understand (i) the origin of the halite thickness variations in the Z2H1 and Z2H2 subunits, (ii) the origin of the anhydrite and polyhalite deposits in the Z2H1, (iii) the depositional setting of the K-Mg salts in the Z2, Z3 and Z4. The presented geometries of halite and K-Mg salt bodies might have been controlled by (i) sedimentary processes linked to hydroclimatic conditions or eustatic variations (e.g., [67,68]), (ii) syn-to post-salt halokinetic deformations (e.g., [30,69]) and/or (iii) basement tectonic activity (e.g., [70]). The alignment of the Z1H unit with NW-SE faults ( Figure 2) and seismic observations already demonstrated the tectonic control on this halite depocenter ( Figure 8A) [4]. We hereafter focus the discussion on the stratigraphic architecture of the Z2, Z3 and Z4 salt units in order to understand (i) the origin of the halite thickness variations in the Z2H1 and Z2H2 subunits, (ii) the origin of the anhydrite and polyhalite deposits in the Z2H1, (iii) the depositional setting of the K-Mg salts in the Z2, Z3 and Z4.

Origin of Thickness Variations in the Z2H1 and Z2H2 Units
The Z2 salt unit filled a large basin trough, having remained relatively sedimentstarved throughout the development of the Z1-Z2 carbonate and sulphate platforms [4]. East of the Texel-IJsselmeer high, the large thickening of the Z2H1 and Z2H2, resulting in their sigmoidal shape under the N-S cross-sections, was interpreted by Geluk and Röhling [34] as a depositional halite platform, having prograded over the previously exiting relief. East of Southern England, in the stratigraphically equivalent Z2 interval, a similar platform geometry was reported by Colter and Reed [71] and Taylor and Colter [72]. There, the authors suggested that the halite platform formed due to a higher precipitation rate in the shelf domain, induced by higher temperatures, evaporation and, thus, salinities. Outside the English and Dutch basins, halite platform-like geometries have never been reported in the Zechstein basin, and neither in other ancient salt giant basins, where halite deposits tend to form aggrading tabular geometries in basin depocenters [24,42,55,67,73,74]. Moreover, recent studies in the 300 m deep Dead Sea have highlighted that, due to downward saline flux transfers in the water mass, the halite precipitation rate in deep basins should be increased in the deep water domain, at the expense of dissolution in the shallow basin margins [75], i.e., the very opposite of what should drive the development of a halite platform. Nevertheless, in the southern side of the Dead Sea, a modern platform-like geometry of halite has recently been recognized [68]. There, it developed because of a salinity gradient (and related precipitation rate gradient) decreasing from the south to the north due to freshwater influxes provided by the Jordan River in the northern side of the lake [68]. Accordingly, by analogy, I suggest that the Z2H1 and Z2H2 units accumulated under a deep-water setting, with the influence of freshwater inflows having reduced the halite precipitation rate in the western offshore domain of Netherlands (and in the basin domain of England), as compared to the onshore one ( Figure 8B). Freshwater influxes could be sourced from continental river water or, more probably, by marine waters coming from the northern North Sea. The Texel-IJsselmeer high probably also acted as a structural barrier, having further confined the eastern onshore area, contributing to higher salinity conditions there.
Finally, following the stratigraphic correlations, the aggradation of the Z2H2 unit over the Z2H1 suggests that the Z2H2 developed during or after a marine recharge having risen the brine level, and probably linked to a slight increase of the eustatic sea level (Figure 7).

Origin of the Anhydrite and Polyhalite Layers in the Z2H1
In the Z2H1 unit, stratigraphic correlations highlight that the thick halite platform of the onshore domain is grading to a thinner halite interval thinly interlayered with polyhalite and anhydrite, to the origins of which must be assessed. Anhydrite can be primary when precipitated from sulfate or halite-saturated conditions, but is more commonly interpreted as resulting from the diagenetic dehydration of primary gypsum beds [76][77][78][79][80]. Similarly, polyhalite can be primary when formed as the first K-Mg salts [81][82][83] or secondary when formed through back-reaction between sulphates and evaporation-derived K-Mg-SO 4 -rich brines [84][85][86]. In the onshore domain of the Netherlands, Biehl et al. [6] argued that the anhydrite and polyhalite interlayered with halite in the Z2H1 were both depositional. However, in a stratigraphically equivalent interval in the English basin, Colter and Reed [71] reported clear petrographic evidence of both primary and secondary polyhalite. Accordingly, in the halite of the Z2H1, the gypsum-anhydrite probably formed during freshwater events, whereas polyhalite formed as a secondary to primary phase during drawdown events, while the halite platform could be emerged.
Considering the deep water setting during the deposition of the Z2H1, it also seems very likely that part of the anhydrite beds reported in the offshore domain mark turbiditic lobe deposits derived from the dismembering of the marginal anhydrite platform (Z2A, e.g., [87][88][89]. Indeed, several studies have already documented these types of basinal deposits in the Zechstein Basin, but only during the development of the Z1, Z2 or Z3 sulfate platforms [42,90,91]. However, the dismembering of the Z2 sulphate platform is not expected to have stopped at the first stages of the basin halite infill, because a topographic relief was still well-developed between the marginal domain and the basin trough. Accordingly, the occurrence of anhydrite turbiditic lobes interlayered in the distal Z2H1 halite deposits is also a very likely scenario (as illustrated on Figure 8B). Finally, such clastic anhydrite could have been lately transformed into secondary polyhalite, similarly to what was recently reported by Shang et al. [83] in halite deposits of the Sichuan Basin, SW China.

Depositional Setting of the Z2, Z3 and Z4 K-Mg Salts Deposits
The intervals of K-Mg salts interlayered with the halite in the Z1 to Z4 salt units record successive drawdown evaporite cycles. Considering the influence of the third-to-fourth order eustatic variations on the large-scale evaporite cycles, the smaller-scale halite-potash alternation hereafter documented are probably linked with fourth-to-fifth order sequences of eustatic variations (Figure 7). Lateral thickness and facies variations reported on the stratigraphic correlations may enable us to decipher the hydrological conditions which favored the deposition of such hyper-saline deposits.
After the deposition of the Z2H1 and Z2H2 halite platform in deep water setting, a probable drop of the eustatic sea-level decreased the seawater influxes, leading to the basin drawdown coupled with the emersion of the halite platform and, ultimately, to the progressive infill of the distal topographic low by the K-Mg-salts of the Z2K1 unit ( Figure 9A). If the K-Mg salt started to precipitate under a relatively deep-water setting (>15-20 m deep), the local occurrence of bischofite in the uppermost part of the potash layer suggests an almost complete desiccation stage reached at the end of the drawdown event. Indeed, bischofite is the last mineral able to precipitate from the evaporation of seawater [9,92], but due to its very high ability to be dissolved under atmospheric conditions, the accumulation and preservation of bischofite must have involved the deepest depocenters devoid of any outflow and with perennial subaqueous conditions [24,93]. Laterally, direct precipitation of carnallite, sylvite or kieserite could occur in ephemeral salt ponds or as subaerial crusts by crystallization of subsurface brines (e.g., [94,95]).
After a new marine flooding rose the brine level and returned the saturation stage to the halite one, the remaining topography was flattened by the deposition of Z2H2 halite unit. A second potash-saturation stage was then reached during the Z2K2 unit because of (i) a new drop of the eustatic level and/or (ii) newly developed hydrological conditions. Indeed, the new very flat configuration of the basin allowed the brine level to lie over an extended surface, promoting shallow-water condition all over the basin, net evaporation, and thus, higher salinities (the so-called "fill and spill" stage, sensu Warren [24]). Accordingly, the Z2K2 unit distributed over a widespread area and bischofite precipitation occurred in the deepest depocenters during the uppermost drawdown maxima ( Figure 9B). This interpretation is in accordance with the depositional setting already inferred for Z2 potash deposits stratigraphically equivalent in Poland and England [72,96].
During the deposition of the Gray Salt Clay, Z3 Carbonate and Z3 Main Anhydrite, the basin trough that had developed by subsidence was quickly filled by the deposition of the Z3H1 halite unit, as evidenced by its thickening offshore-ward. The Z3H1 thus probably started to develop in deep water setting (~100-150 m deep [63]), but ended in shallow water to subaerial conditions, as at the end of the Z2. As a result, the Z3K1, Z3H2, Z3K2 and Z3H3 also developed in a relatively flat and shallow depositional setting, in accordance with that inferred for the Z3 potashes reported in Poland, Germany and England [97][98][99][100].
Around the Groningen High area, Raith et al. [31] demonstrated that the restricted bischofite-rich depocenters of the Z3K1 subunit developed with the influence of basement tectonic movements coupled with early halokinetic deformations of the Z2 salt member. Barabasch et al. [36], Strozyk et al. [101] and Biehl et al. [6] also argued that there was significant tectonic and halokinetic movement during the deposition of the Z3 cycle in the onshore domain. Accordingly, at the scale of the entire Z3K1 depocenters, the limited extension of this K-Mg salt possibly results from similar larger scale tectono-halokinetic movements ( Figure 10A). However, the paleogeography and very shallow depositional setting probably also influenced the basin salinity there. Indeed, the marine influxes coming from the north (i.e., the Barent Sea) were progressively pre-concentrated during their thousand kilometers-long flow among the Zechstein sub-basins [99]. Accordingly, potash and bischofite deposits might have been especially promoted in the onshore domain of the Netherlands because it acted as among the most distal and sheltered depocenter for marine water influxes. Following this hypothesis, the Z3K1 depocenter was laterally equivalent with halite syn-depositionally developed in salt pans or isolated salt ponds. The more regional extension of the Z3K2 unit, as compared to the Z3K1, highlights a basin having become uniformly potash saturated with an overall brine level probably slightly higher than this of the Z3K1. Shallow water conditions coupled with subtle tectonic and halokinetic movements probably promoted the development of a mosaic of depocenters more or less interconnected and with lateral facies and thickness variations ( Figure 10B). Finally, this depositional setting must have remained stable up to the end of the Zechstein, before post-salt dissolution and halokinesis ( Figure 11).

Implications and Recommendation for the Development of Salt Caverns in the Zechstein Group
Using the isopach map of the Zechstein group and the depth map of the top Zechstein, it is possible to constrain the areas where the top of the Zechstein group is not deeper than 1700 m and where the salt is thicker than 200 m, two basic criteria for the development of salt caverns in salt deposits ( Figure 12). In the Zechstein, insoluble materials are dominantly represented by (i) the shale, carbonate and anhydrite layers marking the base of the Z3 unit (Z3G, C and A) and (ii) shale contents in the Z3H3 and Z4 halite deposits. These units being regionally present in the basin, they are expected to compose almost all the eligible salt structures. Regarding the K-Mg salts, thanks to the understanding of their spatial repartition, it is possible to assess the risk of encountering them in eligible salt structures. This is illustrated on Figure 12 by superposing the eligible salt structures for salt cavern and the K-Mg salt depocenters in the Z1 to Z3 cycles.  Bedded evaporites and salt structures located in the southeastern corner of the study area only involve the Z1 salt unit. This unit is halite-dominated and K-Mg salt layers are limited spatially. However, the occurrence of K-Mg salts in the middle part of a halite unit, which is only few hundred meters-thick, remains quite problematic for the development of a salt cavern [103] and potash deposits can still be involved in complex layer-parallel folds in bedded halite deposits [104]. Accordingly, developing salt cavern in the Z1 salts should be relatively safe outside the K-Mg depocenters only.
Further north, eligible salt structures involve the Z2 to Z4 salt units, each of them possibly involving up to two intervals of K-Mg salts. Assessing the risk of encountering K-Mg salts will thus depend on the (i) the spatial localization of the salt structure, (ii) the depth localization of the salt cavern and (ii) the type of salt structure.
Salt pillows dominantly involve salt inflation, brecciation, boudinage and large-scale folding [105,106]. As a result, intra-salt deformations may relatively preserve the primary stratigraphic architecture of the Zechstein in pillow structures. Because of this relative continuity of the bedding, thick K-Mg-rich intervals and the Z3 carbonate-anhydrite layer usually remain well identifiable under seismic images [7,31,102,105,107], even if high amplitude folds develop [101] ( Figure 13A). Accordingly, with good seismic data and the repartition maps of the K-Mg salts, it should be possible to anticipate, with a certain amount of confidence, the safest depth intervals and spatial locations for the development of salt caverns in pillow structures.
Intrasalt deformations in vertical diapiric structures are far more difficult to anticipate, especially due to seismic images remaining blind with the occurrence of sub vertical dips and complex three-dimensional curtain folds [108]. In these structures, intrasalt deformations will depend on several factors, including the stress fields during the different phase of salt flow, the primary stratigraphic architecture of the salt section (e.g., thickness and number of K-Mg rich intervals), the composition of the evaporites, the sedimentation rate during diapir growth or the shape of the diapir. Without performing numerous exploration wells, defining the best location for the creation of a salt cavern in a diapiric structure is, thus, a serious challenge for geoscientists. However, field observations in mined salt domes and numerical to analogical models highlight that, as a general rule of thumb, in diapirs remaining close from an anticlinal shape, the lowermost salt layers tend to rise upward and push outward the uppermost salt layers [27,101,104,109,110]. This deformation pattern was notably observed in mined German diapirs of the Zechstein ( Figure 13B-F) which are cored by the lowermost Z2 halite unit, surrounded by Z3 salts, whereas the Z4 salts are restrained to the outer limit of the diapir.
The Z2 salt involves the Z2K1 and the Z2K2 K-Mg salt units, but the Z2K1 is the most problematic: as it is relatively thick and located in the middle part of the Z2, its occurrence in the core of diapirs is very likely and its ability to develop complex folds with the halite is very high. However, as the Z2K1 is spatially restricted to the offshore domain ( Figure 12A), the risk of encountering this potash layer concerns the offshore diapirs only.
The Z2K2 is further regional but has the advantage of being thinner and at the uppermost part of the Z2, topped by the Z3 carbonate-anhydrite. Because the Z3 carbonateanhydrite has a limited ability to fold and mix with other salts, its occurrence atop the Z2K2 may limit the ability of the K-Mg salts to develop a complex folding pattern with the Z2 of Z3 halite, as observed in the German Gorleben salt dome ( Figure 13E). As a result, the Z2K2 is not expected to be a major risk in the central part of diapirs conventionally cored by the Z2 halite.
In the Z3 cycle, the Z3H1 should be the thickest halite unit and, thus, the best interval to safely develop a salt cavern. However, in diapirs cored by the Z2 salt, the Z3 salt unit might be less inflated than the Z2 salt and complex folds might develop with the two Z3K1 and Z3K2 intervals interbedded with halite in the half upper part of the Z3 cycle. Finally, the Z4 salts and its shale and K-Mg salt content should not be problematic as commonly restricted to the outer borders of the diapirs structures.  Unfortunately, all diapiric stems are not conventionally cored by the lower-most halite unit. Indeed, there are also numerous examples of large Zechstein diapirs for which the flowage and folding have been so complexed that the internal structures are not respecting any first order deformation pattern ( Figure 14). For this structure, the use of several exploration wells and geophysical tools become mandatary to target the best location for the development of a salt cavern. As a first order approximation, we might expect highly complex internal deformations for (i) diapirs with several km high vertical extensions, commonly bounded by welded salt stocks, (ii) several tens of km large diapirs with complex external shapes under map view, (iii) diapirs involved in compressional thrust sheet and (iii) vertical diapirs which, in the cross-section view, display a mushroom shape, salt sheets, a welded stem or evidence of collapse structures.
Unfortunately, all diapiric stems are not conventionally cored by the lower-most halite unit. Indeed, there are also numerous examples of large Zechstein diapirs for which the flowage and folding have been so complexed that the internal structures are not respecting any first order deformation pattern ( Figure 14). For this structure, the use of several exploration wells and geophysical tools become mandatary to target the best location for the development of a salt cavern. As a first order approximation, we might expect highly complex internal deformations for (i) diapirs with several km high vertical extensions, commonly bounded by welded salt stocks, (ii) several tens of km large diapirs with complex external shapes under map view, (iii) diapirs involved in compressional thrust sheet and (iii) vertical diapirs which, in the cross-section view, display a mushroom shape, salt sheets, a welded stem or evidence of collapse structures.   [113]). (B): Cross section of the Bartensleben diapir (after [116]). (C): Cross section of a diapir hosting the Mariglück Mine (after [117]). (D): Cross section of a unnamed diapir (after [118]). (E): Cross section of the Hänigsen-Wathlingen salt dome (after [119]).

Conclusions
In this study, bibliographical, seismic and well data were compiled to perform a geological synthesis of the Zechstein group in the Netherlands, with an emphasis on the K-Mg salts, which until now were poorly constrained in term of thickness and spatial repartition. The results enable us to update the litho-stratigraphy of the Zechstein and contribute towards improving our understanding of its depositional history. The main conclusions of this study can be summed up as follows: -Each Zechstein evaporite cycle (Z1 to Z4) includes one to two intervals rich in K-Mg salts which were named Z1K, Z2K1, Z2K2, Z3K1, Z3K2 and Z4K. - The Z1K was spatially limited to the tectonically controlled deepest depocenters of the Z1 salt unit. -During the first halite precipitation stage of the Z2, deep hydrological conditions coupled with marine freshwater influxes sourced from the northwest induced a higher precipitation rate of halite in the eastern onshore domain of Netherlands, resulting in the development of a 200 m thick halite platform. Basin-ward, the halite was~100 m thick and interlayered with primary to secondary polyhalite and anhydrite, some of which were also turbidite deposits sourced from the marginal Z2 anhydrite platform. - The Z2K1 accumulated from relatively deep to shallow water conditions in the distal topographic low created during the previous halite saturation stage, i.e., north of the halite platform which was emerged during the maximum drawdown. - The Z2K2, Z3K1, Z3K2 and Z4K had a more regional extension as they rather developed on a very flat basin, under shallow water to subaerial conditions marked by a mosaic of depocenters, some of which had been controlled by tectonic activity and related salt movements. -Wireline log data suggest that the four K-Mg-rich intervals in the Z2 and Z3 include local occurrences of bischofite developed in perennial subaqueous and hypersaline conditions, in the deepest depocenters around which complete desiccation probably occurred. -Finally, the constrained spatial repartition of the K-Mg salts in the Zechstein group, coupled with the structural characterization of the halokinetic deformations, enable us to assess the risk of encountering K-Mg rich salts in salt structures eligible in term of depth and thickness for the development of salt caverns. - The updated litho-stratigraphic architecture of the Zechstein group in Netherlands could be useful for other industrial purposes, such as oil and gas drilling operations or for the mining and leaching of K-Mg salts.
Funding: This research was funded by TotalEnergies through R&D projects on H 2 storage.

Data Availability Statement:
All data used in this study are freely available on https://www.nlog.nl/ (11 April 2022).
Acknowledgments: I thank TotalEnergies for their financial support and the supply of working software. I also warmly thank Sidonie Réveillon, Sabine Delahaye and Rene Geluk for their technical support during the study. I finally thank the four anonymous reviewers for their corrections and remarks which improved the final manuscript.

Conflicts of Interest:
The author declares no conflict of interest. The funders had a role in the design of the study. The funders had no role on the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.