Evaluating the impact of the Storegga tsunami on Mesolithic communities in Northumberland

The Holocene Storegga tsunami, 8120–8175 cal a bp, resulted in run‐up heights of up to 3–6 m around mainland UK and coincided with a suggested large population decline in the coastally focused Mesolithic population in Northern Britain. At Howick, Northumberland, the site of a Mesolithic settlement, a nearby sediment deposit may be of tsunamigenic origin, but this is uncertain. Here, a numerical model was used to simulate the Storegga tsunami in Northumberland. Two scenarios of relative sea‐level change, and a third incorporating high tide, were simulated with mortality estimated within the intertidal zone for the Mesolithic sites in the region. The results showed that only with the addition of high tide could the sediment deposit site have been inundated by the tsunami. At Howick, mortality estimates varied but were up to 100% within the resource‐rich intertidal zone. The tsunami inundated a large area and would have led to the loss of key resources such as hazelnuts prior to the winter months. These combined effects would have probably been replicated throughout coastal settlements in Northern Britain, possibly leading to the contemporary population decline estimated to have occurred at this time.


Introduction
Tsunamis are large waves which can be generated by earthquakes, submarine or subaerial landslides, volcanic eruptions and impact events.Recently, public awareness of their impacts has increased drastically due to the devastating Indian Ocean and Tōhoku tsunamis of 2004 and 2011 respectively (Tufekci-Enginar et al., 2021).Unsurprisingly, tsunami research has proliferated since these events with a need to predict future events with greater accuracy so that effective land-use zoning and evacuation plans can be implemented (Cyranoski, 2012).The impacts of large tsunamis are severe but, fortunately, they are comparatively rare events.Therefore, studies of past events are required to provide information on the frequency and severity of past tsunamis in an area, allowing more accurate predictions of the impacts of future events.
Past tsunamis have produced varying impacts depending on the vulnerability and resilience of past populations.During a tsunami, the large water depths generated can cause humans to lose stability which results in an enhanced risk of mortality (Chen et al., 2019).The depth at which this occurs has been quantified using recent observations to construct fragility curves, which correlate water depth to mortality risk (Jonkman et al., 2008).Therefore, for past events where water depth is known these curves can be used to estimate mortality during the event.The mortality and associated impacts of a tsunami will also be affected by the vulnerability of a coastal population, with factors such as a low coastal population density and enhanced awareness and preparedness reducing the impacts of past events (Tufekci-Enginar et al., 2021).
Tsunamis also cause damage to buildings and loss of resources such as crops leading to starvation and disease.Records of historical tsunamis depict these impacts (Chester, 2001) as well as more ancient events such as the tsunami which followed the volcanic eruption of Thera in 3600 BC (Minoura et al., 2000).At the Minoan city of Palaikastro on nearby Crete, 'destruction layers' have been attributed to this tsunami (Bruins et al., 2008).However, other than this event, very little is known about the impacts of ancient tsunamis on contemporary populations (McFadgen and Goff, 2007).
Palaeotsunamis are defined as events that occurred 'prior to the historical period or for which there are no written observations' (IOC, 2019: p. 8).The primary source of evidence for these events comes from sediment deposits formed during the tsunami, which typically form in topographic lows such as marshes, lakes, estuaries and lagoons (Bondevik, 2019).An incoming tsunami wave both erodes, deposits and transports sediment as it inundates land (Goff et al., 2012).Further erosion and deposition occur during the backwash and subsequent tsunami waves, which can leave distinct layers within deposits (Nanayama et al., 2000).However, erosion, bioturbation and anthropogenic activity remove these deposits, with only 50% of deposits from the 2004 tsunami preserved after 4 years (Szczuciński, 2012).Where they survive, deposits can be traced inland, where they thin, to a maximum altitude or distance inland that is the run-up or run-in of the tsunami (Costa and Andrade, 2020).However, these estimates are minimum values as sediment deposits are only found for around 60% of the inundation distance (Goto et al., 2011;Abe et al., 2012).This information can be used to infer the magnitude and frequency of past events.
The run-up heights of the tsunami (Figure 1) are greatest (>20 m) in areas such as Shetland where the narrow valleys amplify the tsunami height (Bondevik et al., 2003) and lowest in distant areas such as Eastern Scotland (Smith et al., 2004).Numerical models have also been used to estimate run-up heights, beginning with the model of Harbitz (1992).More recent models have incorporated palaeobathymetry (Hill et al., 2014) and inundation on land (Bateman et al., 2021).These models indicate the importance of an accurate relative sea level (RSL) prediction on estimations of run-up heights.However, the difficulty in reconstructing RSL for this period provides considerable uncertainty to both the sedimentological and modelling estimates (Shennan et al., 2018).Furthermore, the influence of tide is uncertain but has been shown to significantly influence tsunami hazards elsewhere (Gao and Niu, 2022).Therefore, the considerable uncertainty in reconstructing the tsunami means that predictions of the impacts on the contemporary population are difficult.
During the Mesolithic period (11 000-6000 BP) it is thought people lived nomadic hunter-gather lifestyles.The earliest Mesolithic site in Britain is Star Carr (11 300-10 480 cal a BP; Figure 2A), whose remains include several structures and extensive stone tool assemblages (Milner et al., 2018).These tools had uses which reflected the nomadic hunter-gatherer lifestyle of following the seasonal availability of resources, including animals such as deer, elk and aurochs.The seasonal movement patterns of Mesolithic people are heavily debated with several interpretative models being speculated (Preston and Kador, 2018).The first models suggested movement between upland areas in the summer to coastal sites in the autumn and winter (Clark, 1972).Alternatively, within resource-rich areas such as estuarine or coastal locations mobility was thought to decrease (Binford, 1980).This is supported by the abundance of coastal resources such as shell middens in the archaeological record (Finlay et al., 2019).Groups of people may also have moved within a defined territory such as a river basin, focusing on these estuarine areas (Preston, 2013).Overall, the coastal and estuarine locations were the focus of seasonal movement.However, the social relationships and networks are poorly understood for this period, with recent advances in analysing ancient DNA, lithic assemblages and isotope analysis only beginning to reveal these social relationships within the European Mesolithic (Cucart-Mora et al., 2022;Kjallquist and Price, 2019).Nearly all the major Mesolithic sites that contain structures are found in coastal locations (Figure 2A).This includes five larger 'pithouses' dating to 10 300-9600 cal a BP (Waddington and Bonsall, 2016).These structures are around 30 m 2 in size, with the significant investment in time and resources required in construction suggesting a degree of sedentism (Gooder, 2007).Similar house remains are also found in Northern (Damm et al., 2020) and Western Norway (Astveit and Tossebro, 2023), with the latter associated with a sedentary behaviour with frequent occupation over millennia.Extensive hazelnuts deposits in the pit-houses (Robertson et al., 2014) suggest that occupation probably occurred during the autumn and winter months.
The Storegga tsunami coincided with the 8.2-ka event that led to a rapid sustained decrease in temperatures across parts of Northwest Europe (Alley et al., 1997;Alley and Ágústsdóttir, 2005).However, as with the tsunami, the cold event had spatially varying impacts with some coastal communities experiencing little impact such as those in northern (Blankholm, 2020) and southern Norway (Breivik et al., 2018).Therefore, the spatial interaction between the two events is complex, with areas experiencing differing impact levels from either event.Using proxies of the number of radiocarbondated sites, studies have shown there was a sharp population decline in Northern Britain immediately after the concurrent events of 8200 BP (Figure 2B; Wicks and Mithen, 2014;Waddington and Wicks, 2017;Mithen and Wicks, 2021).Similar decreases were also seen in Norway (Bergsvik et al., 2021), but not in areas of continental Europe distal from the tsunami (Griffiths and Robinson, 2018;Van Maldegem et al., 2021).The population decline in Northern Britain was estimated using the method of summed calibrated probability distribution (SCPD) (Shennan and Edinborough, 2007;Porčić et al., 2016).Whilst useful for obtaining population patterns for which there is little other evidence, this method is subject to several biases (Contreras and Meadows, 2014;Crombé and Robinson, 2014).These include a proliferation of sedentary sites where more radiocarbon dates can be obtained, an absence of sites that have not been preserved such as those under current sea level and plateaux in the radiocarbon curve, where calibrated ages have a larger range and greater uncertainty (Mithen and Wicks, 2021).Moreover, the tsunami may have removed traces of Mesolithic sites from coastal areas (Nyland et al., 2021).These biases were largely addressed during the most recent study by Mithen and Wicks (2021) with similar patterns of population emerging.However, several sites were occupied on either side of the period of the Storegga tsunami, including Castlandhill (Figure 2A; Robertson et al., 2014).There was also no marked change in behaviour during this period (Nyland et al., 2021).

