Tsunami waves extensively resurfaced the shorelines of an early Martian ocean

It has been proposed that ~3.4 billion years ago an ocean fed by enormous catastrophic floods covered most of the Martian northern lowlands. However, a persistent problem with this hypothesis is the lack of definitive paleoshoreline features. Here, based on geomorphic and thermal image mapping in the circum-Chryse and northwestern Arabia Terra regions of the northern plains, in combination with numerical analyses, we show evidence for two enormous tsunami events possibly triggered by bolide impacts, resulting in craters ~30 km in diameter and occurring perhaps a few million years apart. The tsunamis produced widespread littoral landforms, including run-up water-ice-rich and bouldery lobes, which extended tens to hundreds of kilometers over gently sloping plains and boundary cratered highlands, as well as backwash channels where wave retreat occurred on highland-boundary surfaces. The ice-rich lobes formed in association with the younger tsunami, showing that their emplacement took place following a transition into a colder global climatic regime that occurred after the older tsunami event. We conclude that, on early Mars, tsunamis played a major role in generating and resurfacing coastal terrains.

It has been proposed that ~3.4 billion years ago an ocean fed by enormous catastrophic floods covered most of the Martian northern lowlands. However, a persistent problem with this hypothesis is the lack of definitive paleoshoreline features. Here, based on geomorphic and thermal image mapping in the circum-Chryse and northwestern Arabia Terra regions of the northern plains, in combination with numerical analyses, we show evidence for two enormous tsunami events possibly triggered by bolide impacts, resulting in craters ~30 km in diameter and occurring perhaps a few million years apart. The tsunamis produced widespread littoral landforms, including run-up water-ice-rich and bouldery lobes, which extended tens to hundreds of kilometers over gently sloping plains and boundary cratered highlands, as well as backwash channels where wave retreat occurred on highland-boundary surfaces. The ice-rich lobes formed in association with the younger tsunami, showing that their emplacement took place following a transition into a colder global climatic regime that occurred after the older tsunami event. We conclude that, on early Mars, tsunamis played a major role in generating and resurfacing coastal terrains.
The existence of an early Mars northern ocean 1-7 remains a fundamental mystery 8,9 . During the Hesperian Period (~3.71 to 3.37 Ga; ages herein based on Neukum chronology as given in Michael) 10 , Mars' ancient hydrosphere was apparently cold-trapped within vast systems of subsurface aquifers underneath a thick, ice-rich permafrost zone 7 . Groundwater outbursts at the end of the Hesperian may have generated catastrophic floods that produced an ocean in the northern lowlands, as evidenced by a deposit that covers most of this region and generally exhibits a roughly topographically equipotential margin [1][2][3][4][5][6][7]11,12 . Radar-sounding data are consistent with the deposit being comprised of mostly water-ice 13 . This deposit is identified as the Late Hesperian lowland unit (lHl) on the latest geologic map of Mars 14 . However, until now, the lack of wave-cut paleoshoreline features 9 and the presence of lobate margins 8,12 appeared to be inconsistent with the Late Hesperian paleo-ocean hypothesis. Our new geologic mapping in Chryse Planitia and northwestern Arabia Terra regions reveals previously undistinguished, older and younger members of the unit (lHl 1 and lHl 2 , respectively, Fig. 1A). Both members are bounded by south-facing lobes that are typically tens of kilometers in length and width; however, in Chryse Planitia these dimensions reach a few hundred kilometers in scale (Fig. 1B, Fig. S1). The lobes reach upland boundary surfaces distributed between approximately − 4087 m and − 3191 m of elevation Scientific RepoRts | 6:25106 | DOI: 10.1038/srep25106 (Fig. S2). These deposits embay dozens of streamlined promontories scattered over a surface area of ~570,000 km 2 (Fig. 1B,C).
In THEMIS night-time infrared images, the upper reaches of the older deposit that were emplaced along Arabia and Tempe Terrae (member lHl 1 (Fig. 1A)) appear thermally bright (i.e., rocky exposures) 15 and abruptly transition upland-ward into thermally dark (i.e., fine-grained sediments) 15 surfaces (e.g., Figs 2A,B and 3B). Close-up views show that the bright surfaces consist of boulder deposits, with individual boulders typically meters in diameter (Figs 2C and 3E, Fig. S3D). Exhumation of the boulder deposit from beneath ejecta blanket materials along impact crater rims (black arrows in Fig. 2B,C), as well as distinct onlapping contacts (e.g., Fig. S3C), show that the deposit overlies the thermally dark surfaces consisting of finer-grained materials (e.g., Figs 2A,B and 3B). Throughout spatially disconnected locations in the eastern part of northwestern Arabia Terra, the marginal parts of member lHl 1 cover low-slope ramps that are extensively dissected by NNW-trending (Fig. S4A) sets of aligned channels (e.g., Fig. 3A-D). These channels were first identified in Viking data (but only locally along Arabia Terra in association with an older "lowland unit A") 1 .
Upslope flows leading to the emplacement of the lHl unit are implied by the highland-facing orientation of the deposits' lobes as well as their relief gains, which commonly are a few hundred meters (e.g., Fig. 1A,B, Figs. S5, S6). These characteristics rule out emplacement by gravity-driven downslope moving flows such as debris, flood, glacier and lava flows. Uphill unidirectional winds can generate elongate aeolian deposits known as wind streaks. However, these deposits are largely composed of saltating sand-sized lithic particles that are deposited in scattered patches on the lee sides of topographic obstacles (typically impact craters), exhibit surface bedforms, generally cover hills and mesas situated along their paths, and mostly have length-to-width ratios >1 (ref. 16). In contrast, the lobes of member lHl 1 include boulders several meters in diameter (Figs 2C and 3E, Fig. S3D), and those of member lHl 2 appear to be mostly composed of water-ice 6,[12][13][14] . In addition, the lobes in both members diverge around numerous mesas (e.g., Fig. 1C) as well as broad rises (e.g., Fig. S3), and have length-to-width ratios mostly <1 (Fig. S1) (which is consistent with uphill flow along with substantial lateral spreading). Therefore, we propose that the two unit lHl members represent deposits emplaced by highly energetic, sediment-rich tsunami waves that originated from a Late Hesperian paleo-ocean.
In Deuteronilus Mensae, extensive troughs cut the boundary scarps covered by member lHl 1 . The troughs are locally intruded by member lHl 2 run-up lobes (e.g., Fig. 2D), indicating that they formed during the time interval separating the two tsunami events. Active resurfacing leading to the formation of these troughs likely lasted a few million years and could have been the result of Late Hesperian glacial erosion 17 . Crater-count statistics show that,  The boulder deposits of member lHl 1 drape over, and therefore postdate, the incision of adjoining highland channels (e.g., orange arrow in Fig. S3A), ruling out upland fluvial systems as possible discharge sources. Highly energetic, boulder-rich tsunami fronts on Earth show diversion around topographic obstacles as they propagate onshore 18 . Similarly, member lHl 1 boulder deposits exhibit well-defined landward lobate margins around broad promontories (Fig. S3A-C). Member lHl 1 boulders range from rounded to angular and are as much as ~10 m in diameter (Figs 2C and 3E, Fig. S3D), which are also characteristics of some terrestrial tsunami deposits 18 . Thus, we interpret the member lHl 1 lobes as made up of lowland and boundary clastic materials that were captured and transported by a tsunami wave, then beached farther inland as the wave lost its momentum.
Subsequently, we suggest that rapid gravity-forced backwash of the tsunami wave into the paleo-ocean dissected the channel systems on the marginal parts of member lHl 1 in the eastern part of northwestern Arabia Terra (Fig. 3A-D). These channels have remarkable similarities to terrestrial tsunami backwash channels; including the presence of aligned channel heads 19 (black arrows in Fig. 3C), perpendicular orientations to the reconstructed paleoshoreline 19 (Fig. S4A), streamlined bars composed of reworked boulders 20,21 (Fig. 3D,E), and widths ranging between ~50 and ~200 m (refs 19,22) (Fig. 3C,D). Parker et al. 23 observed a few of these parallel channel systems in Arabia Terra using lower-resolution image data, and they also interpreted them as tsunami backwash channels.
The lower terminations of the proposed backwash channels are generally truncated by younger scarps (Figs 2D and 3C, Fig. S4A). However, the identification of a possibly subaqueously emplaced sedimentary lobe adjoining the lower reaches of a set of these channels located at ~− 3795 m in elevation (Fig. S4B) provides an approximate upper boundary to the paleoshoreline from which the older tsunami propagated (Fig. 4A). The lowest margins of the mapped lHl 2 lobes are at ~− 4100 m in elevation (Fig. S2), which we have used as an upper bound to the paleoshoreline elevation from which the younger tsunami propagated (Fig. 4B). The elevation difference between the two paleoshorelines implies a decrease in ocean level of ~300 m, which could have taken place via evaporation/sublimation within several million years 6 .
Based on these paleo-oceanographic reconstructions, we estimate that the areas inundated by the older and younger tsunamis within the study region were ~8 × 10 5 km 2 and ~1 × 10 6 km 2 , respectively (Fig. S5). Measured typical run-up distances are tens to a few hundred kilometers for both the older and younger tsunamis, and their respective maxima reach ~529 km and ~650 km (Fig. S6). Overall, the morphometric characterizations of both tsunamis are strikingly similar. The slightly larger inundation area that was apparently covered during the younger event is consistent with the tsunami extending from a lower shoreline, and therefore, flowing over relatively smooth, older ocean and tsunami deposits. These run-up distances and inundation areas are enormous by terrestrial standards, which explain why the backwash channels exhibit lengths of ~20 km, while some terrestrial examples of backwash channel lengths produced by much smaller tsunamis range between ~200 and ~300 m in length 22 .
Our mapping (Fig. 1, Fig. S6) shows comparatively shorter run-up distances along the rougher and steeper cratered topography of the Arabia Terra boundary terrains, indicative of relatively lower wave heights and velocities, as predicted by tsunami numerical simulations 24 . These simulations also indicate that as the waves overflowed the Arabia Terra cratered boundary, their velocities would have abruptly dropped below the ~1 m/s threshold required to move multi-meter-scale boulders, explaining the occurrence of the boulder deposits in the region (Figs 2C and 3E, Fig. S3D). On the other hand, the more gentle slopes in Chryse Planitia would have resulted in a more gradual decrease in wave velocity, leading to the emplacement of more sorted sedimentary lobes, with their distal-most areas primarily consisting of finer-grained sediments. In addition, prior to their inundation by tsunami waves, the highland boundary surfaces were likely covered by extensive boulder fields, which would have been captured and redistributed by the waves, which is also another important factor accounting for the regional prevalence of boulder-rich lobes. In Chryse Planitia the tsunamis would have mostly propagated over gently-sloping plains that were largely made up of less bouldery outflow channel sedimentary deposits 14 .
The simulations also show that bolide impacts causing craters ~30 km in diameter would have generated tsunami waves with typical onshore heights of ~50 m and local variations from ~10 m to as much as ~120 m (ref. 24). Using run-up distances measured in 71 topographic profiles (Fig. S6), we have calculated the tsunami wave heights and find that they reasonably match the simulations' predicted ranges 24 (see supplementary calculations). In addition, whereas the simulations do not describe the hydrodynamic behavior of the backwash stage, the formation of several marine impact craters on Earth has also resulted in documented tsunami backwash channels 25 .
Using the surface area of the paleo-ocean's region included in the numerical simulation by Iijima et al. 24 (i.e., ~4.5 × 10 6 km 2 ) and the crater production function of Ivanov 26 , we find that ~23 marine impact craters ≥ 30 km in diameter would have formed within this part of Mars during the Late Hesperian Epoch (3.61-3.37 Ga) 10,27 . Of these, 7 fit in the diameter range of 30-35 km, which was used in the tsunami simulations 24 . The prediction is that, within the particular region of the ocean analyzed here, on average about 2 impact craters ~30 km in diameter formed every 30 million years during this time period. Therefore, within statistical constraints for the deposits' surface ages and for crater production rates, impacts can account for generation of both lHl members as tsunami deposits (see supplementary crater statistics).
Briny aqueous chemistry models show that the ocean could have remained in liquid form over millions of years, and consequently mostly free of an ice cover even during cryogenic climatic conditions 28 . Another geologic scenario invokes the formation of an ice-covered ocean soon after the ocean's emplacement 6 . However, no numerical simulations have been performed to detail the behavior of impact-related tsunamis 24 on these types of Martian marine environments.
High rates of marine, and subsequent periglacial 6,12,14 resurfacing, likely reduced the topography of the tsunami-generating crater structures. Such resurfacing can also explain the lack of well-preserved impact craters predating the Amazonian Period in the northern lowlands 12 . The frequency rate of ~30 km in diameter impact craters for the entire ocean's surface area (~24 × 10 6 km 2 , as determined by Head et al. 3 ) is one every 2.7 million years during the Late Hesperian. Although we have only identified evidence for two tsunami events in our study area, other regions in the northern plains likely experienced similar tsunami-related coastal resurfacing, perhaps associated with other impacts, huge landslides, or large marsquakes. Older but less extensive tsunami deposits may have been completely resurfaced by more recent events with similar run-up distances. Thus, the mapped tsunami margins comprise only the largest magnitude tsunami events located at the highest elevations.
Many of the lHl 1 lobes are mostly made up of lithic deposits and exhibit backwash modifications. In contrast, the landward-facing lobate termini of unit lHl 2 lack evidence indicative of a backwash phase subsequent to their emplacement. Like on Earth, the absence of backwash features associated with these flows could have been the result of the waves transitioning into sub-aerial sediment-laden slurry flows extending over low gradient surfaces 29,30 (supplementary calculations), which can also explain the presence of possible contractional folds along the margins of some of the member's lobes (e.g., black arrow in Fig. 2D). However, the lHl 2 lobes appear to be mostly composed of water-ice 6,12-14 , suggesting that the transition into slurry likely involved the formation and incorporation of a significant proportion of ice particles. In May 2013, the Saskatchewan Water Security Agency filmed an ice surge in the Codette Reservoir near Nipawin, Saskatchewan, Canada. The surge comprises a spectacular terrestrial analog of rarely observed catastrophic ice-rich flows leading to the emplacement of enormous lobate fronts, which are remarkably similar to those of member lHl 2 (video link included in ref. 31).We propose that these morphologic differences might be linked to colder environmental conditions following the first tsunami event.
Our mapping of two unit lHl members as tsunami lobes is consistent with the occurrence of two paleoshoreline levels of a receding Martian northern plains ocean during the Late Hesperian ( Fig. 4, Fig. S5). However, resurfacing by the tsunami waves has obscured the paleoshorelines, thus making rigorous testing of their equipotentiality impossible.

Mapping Methodology
Mapping in this investigation was performed using Esri's ArcGIS ® 10.3 software (http://www.esri.com/software/ arcgis). Embayment and overlapping relationships leading to the recognition of the outer margins of members lHl 1 and lHl 2 involved an integrated analysis of (1) thermal infrared image data (i.e., Mars Odyssey Thermal Emission Imaging System (THEMIS) night-time and day-time infrared image mosaics (100 m per pixel)), (2) visible image data (i.e., Mars Reconnaissance Orbiter Context Camera (CTX, (5.15-5.91 m/pixel)) images, and (3) Mars Global Surveyor Mars Orbital Laser Altimeter (MOLA, ~460 m/pixel horizontal and ~1 m vertical resolution) digital elevation models. In some areas, contacts are buried underneath ejecta blanket materials or are locally resurfaced; we mapped these sections as uncertain contacts ( Fig. 2A).