Tsunami generated by the Late Bronze Age (LBA) eruption of Thera were simulated using synthetic tide records produced for selected nearshore (∼20 m depths) sites of northern Crete, the Cyclades Islands, SW Turkey and Sicily. Inundation distances inland were also calculated along northern Crete. Modelling was performed by incorporating fully non-linear Boussinesq wave theory with two tsunamigenic mechanisms. The first involved the entry of pyroclastic flows into the sea, assuming a thick (55 m; 30 km3) flow entering the sea along the south coast of Thera in three different directions all directed towards northern Crete, then a thin pyroclastic flow (1 m; 1.2 km3) entering the sea along the north coast of Thera directed towards the Cyclades Islands. Flows were modelled as a solid block that slowly decelerates along a horizontal surface. The second mechanism assumed caldera collapse, of 19 km3 and 34 km3 modelled as a dynamic landslide producing rapid vertical displacements. Calculated nearshore wave amplitudes varied from a few metres to 28 m along northern Crete from pyroclastic flows, and up to 19 m from caldera collapse (34 km3 volume). Inundation distances on Crete were 250–450 m. Waves produced by pyroclastic flows were highly focused, however, as a function of sea entry direction. Smaller volume pyroclastic flows produced nearshore wave amplitudes up to 4 m in the Cyclades islands north of Thera. Wave amplitudes in the Cyclades from smaller volume caldera collapse (19 km3) were up to 24 m, whereas in SW Turkey were as low as 2.1 and 0.8 m (Didim and Fethye where LBA tsunami deposits have been found). Wave amplitudes for the larger volume caldera collapse (34 km3) were generally 2.5–3 times larger than those generated by the smaller volume collapse (19 km3). These results provide estimates for understanding possible consequences of tsunami impact in LBA coastal zones, thus providing criteria at archaeological sites for detecting inundation damage, as well as for contemporary hazard assessment; they also provide additional criteria for deciphering homogenite layers in the abyssal stratigraphy of the Ionian and eastern Mediterranean Seas.
The tsunami generated by the 17th century BC eruption of Thera (Santorini) volcano has been suggested as contributing to the demise of the Late Bronze Age (LBA) Minoan culture in the southern Aegean Sea (Fig. 1) region through widespread damage to coastal areas, particularly on Crete (Fig. 1) island (Marinatos 1939; Pichler & Schiering 1977; Yokoyama 1978, 1988; Pararas-Carayannis 1992; Driessen & Macdonald 1997, 2000; McCoy & Heiken 2000a, b). The Minoan culture was involved in extensive maritime activities, especially trade (Hood 1971; Doumas 1983; Dickinson 1994), thus being highly susceptible to damage from tsunami. Given the active Aegean tectonic setting and its consequent seismicity (McKenzie 1972; Papadopoulos et al. 1986; Papazachos 1990; Guidoboni et al. 1994; Meier et al. 2004) we might expect that this and other ancient cultures in the circum-Aegean and eastern Mediterranean areas to have experienced frequent tsunamis, much as modern cultures have (Galanopoulos 1960; Ambraseys; 1962; Papadopoulos & Chalkis 1984; Guidoboni & Comastri 2005; Papadopoulos et al. 2007, 2010; Papadopoulos 2009), and that recovery after a damaging event was readily accomplished.
However, many of the LBA coastal communities in the Aegean were abandoned and not rebuilt following the violent Plinian eruption of Thera. The Minoan culture went into decline and was supplanted within a few generations by the Mycenaean culture from the Greek mainland (Hood 1971; Barber 1987). It would seem, was exceptional in addition to other effects from the eruption. And that would imply that wave amplitudes, with consequent coastal inundation and run-up, also may have been exceptional contrary to studies that suggest otherwise (Dominey-Howes 2004; Pareschi et al. 2006).
That tsunami were generated during the explosive LBA eruption of Santorini seems clear from sedimentary deposits found on land (e.g. Didim and Fethye, SW Turkey and Gouves, north Crete, Minoura et al. 2000) and in the deep-sea (Cita 1997), in addition to observations from historic volcanic eruptions of lesser explositivity. Given the variables controlling onshore wave attributes, such as focusing effects during wave generation, bathymetry, nearshore and onshore geography, and more, wave characteristics along coastal areas may vary from negligible to dramatic. Accordingly, the lack of sedimentary deposits at any location need not discount tsunami generation, as suggested by Dominey-Howes (2004).
Contemporary observations of volcanism as a tsunamigenic mechanism allow a better understanding of how eruptive mechanisms, such as edifice collapses, pyroclastic surges and flows, debris flows and lahars, transfer energy to the ocean for initiating waves. New and refined computer programs can incorporate these observations to create simulations of wave propagation and parameters, and help predict nearshore modification of these waves. This is the objective of this paper with focus on tsunami generation during the LBA eruption, and on the impact of these waves in three directions from Thera: north to the Cyclades island complex which is situated in the near-field domain; south to the north coast of Crete; and east to the localities of Didim and Fethye, SW Turkey (Fig. 1). We use these locations as examples of possible impact to coastal areas during the LBA catastrophe by introducing two independent tsunami generation mechanisms: pyroclastic flows and caldera collapse.
There are additional implications in understanding tsunami related to that eruption. For example, predictions of which archaeological sites would have been most affected by wave inundation, and what level of damage might be inferred to ancient buildings and landscapes. Such information would provide focus for understanding and defining geoarchaeological evidence for the presence and magnitude of damage at sites—there are few geological or archaeological criteria at present for clearly distinguishing tsunami activity from that caused by other natural or anthropogenic activity. In addition, contemporary cultures require better tsunami-hazard assessment resulting from the rare mega-event, such as the LBA eruption, for hazard mitigation along coastal zones.
Prior Studies of LBA ‘Minoan’ Tsunami and Estimated Wave Attributes
Wave attributes based upon inferred tsunami deposits
Sedimentary deposits inferred to be from LBA tsunami inundation have been found (1) in western Turkey in coastal sediments (Minoura et al. 2000); (2) on Crete in archaeological sites at Amnissos (Marinatos 1939; Pichler & Schiering 1977), Pitsidia (Vallianou 1996) and Palaikastro (Bruins et al. 2008); (3) within the LBA volcanic succession on Thera at two locations (McCoy & Heiken 2000b); (4) in Israel near Tel Aviv (Pfannenstiel 1960) and other coastal areas (Neev et al. 1987) but also on the continental shelf off Caesarea Maritima, Israel (Goodman-Tchernov et al. 2009); as well as (5) in Cyprus (Astrom 1978; Meszaros 1978) and (6) possibly further in east Sicily (De Martini et al. 2010). These deposits are defined by layers containing rounded pumice from the LBA eruption, often mixed with marine shells and rounded gravels, and are presumed to represent deposition from tsunami rather than from storm waves that have washed in rafted pumice after the eruption (McCoy & Heiken 2000b). At the Amnissos archaeological site, Marinatos (1939) also attributed displaced building-stones to tsunami impact.
Dramatic evidence at the Palaikastro archaeological site comes from exposures in sea cliffs (McCoy & Papadopoulos 2001; Bruins et al. 2008). Here a damaged building with a distinct layer of broken building stone and pottery fragments extending approximately 40 m inland were observed. Lumps up to a few cm in size of LBA tephra, wetting having created partial cementation of the pozzolana, are incorporated within the matrix of the deposit, and thus the deposit is assumed to represent tsunami damage related to the LBA Thera eruption. Additional deposits are exposed in adjacent sea cliffs, where, as one example, a small ruined building is filled with rounded pebbles and cobbles. Elevations of deposits are 1–3 m above present sea level. Assuming subsidence of about 2 m here since the LBA, these structures were perhaps 3–5 m above sea level; their distance inland from the LBA coastline is estimated to have been about 20 m. Artefacts and pottery shards in these deposits remain unstudied and undated, nor have the buildings been excavated, and they are presumed to be of LBA.
These exposures are exceptional. Preservation of sedimentary deposits on land that can be attributed to large waves is poor because of erosion both during, and after, wave inundation (Dawson & Stewart 2007). It is in the exceptional depositional setting, such as an archaeological site or in quiet estuaries/lagoons, where such deposits might be preserved. This difficulty led Dominey-Howes (2004) to question whether damaging tsunami were, in fact, generated during the LBA eruption of Thera.
Additional evidence for tsunami generation by the LBA eruption has come from the deep ocean with the discovery of reworked seafloor sediments from over 50 cores in the eastern Mediterranean and Ionian Seas (Kastens & Cita 1981; Cita et al. 1984, 1996; Hieke 1984, 2000; Cita & Rimoldi 1997). These deposits, known as homogenites, can be up to 20 m thick, and were presumed to represent slope failures of unconsolidated muds into enclosed basins triggered by near-bottom hydrodynamic shear related to tsunami passage. That they are slump deposits that have collected and settled within these basins comes from their mixed-age micropalaeontological assemblages that equate to missing stratigraphic intervals on basin slopes. An LBA date for the slumping comes from the stratigraphic position of homogenites in the well-studied pelagic stratigraphy of abyssal sediments in the eastern Mediterranean and Ionian Seas—homogenites occupy the same stratigraphic position as the LBA ash layer, and occur only where the ash layer is absent. This stratigraphic relationship has been confirmed by carbon-14 dates (Hieke 1984, 2000; Troelstra 1987).
The notion of homogenites as tsunami deposits has been questioned by some authors (e.g. Pareschi et al. 2006), with respect to the magnitude of near-bottom hydraulic shear forces related to wave transit. Other authors, however, questioned if homogeneites were generated by the LBA tsunami or by tsunami caused later by large earthquakes in the Hellenic Arc and Trench (Papadopoulos 2009). That sedimentary deposits on land are the consequence of tsunami impact has been questioned from the perspective of their scarcity with respect to such a major event, as well as from modelling results, bringing into question if indeed tsunami were generated by the LBA eruption and if their geographic reach was limited to the Aegean Sea (Dominey-Howes 2004; Pareschi et al. 2006). In this paper we assume that these sedimentary deposits—on land and in the deep-sea—are the product of tsunami activity related to the LBA eruption, and that these tsunami indeed propagated throughout the eastern Mediterranean Sea basin.
Estimates of wave amplitudes along coastal areas in the Aegean have varied from the absurd to the conservative. Suggestions, for example, of waves almost as high as the Chrysler skyscraper building in New York City (319 m high) by Pelligrino (2001), or that waves deposited pumice to an elevation of 250 m on the island of Anaphi, east of Thera (Marinos & Melidonis 1971), essentially washing over the top of that island, have no basis in physical reality. Pumice deposits on Anaphi were subsequently shown to be of aeolian origin from an earlier eruption of Thera (Keller 1980). Criteria from deposits mapped in coastal zones that are presumed left by tsunami related to the LBA eruption suggest nearshore wave amplitudes on the order of 7–15 m along coastal Thera, SW Turkey and Crete (Pichler & Schiering 1977; Minoura et al. 2000; McCoy & Heiken 2000b; McCoy et al. 2000; Bruins et al. 2008).
Wave attributes from analytical estimates and modelling of tsunami
Yokoyama (1978) calculated possible wave heights at 200 m water depths around Thera at 50–86 m. Kastens & Cita (1981) calculated open-ocean amplitudes of approximately 2–17 m, applying Yokoyama's estimates from around Thera, in their study of homogenite deposits in the eastern Mediterranean and Ionian Seas. Minoura et al. (2000) and McCoy et al. (2000) estimated wave heights up to 17 m along portions of the northern coastline of Crete from computer modelling, supported by additional modelling of wave parameters off eastern Crete (Bruins et al. 2008). Modelling of tsunami generated during the LBA eruption by Pareschi et al. (2006) contradicts these results and implies tsunami activity constrained only to the southern Aegean Sea with minimal wave activity beyond this area into the eastern Mediterranean Sea.
The LBA Eruption and Tsunami Generation
The LBA eruption of Santorini was enormous, estimated by Dunn & McCoy (2002) and McCoy & Dunn (2004) at a Volcanic Explositivity Index (VEI) of 7+, and certainly the largest known eruption in this region during human antiquity. Like most Plinian eruptions, the LBA event followed a sequence somewhat characteristic of such large eruptions (Bond & Sparks 1976; Heiken & McCoy 1984, 1990; Druitt et al. 1989; Druitt & Sparks 1996; McCoy & Heiken 2000a). The activity started with a precursory phase of magmatic and phreatomagmatic eruptions that left a thin 6 cm-layer of ash on southern Thera (Fig. 2). Following this, within months, came the first major Plinian phase that deposited a thick (up to 9 m) layer of pumice via fallout from a central plume. Partial evacuation of the magma chamber with first phase activity led to the initiation of caldera collapse that allowed entry of seawater into the vent—the eruption went into a second phase of highly explosive activity generating numerous thin pyroclastic surges and flows. Complete collapse of the caldera resulted in third phase activity and the generation of one or two, massively thick (up to 55 m), pyroclastic flows; this was accompanied by collapse of two peripheral grabens, to the north and west, opening wide connections between the caldera and the Aegean Sea. Final caldera collapse and degassing of the magma chamber resulted in additional pyroclastic flows during the fourth phase of activity. Rainfall during this final stage resulted in significant remobilization of tephra-laden slopes into debris flows and lahars.
We infer tsunami generation during the LBA eruption from observations of tsunamigenic mechanisms documented during historic volcanic eruptions: (1) lava bench collapse, (2) seismic activity, (3) submarine eruptions, (4) phreatomagmatic eruptions, (5) pyroclastic surges and flows, (6) caldera collapse, (7) landslides, (8) lahars and debris flows and (9) airwaves produced during explosions (Latter 1981; Self & Rampino 1981; Beget 2000). During the LBA eruption, all of the latter eight mechanisms likely produced sea-surface disturbances to trigger tsunami in all directions from the island. Pyroclastic flows occurring in the second through fourth volcanic phases and debris flows/lahars occurring in the third and fourth phases entering the sea would have generated numerous wave sets in all directions, including north towards the Cyclades island complex, south towards Crete and east towards SW Turkey. Caldera collapse and resulting tsunami would have been directed through the collapsed grabens towards the north and west/southwest in the second through fourth phases, the latter passage accounting for the homogenite deposits found on the abyssal seafloor of the eastern Mediterranean and Ionian Seas.
For this study, we simplify these multiple generative mechanisms into two tsunamigenic triggers. For the first mechanism we treat pyroclastic, debris flows, and lahars as a single type of triggering event upon their entry into the sea. A second triggering mechanism considers caldera collapse to create tsunami. We focus on, and model, waves propagated into the northern, eastern and southern Aegean Sea.
Modelling Tsunami Produced by Pyroclastic Flows, Debris Avalanches/Flows, Mudflows, Lahars and Landslides
Tsunamigenic mechanisms implied from historic observations
Pyroclastic flows result from gravitational collapse of eruption plumes that produce hot gas-rich flows of volcanic particles that sweep across the ground at high speeds (e.g. Druitt 1998). The main body of a pyroclastic flow consists of a massive, clast-rich, poorly-sorted assemblage of pumaceous pyroclasts, lithic fragments and ash, and thus may be considered analogous to a debris flow (Druitt 1998; Watts & Waythomas 2003).
Tsunamis can be produced by three mechanisms related to the entry of pyroclastic flows into the sea: (i) phreatomagmatic explosions upon encountering seawater, (ii) physical impact of the flow with the ocean and (iii) plume shear mechanisms. Of these, the first is the most energetic as a source mechanism and is several orders of magnitude greater in efficiently generating waves, thus has the greatest potential for tsunami generation (Sigurdsson et al. 1991a; Nomanbhoy & Satake 1995; Freundt & Bursik 1998; Dufek et al. 2007). Upon encountering the ocean, a pyroclastic flow may ingest water and transform partially or completely into a debris flow, and continue moving along the seafloor as a subaqueous debris flow; the denser portion of the flow may also detach from the less dense portion and flow over the surface of the ocean for long distances (Sigurdsson et al. 1991b; De Lange et al. 2001; Watts & Waythomas 2003). That pyroclastic flows can generate tsunami by entry into the sea was adequately demonstrated during historic eruptions, for example: Awu, 1856, 1892 (Indonesia, tsunami of few metres high; Caruso & Pareschi 1993); Ruang, 1871 (Indonesia, tsunami about 25 m high; Neumann Van Padang 1959); Krakatau, 1883 (Indonesia, tsunami up to 31.4 high; Simpkin 1969; Rampino & Self 1982); Montserrat, 1997, 2003 (Leeward islands, Caribbean Sea, tsunami up to 4 m high; Calder et al. 1997; Pelinovsky et al. 2004).
Production of tsunami by debris flows/avalanches, mudflows, lahars and landslides are also well documented during and after volcanic activity. Water displacement upon rapid entry of these flows and slides produce tsunami as consequence of impact much like pyroclastic flows, and if still hot might also contribute to tsunami production through phreatomagmatic explosions. Historical examples include: Maruyama, 1792 (Japan; landslide with tsunami up to 30 m high; Beget 2000); Krakatau, 1883 (landslides accompanying caldera collapse with tsunami up to 40 m high; Latter 1982); Augustine (Alaska, 1883, debris flow with tsunami up to 8–9 m high, Kienle et al. 1987, and 1964 landslide with tsunami up to 9 m high, Kienle et al. 1987; Beget & Kienle 1992; Beget & Kowalik 2006); Soufriere Hills, 2003 (Montserrat, debris avalanche with tsunami up to 4 m high; Herd et al. 2005,); Stromboli, 2000 (Italy, landslide with tsunami up to 10 m high; Tinti et al. 2000). In the volcanic complex of Thera, strong extra-caldera activity in the submarine volcanic bank Columbo on 1650 October 10 (New Style) triggered a large tsunami up to 50 m high very possibly because of collapse of the submarine volcanic cone (Dominey-Howes et al. 2004).
Assumptions and techniques for modelling tsunami due to pyroclastic flows
Numerical models of tsunami generally use a source mechanism that produces an elevated nearshore dome of water, such as the one that would result from a pyroclastic flow of large volume entering the sea (Nomanbhoy & Satake 1995; De Lange et al. 2001). In this study we used the software package GEOWAVE (Watts et al. 2003), which is a combination of TOPICS (Tsunami Open and Progressive Initial Conditions System) and FUNWAVE. TOPICS uses a variety of curve fitting techniques and was designed (Grilli & Watts 1999) as an approximate simulation tool that provides surface elevations and water velocities as initial conditions for tsunami propagation. The numerical model FUNWAVE (Wei & Kirby 1995; Wei et al. 1995; Kennedy et al. 2000; Kirby et al. 1998, 2003) performs wave propagation simulation, based on the fully non-linear Boussinesq theory using a predictor–corrector scheme.
GEOWAVE allows simulating a wide range of wavelengths and is capable also to model the physics of breaking waves and resulting run-up (Ioualalen et al. 2006, 2007; Grilli et al. 2007). In addition, GEOWAVE accurately captures a wide range of wave parameters and characteristics, such as non-linear shoaling, refraction, diffraction, and wave-wave interactions.
The important parameters of pyroclastic flows that control the interaction between flows and water body appears to be bulk density, velocity, discharge rate of the flow as well as the angle of incidence between the flow and the water surface (Cas & Wright 1991; De Lange et al. 2001). For example, the 1883 Krakatau eruption involved low bulk density pyroclastic flows with high velocities and a high angle of incidence into the sea. In this regard, we assume that the LBA eruption of Thera was similar to that of Krakatau, and that steep slopes both within and surrounding an internal caldera characterized the pre-eruption LBA physiography of Thera (Heiken & McCoy 1984; McCoy 2003), thus representing a similar interaction between pyroclastic flows and water.
In terms of dynamic water surface displacement, De Lange et al. (2001) noted that their dynamic model had a non-linear relationship between flow volume and the maximum wave height observed at any location. Thus, flow volume needs to be well defined. For example, increasing the displacement volume from 1 to 10 km3 increases the wave heights by a factor of 10 for an initial static water surface displacement. However, under dynamic water surface displacement the same increase in volume produces even larger increase in wave amplitudes. Single flows display a strong directivity (focusing) in the resulting wave height distribution, so that the largest waves occur along the linear axis of flow direction. The directivity effect, however, weakens when multiple flows with a large radial spread occur (De Lange et al. 2001). A small eruption column tends to produce single flows, but as the magma discharge rate increases it is possible to produce multiple flows simultaneously. High discharge rates, such as those estimated for the LBA Thera eruption (Sigurdsson et al. 1990) are also associated with high flow velocities (Walker 1979).
For calculating pyroclastic flows as tsunami source mechanisms, following the suggestions of Watts & Waythomas (2003) and Walder et al. (2003), we use a pyroclastic flow (Fig. 3a) with bulk density ρb with physical dimensions of width W, thickness T and length L, entering water of density ρw and constant depth. The flow is represented as a solid block that slowly decelerates as it travels along a horizontal surface with a given initial velocity. The differential equation governing the centre of mass motion along x axis with its origin at the shoreline is approximated byWatts & Waythomas (2003), which is based on experimental results, we choose the values of dynamic coefficients as follows: Cm= 1.0, Cd= 1.0, Cn= 0.05. The characteristic distance, time and velocity of motion are defined as
When a pyroclastic or debris flow enters water it generates a single, coherent wave profile (Walder et al. 2003). Deformation of the pyroclastic/debris flow during transformation and subaqueous translation is assumed to have minimal effect on tsunami generation and as a result it can be neglected (Watts & Grilli 2003). Only the submerged portion of the pyroclastic flow is considered in our calculation, assuming that a density separation of the flow would occur at a short distance from the shoreline if not within the splash zone. However, the tsunami is produced by the submerged portion which is assumed the largest portion of the initial flow.
Parameters of the largest pyroclastic flow for the LBA eruption (Table 1) were adopted by taking into account two attributes. The first concerns the kinematic parameters of the pyroclastic flow (Table 1) which were adopted from the analytical techniques described by Watts & Waythomas (2003). The second attribute used comparative studies from other large caldera forming eruptions, such the Tambora 1815 and Krakatau 1883 for estimating the volume of the pyroclastic flow. At Tambora, according to Scarth (1994) 40–80 km3 of tephra was expelled as nuées ardentes (pyroclastic flows) and ashflows. At Krakatau, the total bulk volume of pyroclastic deposits, including co-ignibrite ash, was estimated to be 18–21 km3 (Self & Rampino 1981). These estimates were compared (Table 1) with the volume of pyroclastic material estimated for the Thera eruption. Pyroclastic/debris flows entering the sea during the LBA eruption along the southern coastline of Thera, capable to produce tsunami, were developed during the second through fourth stages of eruption but mainly during the third stage (McCoy & Heiken 2000a,b). To estimate potential maximum and minimum wave amplitudes, we generated a tsunami by a pyroclastic flow 55 m thick and 5 km wide at the shore during the third eruption phase, and 1 m thick and 1.6 km wide (average thickness from field data) during the second and fourth phases (Heiken & McCoy 1984). The most recently mapped marine pyroclastic deposits around Thera volcano were estimated on the order of 60 km3 (DRE-dry rock equivalent; Sigurdsson et al. 2006). Although our estimates lead to even larger volume of pyroclastic deposits, up to 110 km3, in our final simulation we considered tsunami generation assuming pyroclastic/debris flow with a conservative volume of 30 km3. Three different options for the pyroclastic flow directivity (focussing) were selected: 120o, 145o and 200o measured from north counter-clockwise. Details on the computational domain are explained in the next section.
Modelling Tsunami Produced by Caldera Collapse
Mechanisms implied from historic observations
Production of tsunami by caldera collapse is well documented in prehistoric and historic eruptions such as in Kikai, Japan, where 7.3 ka BP a tsunami of about 20 m was generated following volcano collapse which formed caldera 19 km in diameter (Maeno et al. 2006); in Krakatau, Indonesia, 1883, tsunami up to 40 m was generated following volcano collapse forming a caldera of 7 km in diameter (Latter 1981; Self & Rampino 1981; Sigurdsson et al. 1991a; Carey et al. 2001); in Ritter Island, Papua-New Guinea, 1888, tsunami up to 15 m was generated following volcano collapse forming a caldera of 2.5 km in diameter (Johnson 1987). Sudden and progressive caldera at free fall usually is of 20 s in duration (Maeno et al. 2006).
During the LBA eruption of Santorini, caldera collapse occurred during the second through fourth phases of the eruption, apparently with rapid ingress of the sea following collapse of two peripherals connecting the central caldera to the sea during the third eruption phase (Heiken & McCoy 1984, 2001a; Druitt et al. 1992, 1996; Sakellariou et al. 2006). A pre-existing caldera was expanded to a one of 8 km × 9 km in dimension, with maximum contemporary water depths of 380 m (Perissoratis 1995); collapse volume is estimated at 19–34 km3 (Heiken & McCoy 1984, 2001a).
Summarizing the historical examples described above we can conclude that the sudden and progressive caldera collapse at near free-fall creates high potential energy for tsunami generation. Because of consequent enormous sea level changes, the resulting wave forms one of the most hazardous tsunamis certainly in the near-field.
Assumptions and techniques used for modelling tsunami due to caldera collapse
We have adopted the simulation model proposed by Walder et al. (2003) for tsunami generation by subaerial mass flows, and consider the case of massive slide released at a height, H, above water level. The slide accelerates down a steep slope or falls directly into the sea, decelerating only when it reaches the seafloor (Fig. 3b). For the caldera collapse model, the initial water heights, H, correspond to the depth of the caldera formed. We again apply TOPICS (Watts & Waythomas 2003) to estimate surface elevations and water velocities as initial conditions for the tsunami propagation model performed further with FUNWAVE. To simulate the LBA Thera tsunami we made near-field estimations of tsunami wave amplitudes produced by caldera collapse using an initial minimum collapse volume of 19 km3 (case 1) and then we repeated calculation for maximum volume of 34 km3 (case 2).
Description of the computational domain
Bathymetric 1 arcmin data set from the General Bathymetric Chart of the Oceans (GEBCO) were used for the simulation grid for eastern Mediterranean and southern Aegean Sea regions. In coastal regions of Cyclades islands as well as around Crete, however, supplementary bathymetric data from the Hydrographic Service, Greek Navy, on a grid of 0.5 arcmin, were utilized. For the Aegean Sea computations, based on GEBCO spatial domain we reconstructed a grid with 500 m uniform spacing with water depth measured from mean lower water level. For the caldera collapse mechanism, additional constraints concerning bathymetric variability in and around Thera island varies from 380 m in the source of tsunami generation to 100–400 m in the area surrounding Thera.
Tsunami produced by pyroclastic/debris flows entering the sea
Wave amplitudes offshore, mainly at shallow water depths of about 20 m, were calculated at the locations of twelve synthetic tide-gauges along the northern coast of Crete (Figs 1 and 4). The coordinates and water depths below synthetic tide-gauges are listed in Table 2. These locations include the archaeological sites at Gouves (location 7), as well as two sites near ancient Amnissos (location 6) and Palaikastro (location 12) where tsunami deposits were identified by Minoura et al. (2000), McCoy & Papadopoulos (2001) and Bruins et al. (2008). Our simulations indicate that tsunami generated by pyroclastic/debris flows entering the ocean are highly focusing along the axis of flow direction, as noted also by De Lange et al. (2001).
Maximum wave amplitudes on northern Crete from sea entry along the south coast of Thera of 55-m-thick pyroclastic/debris flows are shown in Fig. 4. The first tsunami wave arrival in north Crete is about 13 min after its generation in Thera. Coloured pictures at the left-hand side of Fig. 4 portray the distribution of tsunami wave amplitudes offshore of the northern coast of Crete. The right-hand side of Fig. 4 shows the tsunami waveforms at the offshore locations of northern coast of Crete (tide gauges 1–12) as well as at the near field islands of Cyclades (tide gauges 13–15), two locations of SW Turkey (tide gauges 16 and 17). In particular, for locations, where tsunami deposits were found such as Gouves (7), Palaikastro (12), Didim (16) and Fethye (17) at the directions of their maximum impact (200o), wave amplitudes range from –19.2 to 19.9 m; from –11.3 to 8.8 m; from –2.49 to 2.6 m and from –1.19 to 1.54 m. Fig. 5 (right-hand site) demonstrates tsunami waveforms at the locality of Augusta Bay (18), where tsunami deposits were found as well (De Martini et al. 2010). The maximum impact of this wave is at 120o with amplitudes ranging from –0.6 to 1.2 m.
Maximum amplitude of waves generated by thin (∼1 m) pyroclastic /debris flows entering the sea along the north side of Thera are on the order of 4 m at the islands of Ios, Folegandros and Sikinos located to the north of the volcanic source (Fig. 6b). We would expect similar near-field amplitudes wherever such flows entered the ocean, which occurred to all directions peripheral to Thera as suggested by McCoy & Heiken (2000b). Such amplitudes are consistent with estimates from exposures of sedimentary deposits (Fig. 6a) presumed left by tsunami impact on Thera during the fourth phase of eruptive activity (McCoy & Heiken 2000b).
Estimates of maximum inundation in coastal areas along northern Crete vary from 250 m at location 3, to 370 m at location 6 and to 450 m at location 8. As for the calculation of inundation and vertical run-up GEOWAVE assumes a very rough shoreline and onshore topographic condition. However, GEOWAVE does not explicitly take into account the presence of obstacles such as trees, houses, boulders, etc., but it runs as if such roughness is present. This is expressed with an 8% land porosity that diminishes onshore wave motion. Finally, it should be noted that there is significant sensitivity to grid resolution for wave run-up elevations computed at a particular site (Titov & Sinolakis 1997; Ioualalen et al. 2007). Detailed run-up simulations need a fine grid, with accurate bathymetric and topographic data. Thus, the inundation estimations here as well as estimates of run-up that can be easily obtained from offshore wave heights, using formulas such as those proposed by Chesley & Ward (2006) and Ward & Day (2008), should be considered as approximate.
Comparison of our results for tsunami generation by pyroclastic flows with that of Bruins et al. (2008) at Palaikastro, where authors suggested 9 m as a minimum run up and about 300 m inland, based upon tsunami deposits found there, led us to conclude that our mechanism of tsunami generation by pyroclastic flow was constrained correctly. The values of the initial wave height in the generation region at Santorini are also in good correlation: 25 m is of Bruins et al. (2008) and 26 m is of our modelling. We should emphasize, however, that 55-m thick pyroclastic flow can be represented by two flows, as earlier suggested by Heiken & McCoy (1984). This would allow a better correlation between the size of tsunami waves at the locations and the tsunami deposits found there.
Tsunami produced by caldera collapse
Calculations for offshore tsunami wave heights were made for five locations. Three of the locations are situated in islands surrounding Thera (see Fig. 1): the southern coasts of Amorgos Island (location 13) and Ios Island (location 14) both situated directly to the north and northeast of Thera; and location 15 in Naxos Island which is situated to the northwest of Thera. The other two locations are Didim (location 16) and Fethye (location 17) in south-western coast of Turkey where traces of tephra of the Thera eruption and tsunami deposits were found and studied by Minoura et al. (2000).
Results of simulation for a volume of collapsed material of 19 km3 are shown in Fig. 7. At the southern coasts of Amorgos (location 13, Fig. 7a) amplitude of about 11 m is calculated. At the other two locations of Cyclades, that is in Ios and Naxos islands (locations 14 and 15, Figs 7b and c) the wave would have amplitude of 24 and 19 m, respectively. Perhaps it was the propagation of such tsunami waves to the west and southwest, into the eastern Mediterranean and Ionian Seas that triggered deposition of the homogenite layers found in the abyssal pelagic stratigraphy.
In this modelled scenario, the tsunami wave reached the locations of Didim and Fethye, SW Turkey, about 1.2 and 1.6 hr after the caldera collapse with offshore amplitude equal to 2.1 and 0.8 m, respectively (Fig. 7d). Based on the above results we calculated that increasing the volume of the collapsed material to 34 km3, the tsunami wave height in Cyclades increases up to 50 m. These results indicate that the very complicated process of caldera collapse concludes in a series of chaotic waves with very high amplitudes. Fig. 8 illustrates results obtained for the northern coast of Crete (upper panels) and locality of Augusta Bay, Sicily (lower panels) from the simulated model of tsunami wave produced by sudden caldera collapse with total volume of either 19 or 34 km3. The comparison of wave profiles received from synthetic tide-gauges shows that in the case of caldera collapse involving the highest limit of the material volume (34 km3), the wave amplitude will be 2.5–3 times larger than that for material volume of 19 km3 involved in the collapse.
Figs 7 and 8 suggest that the case of low-volume caldera collapse (19 km3), can explain tsunami wave amplitudes at Gouves (7) and Fethiye (17), consistent with tsuanmi deposits found at these two archaeological sites. We note again that wave amplitudes at Palaikastro, however, need to invoke a high-volume caldera-collapse mechanism.
Minoura et al. (2000, 2003) and McCoy et al. (2000) estimated that along portions of the northern coastline of Crete maximum wave height from caldera collapse ranged between 6 and 11 m, which is approximately that of our model, which yields offshore amplitudes ranging between 4 m and 18 m depending on the volume of the material involved in the caldera collapse. Pareshi et al. (2006) modelled extreme scenarios for the LBA eruption and estimated amplitudes from 10.6 to 14.0 m for the locations of Amnissos and Gouves, north coast of Crete. However, they calculated run-ups of only a few metres in the same locations from triggers of smaller size.
We have considered pyroclastic flows and caldera collapse as two dominant tsunamigenic mechanisms, characteristic of the large Plinian-type, caldera-forming processes during the LBA Thera eruption. We then modelled tsunamigenic consequences for selected coastal sites in the southern Aegean region. In this study we combine the sea entry of both pyroclastic flows and debris flows for modelling purposes, and consider two extreme cases during the LBA eruption derived from field mapping on Thera: massive flows (∼55 m thick) and thin flows (∼1 m in thickness). The caldera collapse mechanism is modelled as a series of large landslides that produce rapid vertical displacement, assuming caldera collapse volumes of 19 and 34 km3. To better understand coastal impacts, applying these two tsunamigenic mechanisms we first modelled tsunami generated by pyroclastic flows entering the sea along the southern coast of Thera to constrain a range of values for tsunami amplitudes along the northern coasts of Crete. Further we considered tsunami generated by caldera collapse with waves propagating to the west and north, the former as a possible source mechanism for the generation of homogenite layers recovered in abyssal sediments of the eastern Mediterranean and Ionian Seas, and the latter as an estimate of wave impact upon islands north of Thera.
Along the north coast of Crete, calculated wave amplitudes were highly variable (Fig. 4), ranging from negligible to 28 m, due to the entry of pyroclastic flows of varying thickness entering the ocean along the southern coast of Thera. We found significant focusing of wave energy where such flows have produced tsunami. Energy focussing combined with nearshore bathymetric variability may account for the wide range in calculated wave amplitudes. Inundation distances on land of between 250 and 450 m were calculated assuming an 8% land porosity.
Our calculated wave amplitudes for tsunami along the north coast of Crete generated during the LBA eruption of Thera are in agreement with estimates by others based upon modelling and/or field criteria (Pichler & Schiering 1977; McCoy et al. 2000; Minoura et al. 2000; McCoy & Papadopoulos 2001; Bruins et al. 2008).
As discussed above our results for both scenarios contradict partly those of Pareschi et al. (2006) who suggested runups of a few tens of metres on the northern coasts of Crete. Rather we show wave heights to have been highly variable along this coastline, from negligible to large, thereby accounting for the discrepancy with the conclusions from these two publications. Discrepancies between conclusions presented with those of Pareschi et al. (2006) we attribute to differences in assumed physical parameters determining tsunamigenic factors during the eruption that were applied to computer modelling. For example, Pareschi et al. (2006) modelled the tsunami produced by caldera collapse, using a simplified framework that includes the total volume of collapsed material which is equivalent to the total volume of water displaced. Their assumption was that the initial negative water displacement is instantaneous. In the present study, however, we employed an approach that considers dynamics of the massive slide masses.
Modelled tsunami wave amplitudes for nearshore areas along the southern coast of Ios due to caldera collapse are on the order of 24 m; it may be assumed that it is the propagation of these waves to the west and southwest, into the eastern Mediterranean and Ionian Seas, that triggered deposition of the homogenite layers found in the abyssal pelagic stratigraphy.
These results are in good agreement with those of Bruins et al. (2008), suggesting that tsunami waves at Palaikastro had coastal amplitudes of about 9 m high. We also emphasize that for the northern coast of Crete, such a wave was produced by pyroclastic flow entering the sea at azimuth of 200o. Very likely the caldera collapse source could only be responsible for tsunami waves propagating north and northwest from Thera.
Finally, it is hoped that by providing numerical model of tsunami generation during the LBA eruption of Thera we might assist others to search for residues of tsunami impact in sedimentary deposits or ancient cultural debris in the southern Aegean region. We demonstrate that the production of tsunami is not only entirely likely, but could have resulted in significant wave amplitudes and run-up on coastal areas north and south of Thera, based upon physical models employing criteria from observations of tsunamigenic processes during historic volcanic eruptions. These conclusions have significance for interpreting the consequences of tsunami impact on the ancient, Late Bronze Age, cultures that inhabited these coastal areas.
This work was done while the senior author (T.N.) had a Research Associateship at the Malcolm H. Wiener Laboratory at the American School of Classical Studies at Athens (September 2004-January 2005), and the third author (F.W.M.) was the Visiting Malcolm H. Wiener Professor at the same laboratory (2007–2008). T.N. and F.W.M. acknowledge the financial and scholarly assistance offered by this laboratory and its director Dr. S. Fox. Laboratory facilities at the University of Hawaii are also acknowledged. Portions of this research were done while F.W.M. was also a Senior Fulbright Scholar to Greece (2007–2008) and an Alexander Onassis Scholar (2008). G.A.P.'s contribution was funded by the EC research project TRANSFER, contract n. 037058, FP6–2005-Global-4, Reduction of tsunami risks. The authors are thankful to Dr. Ph. Watts, Dr. Ah. Yalciner and Dr. G. Heiken for their kindness in providing helpful suggestions as regards tsunami modelling the first two and volcanic processes the third. Authors are grateful to two anonymous reviewers, whose thoughtful comments improved substantially the initial version of manuscript.
Additional Supporting Information may be found in the online version of this article:
Movie S1. This animation corresponds to the case,when 55-m-thick pyroclastic flowentered the sea from the south coast of Thera with an azimuth of 145° and the produced tsunami wave is approaching Crete island. To get an idea of the values of the tsunami wave amplitude along the coast of Crete islandwe refer readers to Fig. 4(b) of the manuscript.
Please note: Wiley-Blackwell are not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the article.