Howick
Howick is located on the Northumberland coast and is one of the most important Mesolithic sites in Britain (Figure 3).The site contains a pit-house occupied for around 100-300 years from 9800 cal a BP, sited within a resource-rich location on a river estuary a few hundred metres from the coast (Bayliss, 2007).Extensive hazelnut accumulations suggest at least autumn and winter occupation, although year-round living may have been possible (Waddington, 2007).Howick shares similar material remains and position relative to the coastline with other Mesolithic sites in the wider region including Low Hauxley in Northumberland and East Barns, Castlandhill and Echline Fields in southeastern Scotland (Mithen and Wicks, 2018;Figure 2A).Mesolithic people at these sites probably occupied similar environments and utilized similar resources to those at Howick.
Occupation at Howick is dated to before the Storegga tsunami.However, nearby sites were later reoccupied such as the site of East Barns over 1000 years after the initial occupation (Engl et al., 2021).The site at Howick was also later reoccupied during the Neolithic, Bronze and Iron ages (Figure 3; Waddington et al., 2005).This shows that the resource-rich location was a favoured site of occupation for later societies.Therefore, Howick was likely to have been occupied when the tsunami occurred in 8120-8175 cal a BP.
A sediment core taken from Howick burn (Core A; Figure 3) identified a discontinuity due to a 'short-lived erosional event' occurring before 8100 cal a BP (Boomer et al., 2007).This event removed several thousand years of sediment from the stratigraphic layers below and the deposit is repeated in another undated core further upstream (Core B).This evidence suggests the Storegga tsunami produced this deposit.However, these deposits were much coarser than the typical Storegga tsunami deposit from Britain.Moreover, a different highenergy depositional event may have produced the deposit such as a storm surge or a landslide of the adjacent riverbank (Boomer et al., 2007).Therefore, the tsunami provenance is uncertain.Potential Storegga tsunami deposits have been found in Northumberland at Broomhouse Farm around 25 miles north of Howick near Lindisfarne (Shennan et al., 2000).However, evidence of the Storegga tsunami at Mesolithic sites across northwest Europe is currently sparse, only being found at one site in Inverness (Dawson et al., 1990) and another in western Norway (Bondevik, 2003).Neither of these sites was occupied at the time of the tsunami.Therefore, the deposit at Howick could represent rare evidence of the direct impacts of the tsunami on Mesolithic people.2023) is used to simulate the Storegga tsunami at Howick and the surrounding areas.A spatially varying resolution mesh is used with multiple sealevel reconstructions to assess the impact of the Storegga tsunami on this coastline.This will provide evidence as to whether the sediment deposit probably comes from the tsunami and provide an understanding of the impact on Mesolithic people at Howick and the wider region.

Numerical modelling
Numerical models that represent coastal ocean regions can be used to predict the altered hydrodynamics caused by the tsunami wave.Here Thetis (Kärnä et al., 2018) was employed as a 2D flow solver for simulating coastal flows (Angeloudis et al., 2018;Goss et al., 2019;Bateman et al., 2021;Hill et al., 2023), which is implemented using the Firedrake finite element Partial Differential Equation (PDE) solver framework (Rathgeber et al., 2016).Thetis solves the non-conservative non-linear shallow water equations: where u is the depth average velocity vector, H is the total water depth, η is the free surface height and v is the kinematic viscosity of the fluid.The Coriolis term is represented by   u f , where   u is the velocity vector rotated anticlockwise through 90°.Moreover, f = 2Q sin c with Q designating the angular frequency of the Earth's rotation and c the latitude.Manning's n formulation is represented through the fluid density ρ, the gravitational attraction g and the bed shear stress τ b : where n is Manning's friction coefficient of 0.025, which represents medium-grained sand.Thetis incorporates wetting and drying inundation processes implemented by the algorithm of Kärnä et al. (2011).Here models were created using a discontinuous finite element discretization (DG-FEM) using the PIDG-PIDG velocity-elevation finite element pair.Newton non-linear solver algorithms were used through the PETSc library to solve the discretized equations (Balay et al., 2016).

Mesh generation
For the three main models, QGIS 3.26 (QGIS, 2023) was used to define the domain for the mesh which was then converted into a format compatible with meshing via qmesh (Avdis et al., 2018).The landscape for the model was constructed using OS Terrain 5 DTM topography (5-m resolution) and Marine DEM 1 arc second bathymetry (~30-m resolution) (Digimap, 2022a(Digimap, , 2022b)).These were adjusted to ordnance datum.To adjust these data to incorporate tidal dynamics and RSL change over the past 8100 years, three scenarios were produced.Two separate models were created to account for RSL change.First, the RSL from the glacio-isostatic adjustment (GIA) model reconstruction of Bradley et al. (2011) was subtracted from the modern-day bathymetry/topography.Similarly, the relative height of the Sea Level Index Point (SLIP) at Howick (3.23 m) was subtracted from the modernday bathymetry/topography (Shennan et al., 2018).This created two different models of RSL change labelled GIA and SLIP (Figure 4B, C).Next, the influence of tides was considered by using the modern-day tidal range to artificially increase the tsunami height.This involved subtracting half the tidal range at Howick (4.2 m; Bicket et al., 2017) from the topography/bathymetry of the SLIP model, which is representative of the wave having occurred at high tide (Gao and Niu, 2022).
From each of these scenarios, a 20-m contour was extracted in QGIS to be used as the landward boundary to allow for wetting and drying to take place within the whole area that the tsunami may inundate.Narrow inlets with a width of less than 200 m or where they persisted for more than a few kilometres outside of the high-resolution region were removed.These would cause oscillations potentially leading to instability in the model (Avdis et al., 2018).Small areas of high ground and islands (with circumference <1000 m) were removed from the boundary conditions but still included in the model.At the northern and southern edges of the domain, a smoother arced forced boundary was added.The same forcing boundary was used for each mesh to ensure the boundary conditions remained the same (Figure 4A).
To convert the boundaries for mesh generation qmesh (Avdis et al., 2018) and Gmsh (Geuzaine and Remacle, 2009) were used.Each mesh was produced in the UTM30N coordinate reference system.The mesh resolution was determined via the depth and distance from the coastline and external boundaries.In similar high-resolution studies, a mesh resolution of 100 m was used at the coastline with a linear increase to 5 km distal from this boundary (Bateman et al., 2021).The larger-scale study of Hill et al. ( 2014) used a multiscale mesh starting at 500-m resolution at the coastline.In the present study, the minimum resolution at the coastline was taken as 75 m, with a sigmoidal increase after 500 m away from the coastline.Furthermore, at the external boundary, an element size of 1 km was used which increased linearly after 2 km from the boundary.The largest element size was 5 km at 25 and 50 km from the external and coastline boundaries.Depth (d) was used to control the mesh resolution (R) in a sigmoidal function: (4) This decreases resolution from 75 m at 0 m depth to just below 5 km at around 70 m depth.Finally, a high-resolution region was created for the study area around Howick, using OS Terrain 5-m DTM topography (Digimap, 2022b).A polygon around Howick was created using QGIS and a resolution of 15 m was used for the area inside this polygon.The resolution then decreased linearly with distance (D) from the edge of the polygon: This decreases resolution from 15 m at the polygon edge to 1 km at around 2 km from the boundary with the function only applying for 3 km around the polygon (Figure 4D).Overall, the resolution was taken as the maximum of the coastline, boundary, depth and high-resolution functions, with a minimum resolution of 75 m.This resulted in a maximum mesh size of 427 267 triangular elements and 220 303 nodes.Finally, the bathymetry and topography data were linearly interpreted onto the mesh via HRDS (Hill, 2019).
Each of the simulations used the larger-scale Storegga tsunami model of Hill et al. (2014) as the water elevation on the forced boundary.These forcing data were applied at each timestep, linearly interpolating between the output times of the regional model, beginning at 6.5 h into the Hill et al. (2014) model for a further 8.5 simulated hours.The boundary data were also bi-linearly interpolated onto the forcing boundary node points.Every model generated an output every 60 s and used a timestep of 1 s.The models were run with a varying alpha parameter for wetting and drying set between 0.5 and 75 and a viscosity of 10 m 2 s −1 .Manning's drag coefficient was set at 0.025.Each model outputted the initial bathymetry at the start of the simulation followed by the free surface elevation and velocity at each timestep, with a post-processing script outputting the maximum and average bed shear stress (BSS), velocity and speed and maximum elevation.The maximum elevation and bathymetry were used to show the inundation of the wave.

Grain size analysis
For each of the simulations where the tsunami inundated the core locations, the maximum BSS and velocity were compared to the critical shear stress and velocity needed to deposit the material of different size classes (Table 1).The Core A deposit consisted of coarse sands with pebbles and cobbles incorporated (Boomer et al., 2007).This suggests that coarse gravel was present in the deposit.Therefore, the maximum BSS and velocity were compared to the critical values for the coarse gravels.

Mortality estimation
The intertidal zone was crucial for Mesolithic people for fishing and shellfish collection, and therefore the mortality within this zone was estimated.Three functions were used to represent the relationship between water depth and mortality during flooding.The first was derived from observations from tsunamis (CDMC, 2003;Jonkman et al., 2008) where F d is the probability of mortality (0-1) and H is the wave height of the tsunami.Two other functions were considered which are represented by exponential and log-normal functions (Figure 5 Using QGIS, the maximum elevation of the wave in the intertidal zone at Howick was subtracted from the topography/ bathymetry to obtain the water depth.The mean water depth within the intertidal zone was used to estimate the mortality from Equations 6-8 respectively for the GIA and SLIP tsunami models.This was repeated for the other major Mesolithic site in the domain of Low Hauxley to provide a more regional perspective of mortality.

Wave dynamics
The primary tsunami wave from the submarine landslide is simulated to dissipate before it reaches northern England and is only represented by a very small peak in the northern area of the domain.However, the secondary wave from the initial tsunami is the largest wave observed at the study site, which refracts around the Shetland Islands (Hill et al., 2014).This wave arrives at the domain from the northeast ~2 h into the simulation, being preceded by a negative wave (Figure 6; see animations of Supporting Information Figures S1-4).This is ~8.5 h into the simulation of Hill et al. (2014).The simulated wave first impacts the northern areas of the domain including around Lindisfarne and the site of Nessend.The wave then travels southwest down the coastline, reaching Howick around 20 min later, followed by Low Hauxley in the southern area of the domain.At Howick, subsequent waves occur at 3.5, 5.5, 6.5 and 7.5 h (Figure 7; see animations of Figures S1-4).These waves appear to be caused by oscillations of the primary waves from different areas of the North Sea.The first is caused by oscillations from the largest wave, and the second originates in the northwestern region of the domain and moves south.The final two simulated waves originate from the north and south and the southeast of the domain.
The largest simulated wave heights were generated from the first wave.Maximum wave heights varied between 2.2 and 6.2 m across the domain, with the largest wave heights seen at Berwick (north of Nessend), south of Lindisfarne and a little south of Howick.Similarly, the maximum velocity varied throughout the domain, being greatest in passages between islands, such as the southern end of Lindisfarne and at the Farne Islands, where simulated velocities reached up to 6.7 m s −1 .At Howick the wave height at the two core locations varied between 3.33 and 5.12 m across the three models, with a maximum velocity of between 0.28 and 1.26 m s −1 (Figure 7).There was a large second peak at Howick as the first wave had not yet fully drained from the site, but overall subsequent waves dissipated in height.

RSL change reconstructions
The two different RSL change reconstructions showed different patterns within the results.Neither model inundated the two sediment core sites, with the wave height below the reconstructed sea level throughout the simulation (Figure 7A,  B).However, the SLIP model is far closer to inundating the core site, being 1.68 m lower than the land surface (Figure 7B).The subsequent waves gradually decrease in height with five distinct peaks in the SLIP model and only four in the GIA model.The SLIP model inundates a far larger area across the domain than the GIA model, as expected from their relative sea level positions (Figure 7).

High tide model
With the addition of a high-tide scenario to the SLIP model, Core A was inundated with a depth of 0.18 m, but Core B was not inundated (Figure 7C).Core A was only inundated for 3 min of the simulation.To ascertain whether a deposit could have formed, the BSS and velocity around Howick were analysed and compared to the values from Table 1.The maximum grain size that could have been deposited or transported around Howick is shown in Figure 8.At Core A, the BSS was high enough to be able to move and deposit sediment up to coarse gravel in size, whereas the velocity indicates a maximum grain size of fine gravel could have been deposited.

Inundation at Howick
The surrounding area around the site of Howick was inundated to differing degrees with each model.The GIA model inundated a very small area of land, being less significant than a high spring tide.On the other hand, the SLIP and Tidal models inundated a far larger area covering from the high tide mark to the modern cliff line (Figure 9A, B).Within the intertidal area, the maximum water depths are large with an average depth of 5.2 m but up to 7.4 m at the high tide mark (Figure 9A).These depths are amplified for the Tidal model, with a mean depth of 7.1 m across the intertidal zone.Finally, within the intertidal zone, the maximum velocity varied between 0.65 and 3.5 m s −1 for the SLIP model and between  1 and 2.8 m s −1 for the Tidal model (GIA: 0.33-3.5 m s −1 ).At Low Hauxley, the other major Mesolithic site within the model domain (Figure 2A), water depths in the intertidal zone were smaller but still significant for the Tidal and SLIP models (Figure 9C, D).Here, the average water depths were 4.1 and 6.2 m respectively for the two models.The inundated land area was larger here, being around 500 m further than the high tide mark.

Mortality risk
The risk of mortality in the intertidal zone was estimated using the three mortality estimations.Across these different estimations, there was a large range of mortality risk.However, there was very little variability between the GIA and SLIP models (Fig. 10).The highest mortality was predicted by the exponential estimation with a median mortality risk of 1 in the intertidal zone for both tsunami models.However, this estimation also had a large range encompassing the full scale of mortality from 0 to 1.The log-normal estimation was less extreme than the exponential estimation with mortality ranging between 0.25 and 0.9, with a median of 0.67 and 0.71 respectively for the two models within the intertidal zone.
Conversely, the CDMC estimation predicted a lower mortality with a median of 0.096 and 0.089 for the SLIP and GIA respectively with a small spread of 0.1 mortality risk.At Low Hauxley, the intertidal zone produced similar mortality estimates but these were generally slightly lower across all three estimations than those at Howick in the intertidal zone (see Supporting Information Fig. S5).The variability within both the exponential and log-normal estimations was also greater than at Howick.

Discussion
Were the core sites at Howick inundated by the tsunami?
The two sediment deposits from cores A and B (Boomer et al., 2007) could be direct evidence of the tsunami at the site of Howick.However, the largest initial wave of the simulated tsunami failed to inundate the sediment core locations for either of the two models incorporating RSL change.Inundation should extend beyond these locations as deposits only form for around 60% of the inundation distance (Goto et al., 2011;Abe et al., 2012).Therefore, the results suggest that either the deposits were not from the tsunami or that the model underestimated run-up heights.
The tsunami provenance of the sedimentary layer at Howick is uncertain.The dating and characteristic erosion of underlying stratigraphic strata suggest the Storegga tsunami could have left the deposit.However, the unusually coarse composition of the deposit casts doubt on this interpretation.Tsunami deposits typically consist of fine-grained material, generally with a maximum grain size of coarse sand (Costa and Andrade, 2020).Deposits from the Storegga tsunami elsewhere in Mainland Britain are fine-to medium-grained sand (Long, 2018).However, in this deposit gravels dominate (Boomer et al., 2007).An alternative source is from a landslide of the adjacent river valley sides that comprise similar sized particles (Boomer et al., 2007).Tsunamis can trigger and increase the probability of a landslide (Gill and Malamud, 2014) so the deposit could be an indirect rather than direct tsunami deposit.A storm surge is not likely to have been able to produce the deposit given that the largest storm during the last 70 years only reached a height of 1.2 m above the high tide mark, significantly lower than the deposit site (Haigh et al., 2015).The second deposit upstream, Core B, is in the same stratigraphic context as the Core A deposit, but the lack of dating of this deposit means the connection cannot be confirmed.Overall, the deposits are not certain to have been from the tsunami.
The underestimation of run-up height in past modelling studies has been most profound in narrow valleys and inlets.Narrow inlets are thought to cause a funnelling effect during tsunamis whereby the water is pushed upwards as it enters the narrow width of the valley (Didenkulova and Pelinovsky, 2011).This results in a larger wave height in valleys such as in Norwegian fjords compared to coastal locations (Vasskog et al., 2013).However, this is extremely difficult to model accurately, leading to large discrepancies in narrow river valleys in Shetland between the model estimates of 13 m runup height and the sedimentological estimate of 31.8 m (Dawson et al., 2020).Furthermore, Bateman et al. (2021) did not predict inundation of the site of Fullerton located in a narrow valley even with a detailed geomorphological reconstruction.Therefore, the narrow 100-m-wide valley at Howick may have experienced this amplification of run-up height that cannot be represented in the model presented here.This could explain the presence of Cores A and B beyond the inundation limits of the tsunami.
Underestimation may also have occurred due to the difficulty in reconstructing sea levels for the Storegga tsunami in Northumberland.The region is located on the 'hinge' of isostatic uplift, meaning it is uncertain if sea level has risen or fallen since the middle of the Holocene (Shennan et al., 2018).Furthermore, at the time of the Storegga tsunami there was a period of rapid sea level rise meaning that there is a larger range of values within the dating range of the tsunami (Bradley et al., 2011).This can be seen in the SLIP data used in the model for Howick which had a two standard deviation range of 1.21-5.25 m (Shennan et al., 2018).This range is greater than the 1.68 m difference between the elevation of the SLIP model and the sediment core.Therefore, the large uncertainty in the sea level estimate could account for the lack of inundation of the site.
Even with an accurate sea level position, the geomorphology will have altered, with features such as dunes, spits and beaches being formed in the last 8000 years.In the Montrose basin, a study using a similar methodology found only one of six tsunami deposit sites was inundated by the model (Bateman et al., 2021).However, a spit blocked the entrance to the river valley, and when this was removed five out of the six locations were inundated.Similarly, Hill et al. (2023) showed that removing more recent geomorphological features at the mouth of the Ythan River meant that the model and sedimentary data matched very well.At Howick, a headland and wave-cut platform block the entrance to the burn from the tsunami wave coming from the northeast, preventing direct inundation into the burn.However, evidence suggests the river valley would have extended further out to sea in the early Holocene with a wider river valley than present today (Bicket et al., 2017).This would allow the wave to penetrate directly into the river valley from the northeast.Further work could investigate the effect this has on inundation at Howick.
An alternative is that the tsunami wave was amplified as it coincided with a high tide.Tides have a major effect on coastal flooding, with the most severe storm surges coinciding with a high tide such as the 1953 North Sea floods (Wadey et al., 2015).The same principle applies to tsunamis.If a tsunami occurs at high tide it will have a greater inundation distance due to the increase in water height (Kowalik and Proshutinsky, 2010).The magnitude of this effect is not certain, but recent work suggests that the effect of tides is roughly equivalent to adding the high tide height to the tsunami run-up height (Gao and Niu, 2022).The tidal range may have been key in inundation of the Storegga tsunami to sites on the west coast of Scotland (Woodroffe et al., 2023).Here, the incorporation of the high tide into the SLIP model enhanced the inundation area, with the initial wave inundating the sediment Core A at Howick.Furthermore, the wave produced a BSS sufficient to transport coarse gravel, which the sedimentological data suggested was present in the deposit (Boomer et al., 2007).However, the velocity was only sufficient to transport coarse gravel a few metres from the core site.Overall, this suggests that the simulated wave could have created the coarse gravel deposit in Core A. However, the tsunami came close to but did not inundate Core B. The area of inundation near to Core B is disconnected from the sea.This is due to the model simulating the wave irrespective of whether it is above or below ground level so a depression inland can be inundated.Overall, the evidence from the tidal model suggests that had the tsunami occurred at high tide the Core A deposit could have been formed by the tsunami.
Numerical simulations of the Storegga tsunami have found an underestimation of run-up height in Shetland and an overestimation in western Norway (Hill et al., 2014;Løvholt et al., 2017).These differences could have in part been caused by a difference in the tidal regime.Tides are highly predictable and occur in semi-diurnal cycles in the North Sea (Vindenes et al., 2018).Therefore, tidal conditions can be predicted at other locations by assuming a high tide at Howick when the first wave impacts the site at around 8.5 h into the simulation of Hill et al. (2014).This wave impacts the Shetlands ~6 h earlier in this simulation.This roughly matches the just under 6-h difference in tidal regimes between Lerwick (Shetland) and North Shields (near Howick) (National Tidal Facility, 2023).Taking this together means that if the tsunami occurred at high tide at Howick, it would have too at Lerwick.This could explain some of the differences in run-up height seen in the Shetland Islands.Future work could use numerical modelling to quantify the influence of high tide on the tsunami run-up in Shetland.

Direct impacts
The intertidal zone during the Mesolithic a key area for marine resources including fish and shellfish.This is evidenced by the extensive and almost ubiquitous shell middens from Scotland during the Mesolithic (Finlay et al., 2019).Boats were used during the Mesolithic (Warren, 2000), as evidenced by the occurrence of bloodstone from an island in the Hebrides to the mainland 40 km away.However, it is unclear whether these boats would have been used for fishing activities.Therefore, fishing probably predominately took place in the intertidal zone along with shellfish and seaweed collection.So this key resource area could have been occupied by Mesolithic people when the tsunami occurred.
A tsunami's sudden onset gives no warning unlike storm surges that are preceded by a period of stormy weather (Khalfaoui et al., 2023).Past hunter fisher gatherer societies in tsunami-prone regions such as the Northern Pacific have shown a degree of resilience to the impacts of tsunamis, mitigating effects by moving to higher ground (Fitzhugh et al., 2016).However, this is unlikely to have been the case for Mesolithic people in Britain where there would have been no living memory of tsunamis (Long, 2018) and no appreciation of the danger as the sea receded.Indeed, the sudden withdrawal of the sea may have attracted people to the coast to collect the exposed shellfish and other resources.Even during the Indian Ocean tsunami of 2004, this behaviour occurred when thousands of people who had not experienced a tsunami before walked out onto the exposed seafloor when the water receded (Gregg et al., 2006).In the intertidal zone at Howick, mortality is predicted to have been variable but high, reaching up to 100% depending on the mortality functions used and assuming people were there at the time of the tsunami.Mortality in a flood event can be estimated from the water depth, but these estimates vary depending on the model that is fitted to the observational data (Jonkman et al., 2008).This variation is largest at high water depths (Ge et al., 2022), where observations are rare.This explains the high variability in mortality estimates within each simulation, which varied by up to 0.65 due to the large water depths.The median mortality risk peaked at 100% for the exponential function for both models but was variable with the median mortality at ~70% for the log-normal function.However, the third mortality function predicted a median mortality of <10%.Therefore, depending on the function used there could have been significant mortality in the intertidal zone.Occupation of this zone would be less prevalent at high tide as shellfish collection is more difficult.Therefore, the tidal model has not been included in our comparison of the mortality in the intertidal zone.Instead, mortality would be directed into the immediate coastal region where people would be foraging.At Low Hauxley, these estimates are echoed with slightly lower mortality overall at a similarly located site only a few hundred metres from the intertidal zone.Together, these mortality estimates represent the potential impact on Mesolithic people in Northumberland and could explain some of the population collapse seen at the time (Mithen and Wicks, 2021).

Indirect impacts
The tsunami would have also caused indirect impacts through damage to key resources that Mesolithic people used to survive.The simulated tsunami wave inundated a large area of resource-rich land beyond the usual tidal range throughout the Northumberland region.This includes the sites of Howick and Low Hauxley where inundation was at several metres depth.As outlined earlier, Mesolithic people lived a nomadic lifestyle focusing on these coastal and estuarine locations.Here, they utilized marine resources including fish, shellfish and seaweed which the erosive tsunami waves would have disrupted (Tanyaros and Crookall, 2011).Furthermore, the coastal area at Howick and Low Hauxley was probably forested, being dominated by hazel trees (Boomer et al., 2007;Waddington and Bonsall, 2016).These trees are often considered to have been heavily exploited by Mesolithic people for hazelnuts with 200 000 fragments being found at Howick (Cotton, 2007).Hazelnuts were a staple food source that could be stored through winter when other sources of food were scarce, although their high preservation potential can introduce some bias (Conneller, 2021).The Storegga tsunami would have produced similar impacts to the 2004 Indian Ocean tsunami where coastal forests were defoliated within the lowest 2.5 m (Hayasaka et al., 2009).This would have removed most of the hazelnuts which were collected during autumn (Bishop et al., 2014) when the tsunami occurred (Rydgren and Bondevik, 2015).Other seeds and nuts collected and stored during autumn would also have been lost.This would have increased the vulnerability of Mesolithic people to food shortages or even famine later in winter.Furthermore, some of this coastal forest would have been permanently lost with around 5% of coastal forests in Thailand in the 2004 tsunami not recovering for at least 5 years (Kamthonkiat et al., 2011).Overall, the tsunami would have caused short-term impacts on vegetation and resource availability, which would have extended through the winter months.
There may also be longer-term effects including lower resilience to the harsher winter months and effects of saline intrusion.The Storegga tsunami occurred during late autumn (Rydgren and Bondevik, 2015) just before winter when resources were hardest to obtain (Mithen and Wicks, 2018).Therefore, Mesolithic people would have been at their most vulnerable to the effects of the loss of resources from the coastal zone.Furthermore, saline intrusion would have caused longerterm damage to plants, with salinity in coastal croplands markedly increased following the 2004 (Kume et al., 2009) and 2011 (Roy et al., 2014) tsunamis.This salinity caused the majority of crops to wither and die in coastal areas after these tsunamis with damage still evident several years later (Nakaya et al., 2010).Whilst Mesolithic people had not yet developed farming, they would still have exploited wild grasses and plants that would have experienced similar effects.The removal of resources would have further lowered their resilience to survival through the winter months.These effects would have been amplified under the high tide model as saline intrusion would have reached farther inland, with inundation at Howick extended a few hundred metres through the river valley and up to 1 km farther inland in other areas of the domain.would have expanded the potential zone of lost resources, further lowering their resilience to food shortages.
The effects of the tsunami would have been unprecedented in living memory, as even the largest storm surges would have been comparatively insignificant.In Northumberland, the largest storm surge in the last 70 years produced an elevation of 1.2 m above the high tide level (Haigh et al., 2015).Such an event would be less likely to occur in the Mesolithic period than today due to recent climate change (Palmer et al., 2018).Even so, this is still significantly smaller than the 3.1 and 5.3 m above the high tide mark for the SLIP and tidal models of the tsunami at Howick.The decreased warning time of the tsunami would have further amplified the potential impact.Therefore, the unprecedented nature of the Storegga tsunami would have reduced the resilience of Mesolithic people to its effects.

Wider implications for Northern Britain
The inferred effects of the Storegga tsunami at Howick and Low Hauxley would probably have been replicated throughout Northern Britain.During the Mesolithic, ethnographic comparisons estimate the population of Northern Britain to be in the region of around 1000 people, so the population of a single site is small but not irrelevant on a wider scale (Mithen and Wicks, 2021).Moreover, settlement activity was predominantly coastal, with several other major sites in nearby southerneastern Scotland located within a few hundred meters of the coastline at East Barns, Castlandhill and Echline fields (Figure 2A; Mithen and Wicks, 2018).People living at these sites would have made use of similar coastal resources such as shellfish in the same manner as at Howick and Low Hauxley.Therefore, both the direct and indirect effects identified in this study would probably have been replicated at these other sites.Taken together these sites would represent a not insignificant proportion of the population at the time, with the combined effects likely to have had a major impact on the Mesolithic population of Northern Britain.Furthermore, recent work suggests that the tsunami also impacted the west coast of Scotland (Woodroffe et al., 2023), where several Mesolithic sites are located including Criet Dubh in an estuarine location (Mithen and Wicks, 2018).This suggests that the effects of the Storegga tsunami would have extended throughout the majority of the coastally orientated population of Northern Britain.This provides evidence to suggest that the tsunami was likely to have been the contributor to the inferred population decline in Northern Britain in the period after 8200 BP when the tsunami occurred (Wicks & Pirie, & Mithen, 2014;Waddington and Wicks, 2017;Mithen and Wicks, 2021).As outlined earlier, the model used to estimate this demographic collapse was imperfect and is difficult to apply at the site level.
However, the inferred effects from this study provide an evaluation of human-scale impacts of the event.Further work at some of these other sites would provide more evidence for this claim and increase confidence in the results of this study.

Conclusion
This study represents one of the first attempts to directly link the Storegga tsunami to its effects on Mesolithic people.Modelling results suggest that the sediment deposit from the Mesolithic site of Howick could have been formed by the tsunami, although only if the tsunami occurred at high tide.However, the tsunami provenance for the deposits is uncertain.The tidal regime would have been similar in Shetland when the initial wave impacted the islands, which could explain some of the underestimation of run-up height seen here.At Howick the tsunami's impact is shown to be severe with potentially catastrophic direct mortality as well as longer-term impacts on resource availability for survivors.Furthermore, a similar magnitude of impact at the nearby site of Low Hauxley suggests the effects were probably replicated throughout the coastally dominated sites of Northern Britain.Therefore, the tsunami was part of the population decline contemporary with the event.The framework of this study could be replicated at other sites or for different palaeotsunamis to better inform present-day coastal managers on the potential impacts of future high-magnitude tsunami events.
for the three different mortality estimates (CDMC, 2003;exponential and log-normal: Jonkman et al., 2008).The mortality risk is on a scale from 1 being a 100% chance of mortality from the incoming wave and 0 being no chance of mortality.

Figure 1 .
Figure 1.The locations of the Storegga slide and Storegga tsunami sediment deposits (red circles) with the run-up heights at each location indicated (produced in QGIS; based on Bateman et al., 2021).[Color figure can be viewed at wileyonlinelibrary.com]

Figure 2 .
Figure 2. (A) Map of the locations of Mesolithic sites that contained structures in Northern Britain (pit-houses and shelters) and the locations of sediment deposits from the Storegga tsunami (red).(B) The population of Northern Britain in the Mesolithic period, using the frequency of activity events as a proxy for population size.The timing of the Storegga tsunami is shown in purple along with two cold events of 8.2 ka (8200 BP) and 9.2 ka (9200 BP) in blue (reproduced after Mithen and Wicks, 2021).[Color figure can be viewed at wileyonlinelibrary.com]

Figure 3 .
Figure 3.The Mesolithic site of Howick with the positions of the Mesolithic site and later Bronze and Iron Age sites as well as the locations of the two sediment cores taken by Boomer et al. (2007).Satellite imagery: Aerial Digimap Service (2022).[Color figure can be viewed at wileyonlinelibrary.com]

Figure 4 .
Figure 4. (A) Overview of the SLIP model domain of Northumberland; the boundaries for the model are shown in black (SLIP 20-m contour) and red (forcing boundary for wave from Hill et al., 2014).(B) The coastlines in the different RSL reconstructions from the central section of the domain and the key sites mentioned in the text.(C) Zoom-in area of B, at Howick.(D) Close-up of the mesh depicting the variation in mesh resolution across the domain.(E) Zoom-in area of D at Howick, indicating where domain includes topography up to the 20-m contour.Topography/bathymetry is from SLIP model.[Color figure can be viewed at wileyonlinelibrary.com]

Figure 5 .
Figure5.The observed data and two differing functions that explain the relationship between mortality and water depth for areas with rapidly rising water, e.g. a tsunami.Mortality is on a scale from 1 meaning 100% chance of mortality from the incoming wave to 0 meaning there is no chance of mortality.After:Jonkman et al. (2008).

Figure 6 .
Figure 6.The simulated wave height of the SLIP model for the entire domain after each hour of the simulation.The topography of the land is shown where inundation has not occurred.Two hours in this simulation is ~8.5 h into the simulation of Hill et al. (2014).[Color figure can be viewed at wileyonlinelibrary.com]

©Figure 7 .Figure 8 .Figure 9 .Figure 10 .
Figure 7. (A-C) Wave height above mean palaeo mean sea level from the GIA (A), SLIP (B) and tidal (C) models at the location of the sediment core (Core A) taken by Boomer et al. (2007).The horizontal line indicates the position of relative sea level in each model, where the wave height is above the line inundation has occurred.(D) Speed of the wave through time at the sediment core location for each of the three models.[Color figure can be viewed at wileyonlinelibrary.com] ): * * d (8)