Computer simulations of tsunamis due to sector collapse at Stromboli, Italy

Journal of Volcanology and Geothermal Research 96 Ž2000. 103–128 www.elsevier.comrlocaterjvolgeores

Computer simulations of tsunamis due to sector collapse at Stromboli, Italy Stefano Tinti a,) , Elisabetta Bortolucci a , Claudia Romagnoli b b

a Dipartimento di Fisica, Settore di Geofisica, UniÕersita` di Bologna, Viale Carlo Berti Pichat 8, 40127 Bologna, Italy Dipartimento di Scienze della Terra e Geologico-Ambientali, UniÕersita` di Bologna, Piazza di Porta S. Donato 1, 40127 Bologna, Italy

Received 15 September 1998; received in revised form 30 August 1999; accepted 10 September 1999

Abstract Stromboli is an island volcano of the Aeolian Volcanic Arc, characterised by persistent activity. The cone rises about 2500–3000 m from its submarine base with very steep slopes; its summit at 924 m above the sea level. The subaerial growth of Stromboli, occurred in the last 100 ka, has been marked by repeated episodes of large gravitational collapses especially affecting the NW flank of the island in the last 13 ka. The last one occurred less than 5000 years ago forming the deep depression on the NW seaward flank, named Sciara del Fuoco ŽSdF., and it produced very likely large water waves. This paper envisages a scenario where a huge mass of volcanic material collapses into the sea in the same sector in which the SdF collapse took place in Holocenic times and it computes possible tsunami evolutions assuming the present-day bathymetry. Numerical simulations are performed by means of two distinct models; one for the mass collapse and one for the tsunami. Slope failure dynamics are calculated with the aid of a Lagrangian model: the landslide is subdivided into blocks, and the motion of each constituent block is calculated by applying the basic principle of mechanical momentum conservation, with block–block and block–ambient interactions being taken into account. Water waves are computed by solving a system of shallow-water equations including a forcing term dependent on the sliding mass motion. The finite-element ŽFE. technique is employed since it permits the use of non-uniform grids, which are adequate to account for marine basins with irregular coastlines. In addition to the sensitivity analysis concerning the main parameters governing the slide motion, two main cases are explored, differing in the slide path followed by the mass. The resulting tsunami is very large, with giant waves as high as several meters Žtens of meters in the worst cases. impinging the coast. Due to the strong wave refraction induced by bathymetry, waves travel around the island, affecting even the island coast opposite the source. q 2000 Elsevier Science B.V. All rights reserved. Keywords: debris avalanche; numerical model; sector collapse; Stromboli; tsunami

1. Introduction Stromboli is a persistently active island volcano of the Aeolian Volcanic Arc, which is located on the )

Corresponding author. Tel.: q0039-51-209-5025; fax: q0039-51-209-5058.

Sicilian–Calabrian northern continental slope in the southeastern Tyrrhenian Sea. The geologic setting is well studied and described in the literature, and is only briefly outlined here. Stromboli rises from water depths ranging between 1200–1500 m, on the east and south flanks, and 1700–2200 m, on the north and west flanks Žsee

0377-0273r00r$ - see front matter q 2000 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 7 - 0 2 7 3 Ž 9 9 . 0 0 1 3 8 - 9

104

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

bathymetry sketched in Fig. 1.. With its summit at 924 m above sea level, it rises up to 3000 m towards the northward-deepening SE Tyrrhenian slope and bathyal plain. The volcano has broadly bilateral symmetry about a NE-trending axis, along which major vents dikes and eruptive fissures are preferentially located, and which reflects the influence of regional tectonism ŽRosi, 1980; Gabbianelli et al., 1993.. From its elongated submarine setting, and from the distribution of marine wave-cut terraces truncating the submerged summit of older NE and SW portions of the volcanic edifice, it appears that the subaerial cone developed along older NE–SW-aligned centres ŽRomagnoli, 1990; Gabbianelli et al., 1993.. The island formed mainly within the last 100 ka and its stratigraphy and petrography are well known ŽRosi, 1980; Francalanci et al., 1989; Gillot and Keller,

Fig. 1. Sketch of the Stromboli island and of the surrounding areas. The island is dark grey with isohypses superimposed Žequidistances 200 m.. Contour lines of the bathymetry are also plotted with labels given in meters. The collapse-scar area Žpale grey., is partly subaerial ŽSdF. and partly submarine. It begins from the cone summit and goes down to a depth of about 700 m, involving the NW flank of the volcano. The points on the island coast and in the sea marked with small solid circles are places for which time histories of the water waves excited by the landslide are computed and graphed.

1993; Hornig-Kjarsgaard et al., 1993.. The subaerial evolution of Stromboli comprises two main stages of construction, followed by major flank and sector collapse events ŽRosi, 1980; Pasquare` et al., 1993.. The related volcanic products occur on opposite sides of the island and are in contact with major discordances along a NE–SW-trending surface. Deposits of the later constructional stages Žlast 13 ka, Gillot and Keller, 1993; Hornig-Kjarsgaard et al., 1993. are mainly confined to the NW side of the island. A NNW-ward migration of the eruptive activity occurred during the evolution of Stromboli, concurrent with the series of large-scale gravitational collapses, which were probably facilitated by tectonically induced zones of weakness such as NE–SWtrending fissures, mainly with NW-dipping planes ŽRosi, 1980; Zanchi and Francalanci, 1989; Pasquare` et al., 1993.. Moreover several interacting factors, including constructional loading and oversteepening, lack of basinward support for the slope along unbuttressed NW and SE flanks, differential stress induced by tectonically-controlled intrusion along the main NE–SW axial zone, seismic and tectonic activities, play a role in promoting the gravitational instability of Stromboli ŽKokelaar and Romagnoli, 1995.. The most recent sector collapse, which took place less than 5000 years ago, generated the Sciara del Fuoco ŽSdF., a steep depression marking the NW seaward flank of Stromboli. It is laterally limited by two steep walls up to 300 m high, which represent the subaerial part of the scar that formed after the collapse. Lavas and pyroclasts produced after the SdF collapse prograde across, and partially fill, the slide surface. The slope of the partial cone, constructed between the lateral walls of the collapse scar, acts as a channelway to the sea for most of the recent and present-day eruptive products. SdF is presently the site of persistent volcanism of Stromboli and most loose material below the vents move gravitationally down its steep Žup to 388. slope by creeping, avalanching and similar processes ŽBarberi et al., 1993; Kokelaar and Romagnoli, 1995.. 2. Morphologic and volumetric constraints of the Sciara del Fuoco collapse An oceanographic survey carried out in 1991 aboard the ship ‘‘Minerva’’ of the Italian National

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

Research Council ŽCNR. obtained the morphological features of the area offshore SdF, determined down to 2200 m below sea level Žbsl. by means of echosounding and high-resolution single-channel seismic profiling, integrated with observation of seabed samples ŽRomagnoli et al., 1993; Kokelaar and Romagnoli, 1995.. From this, it was inferred that the SdF is just the subaerial half-portion of a partially filled sector-collapse scar that extends down to 700 m bsl. The horse shoe shape and dimensions of the whole scar Ž3 km long and up to 2 km wide; 1.5 km of vertical distance between the lowest and the highest points. are consistent with it being the source of a debris avalanche ŽKokelaar and Romagnoli, 1995.. These are characterised by the fast movement of large volumes of rock and debris fragments over great distances also on low angle slopes ŽVoight and Elsworth, 1997. and the related deposits show greater mobility than slumps. From the 1991 and most recent surveys Ž1994, 1995 and 1998. a major fan-shaped feature has been recognised in front of the SdF scar from 700 m bsl ŽKokelaar and Romagnoli, 1995.. It extends far from the island and beyond the submerged area considered in this simulation, having been mapped down to over 2500 m bsl and more than 18 km far from the Stromboli coast. The morphological characters, surface and internal appearance of the deposit Žas revealed by deep-tow Side-Scan Sonar images and high-resolution seismic profiles. agree with a primary debris-avalanche origin, probably accompanied by extended disintegration of the collapse material. Its volume has been estimated, on the area investigated in 1991, as being larger than 4 km3 ŽKokelaar and Romagnoli, 1995.; more precise calculations are presently in progress, considering the farther extended part of the collapse-related material downslope. The morphologic and volumetric constraints used in the present work for the SdF collapse have been obtained by digitising and computing 3D reconstructions of both the emerged and submerged portions of the scar from the 1991 data set ŽKokelaar and Romagnoli, 1995.. Estimates of the original volume of material removed in the collapse have been obtained by assuming an original conical form of the pre-collapse Stromboli edifice. Two possible projected slide surfaces A and B, with U-shaped cross-sections and

105

a smoothly curved long-section profiles, have been hypothesised for the slide surface Žsee Fig. 11 of Kokelaar and Romagnoli, 1995.. Profile A has a gradually decreasing axial slope, downwards from 408 to 208, and the maximum thickness of collapsed material is about 250 m. Profile B has a steep back wall, but it becomes nearly horizontal at the base Žabout 58., and corresponds to a maximum thickness of 700 m for the collapsed material. Profile B should give the likely maximum permissible volume of rock involved in the collapse. The convergence of the submerged walls of the collapse-scar suggests that the lowest level of the slide plane can have been no deeper than 700–750 m bsl. Both hypothesised surfaces A and B have been digitised, and for each, a 100 = 100 m matrix was obtained by computerised bidimensional interpolation of points. Then the volume enclosed between each hypothetical slide surface and an original conical surface was calculated. Despite the uncertainty of locating the slide surface beneath the partial fill of the Sciara scar, the likely volume of material that was involved in the SdF collapse is in the range of 0.97 to 1.81 km3 , respectively obtained for surfaces A and B ŽKokelaar and Romagnoli, 1995.. The large disparity between the volume of the SdF collapse-scar and the volume of the surveyed debris-avalanche deposit should not be a cause of surprise, considering that the latter comprises deposits from more than one large-scale sector collapses that affected the NW flank of Stromboli also before the SdF event, involving volumes of the order of 1–2 km3 each one for a total volume estimated about 3 to 6 km3.

3. Simulation of the landslide Avalanches and collapses are not infrequent at volcanoes, and may involve large masses with consequent disastrous effects. Theoretical modelling of slope failures and mass flow is influenced by our present lack of full understanding of the basic mechanics of the very complex processes involved. Questions concerning the way energy is dissipated, the interaction of the mass with its boundaries, the mass rheology, and how momentum is transported in the flow, are still unanswered ŽPierson, 1995.. De-

106

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

spite this, computer models to study slope failure dynamics have been developed and applied by a number of researchers in the last 20 years Žsee the review papers by Takahasi, 1981, and by Iverson, 1997., to simulate laboratory experiments Žcf. Savage and Hutter, 1991; Hutter et al., 1995. as well as real events, especially so for those relatively few recent cases where experimental field data are available Žsee e.g. the applications to study the 1980 Mount St. Helens, and the 1984 Mount Ontake avalanches ŽMcEwen and Malin, 1989; Voight and Sousa, 1994; Sousa and Voight, 1995... The numerical model used here to simulate the landslide at Stromboli was developed by Tinti et al. Ž1997. with the purpose of providing a numerical tool suitable to a large variety of mass sliding, from rock avalanches to flows of granular material, and with the motivation of investigating processes of tsunami generation by submarine landslides. The model, which was tested against laboratory data Žsee also Tinti et al., 1998., is based on a Lagrangian view, in the sense that the sliding mass is partitioned into a finite number of blocks and then the evolution of each block is calculated by considering all mechanical forces acting on it, viz. gravity, ambient pressure, resistance forces exerted by the surroundings on the block surfaces, and the interaction forces due to the other blocks. All blocks are seen as deformable mechanical bodies that cannot change their volume during their motion. Mathematical details are given in Appendix A, and only some general considerations are made here. It is remarked that this approach can be interpreted as an extension of the centre-of-mass ŽCoM. models, which was once extensively applied to study the motion of landslides and snow avalanches. CoM models are essentially 1D models, where the landslide is viewed as a single block moving on a prescribed curve Žcalled the flow line., and the main interest is focused on the motion of its barycentre Žsee Korner, 1976; Perla et al., ¨ 1980; Hutchinson, 1986; McEwen and Malin, 1989.. There are two main differences between our model and CoM. First, considering a multiblock mass instead of a single block implies the need to model interaction forces between adjacent blocks. In our model, blocks are aligned along a chain that is assumed to be oriented in the direction of sliding. Block k interacts with the preceding block k y 1

and with the following block k q 1. The interaction is modelled as a sort of a collision process able to modify the interacting block velocities, by preserving however the total momentum of the system. The interaction law is described in terms of interaction coefficients with values depending upon the characteristics of the material, and upon the instantaneous distance of the barycentres of the interacting blocks. This means that their values depend on the slide dynamics, since they take different values in different parts of the slide and at different times. The second relevant difference concerns the resistance forces due to the ambient fluid. These are very important for masses moving in water and able to generate water waves, as in the case of interest in this paper, although they are unimportant in the case of subaerial slides, where they can be neglected. Since these forces depend on the extent of the lateral and upper surfaces of the slide, our model, different from CoM models, cannot be 1D, but it is required to compute the surface of the sliding blocks at any time. Therefore, besides the flow line on which block barycentres are supposed to move, even the basal surface swept by the blocks during their motion must be specified here. In the following, the flow line and the basal surface will be denoted as slide path and as slide surface respectively. The landslide we have simulated for the W flank of Stromboli has the mass corresponding to the hypothesis A of Kokelaar and Romagnoli Ž1995., with a volume slightly smaller than 1 km3. The scar surface is shown in plan in Fig. 1, and contoured in Fig. 2, where also the assumed top surface of the slide body is plotted. The slide path and the slide surface relative to the landslide motion must be provided as input data to the model. Since they are determined by the bathymetry and the slide dynamics, in principle they are not known a priori, which seems to be inconsistent with the model requirement. In practice, however, in many instances they can be given a viable estimate from field observations and the requirement of a priori knowledge is not so great. In our case, since the simulations refer to a predictive tsunami scenario considered in an actual frame, although inspired to the Holocenic sector collapse, the present-day bathymetric setting offshore SdF has been taken into account ŽFigs. 1 and 3.. This shows the presence of an active canyon, which, starting

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

107

canyon, that is, the eastern edge of the collapse deposits. The wider surface ŽFig. 10b. has been enlarged further to the west, considering the location of the debris-avalanche deposits attributed to the SdF and previous collapse events and the fact that overbank currents and other unconfined gravity-driven flows appear to have widely covered them with further deposits ŽKokelaar and Romagnoli, 1995.. The slide paths corresponding to the two hypotheses are depicted in Figs. 10a and 10b and do not differ too much from one another. Both slide surfaces are arbitrarily interrupted around a depth of 2000 m, since the most important tsunami generation occurs where the mass moves in shallow water, as we will see later. In the following, the narrower and the wider surfaces will be more conveniently referred to as surfaces I and II respectively. On surface

Fig. 2. Contour plot Ždashed lines. of the slide surface assumed in this paper for the collapsing mass Žhypothesis A of Kokelaar and Romagnoli, 1995.. The assumed top surface for the mass is also contoured Žsolid lines.. Slide and top surfaces are known over a regular grid of squared cells with 100-m long sides. The solid triangles show boundary nodes of the mass in the discrete grid.

from the lower limit of the collapse scar at 700 m bsl, and being gradually deflected northward, marks the eastern edge of the collapse deposits downward to its confluence with the Stromboli Canyon. The canyon is presently the major pathway for volcanoclastic density currents generated at or near the SdF shoreline during eruptions, storms and following slope-failure events, that are then mainly constrained to follow the canyon axis downslope ŽKokelaar and Romagnoli, 1995.. Two different hypotheses have been explored. The two possible slide surfaces, shown in Figs. 10a and 10b, embrace the same scar surfaces of the Sciara collapse; both are centred on the active canyon which extends from the SdF scar and have the same eastern boundary, which is topographically constrained by the NW submerged base of the volcanic edifice. On the contrary, the two surfaces differ for their western boundary, which for the narrower surface ŽFig. 10a. corresponds to the left side of the

Fig. 3. Computed shaded-relief model of Stromboli and its submerged NW flank in front of SdF with the slide surface used in the simulation Žfrom q800 m asl to y800 m bsl. and the present-day bathymetry Žfrom y800 m to y2200 m..

108

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

I, numerous computer runs have been performed with the purpose of testing the sensitivity of the code to the parameters governing the motion of the slide. The main reason resides in the total lack of experimental data in the field, with the consequent difficulty of selecting specific and reliable values for the model parameters. To make the illustration of the numerical results easier, it is convenient to present a specific reference case first, for which the sliding mass has been partitioned in a chain of 10 blocks that have similar volumes. The parameters used in the computer runs are listed in Table 1, and the results of the computations are shown in Fig. 4. The velocities of the barycentres of the blocks and the heights of the blocks are plotted respectively in the upper and in the lower panels of the figure. As expected, velocities increase quite rapidly in the first 50 s after slide initiation and then decrease regularly and slowly, peak velocities attaining values around 60–65 mrs. Block velocities differ from one another, especially at the beginning, since blocks in the front of the landslide are much slower. The reason is that they slide on a slope that is less steep than that run by the blocks at higher positions and, moreover, being underwater, they are affected by the resistance exerted by the sea water since the beginning of their motion. This causes the landslide to become shorter in the direction of sliding. Around 30 s and later, the first block in the landslide front Žblock no. 10. becomes the fastest, while all others advance with Table 1 Parameters used in the computer simulations for the landslide model C1a s 0 C2a s 0 C1w s 0.008 C2w s 2 pi k s 0 ri k s 0 g s 0.5 ma s tan da s 0.09 m w s tan d w s 0.02 l s 0.01 r s 2500 kgrm3 r X s1000 kgrm3 s s 0.2 s 1 s 0.8 max D t s 0.5 s

Fig. 4. Case I: narrower slide surface. Barycentre velocities of the 10 blocks into which the landslide has been subdivided Žtop panels. and block heights Žbottom panels.. Blocks have been labelled from the rear Ž1. to the front Ž10. of the slide body.

similar velocities. The explanation is that this block is decelerated less than the blocks following it, according to the water resistance term, since the height of its frontal surface is the smallest. The consequence is that the slide tends to become longer. About 160 s after the slide initiation, the mass reaches the end of the assumed slide surface, which limits the domain used for computations, and then the run stops, though the slide still possesses a considerable velocity Žover 40 mrs.. The lower panels show that the heights of the blocks change considerably during the motion. Their initial values range between 30 and 150 m, with the highest blocks being in the central part of the landslide. The frontal block tends to become progressively less high, while all others have more complex histories including local thickness maxima. It is remarkable that the largest values for the height of the central blocks Žnos. 6 and 7. around 50 s exceed 260 m, which means a thickness increase larger than 70% over their initial value. The case corresponding to the mass sliding on the wider surface ŽII. leads to similar results. All qualitative considerations made previously for block veloci-

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

ties and heights hold even for this case. It is worthwhile noting, however, that block heights are smaller. This can be clarified by considering that block k in case I and the corresponding block k in case II have the same size Ži.e., same length, width and height. at the initial time when the landslide lies completely on the scar surface Žobserve on Figs. 10a and 10b that surface I is identical to surface II in this zone.. But, during the motion, blocks sliding on surface II tend to become wider than on surface I, because they must fit a broader basal surface, and accordingly they have smaller heights.

4. Sensitivity analysis and discussion The case with slide surface II has been taken into account to explore the dependence of the solution upon some physical parameters of the model. These can be grossly divided in three classes: discretisation parameters concerning the algorithm of transforming the sliding mass into a discrete set of blocks, resistance parameters regarding the interaction of the mass with the ambient, namely ocean bottom and ambient fluid ŽEqs. 8a–9 in Appendix A., and internal interaction parameters involved in the laws describing the block–block interaction ŽEqs. 10a–12c, in Appendix A.. The sensitivity analysis has been focused on the parameters of the second class, that is, the bottom friction coefficient m x and the water resistance coefficients C1 x and C2 x , chiefly because sensitivity of the model to the parameters of the other classes was already investigated in previous papers ŽTinti et al., 1998, 1999.. For the sake of completeness, before coming to the results of the sensitivity test it is however convenient to make some considerations concerning also the main parameters not explicitly investigated in this paper. The most relevant parameters of the first class are the number of blocks N into which the mass is sliced and the characteristic inter-node distance Ž D x . of the grid used to discretise the initial bottom and top surfaces of the mass, as well as the basal slide surfaces I and II. It can be shown that for a given D x, convergence of the solution is reached quite quickly as N grows larger ŽTinti et al., 1999.. In our case, the initial geometry of the mass was available on a regular square-cell grid with cell step D x s 100

109

m, and substantial convergence has been attained as soon as N s 10, which entails a mean block length of about 450–550 m with smallest lengths of about 150 m, and which makes it useless to use larger numbers of blocks. The third class of parameters is involved in the interaction between contiguous blocks and governs the strength of the internal forces Žfirst of all the interaction coefficient e i k , the intensity coefficient l, the shape coefficient g , and the deformability coefficients s , s 1 .: according to their values, the mass can behave as a body with no deformation allowed along the path, or, on the other extreme, as a largely deformable body with possible great changes of its axial length Žsee Tinti et al., 1997, 1998 for a discussion.. In this study we have adopted the reference values listed in Table 1 and already used to simulate the recent tsunamigenic landslide event of Vulcano in 1988 ŽTinti et al., 1999.. In our case, this set of parameters allows substantial variations of the mass length Žup to about 150%. that are smaller than, though in the same order as, the changes of the mass width due to the widening of the slide surface. Practically, it is assumed that the slide moves as a deformable mass rather than flowing as a turbidity current. It should be stressed that our model is not adequate in principle to describe the dynamics of turbidity currents that may originate from the possible incorporation of the ambient fluid into the sliding mass and its consequent extreme dilution. Here we assume that the mass is partitioned in constantvolume blocks, however deformable they can be. This limitation however is not important in our study, in view of the very huge volume of volcanic rocks involved in the collapse ŽKokelaar and Romagnoli, 1995., and to the very short interval of time after the motion initiation that is of interest here in the view of tsunami generation, as will be pointed out in later sections. The idea is that, though processes of mass fluidification cannot be ruled out, they may have involved at most a small marginal fraction of rocks during the tsunamigenic phase of slide motion, without turning it into a full turbidity flow. It is, however, pointed out that fluidification of the sliding mass, which was the prevalent traditional explanation for the observed large runout of several submarine slides, might not be the only possible explication. Alternative mechanisms, such as hy-

110

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

droplaning, have been recently recognised as being equally effective. This phenomenon, which has been observed in laboratory experiments for the front of subaqueous flows of incoherent granular materials ŽMohrig et al., 1998., can be easily accounted for by our model by adequately reducing the bottom friction coefficients. Furthermore, it should be born in mind that fluidification is a very difficult process to model properly in the framework of continuum mechanics, which is the viewpoint adopted here. However, it is even more complicated in the frame of the approach of discrete particle theory, which was not frequently used as an alternative tool to study the motion of incoherent material. Indeed, the application of this latter view gave rise to computer simulations where the slide is considered as formed by a very large number Žin the order of 10 6 . of discrete interacting particles ŽHopkins, 1987; Campbell, 1990., and has been used so far only for subaerial mass slides, where the interaction with the ambient fluid may be neglected. These models can certainly account for separation of grains Žwith grain–grain collision causing grains to bounce and leap at the mass surface.. The grain dynamics however is often influenced by the assumed shape and size of the particles, which introduces a strong element of arbitrariness. Moreover, the number of the involved equations is generally so high to challenge the present state of computer power and capacity and to permit the computation of the solution at most over vertical cross-sections, involving one horizontal and one vertical coordinate ŽCampbell et al., 1995.. The latter is a limitation that cannot be accepted in the present study, where the slide is the source of the tsunami, and its evolution on the 2D ocean bottom surface is of maximum interest. Finally, we remind that the model has been successfully tested against laboratory experiments performed by Hutter et al. Ž1995. concerning subaerial sliding masses formed by cohesionless grains ŽTinti et al., 1997, 1998.. Thus, we are confident that this model, provided that the mass friction with the ocean bottom is reduced suitably and a very high degree of strain for the constituent blocks is allowed can produce acceptable results even for underwater events involving incoherent materials, such as debris avalanches, though under the restrictions that the mass motion does not turn into pure turbidity flow.

Let us now turn the attention to the results of the sensitivity analysis, starting with the bottom friction coefficients. The range of values of m x that can be found in the literature is very large and they are derived mostly from laboratory experiments ŽHutter et al., 1995; Watts, 1998. or from geological observations of runout and subsequent data inversion ŽLabazuy, 1991; Hungr, 1995; Iverson, 1997; Voight and Elsworth, 1997.. Distinction must be made between the friction produced at the basal surface outside Ž ma . and inside Ž m w . the sea water, since underwater friction is expected to be less effective due to possible mass hydroplanning, but mostly to the minor capacity of submarine loose sediments to produce large resistance stresses. In his study of volcanic landslides at Piton de La Fournaise at the French oceanic island named La Reunion, Labazuy ´ Ž1991. explored values of ma in the range 0.03–0.25 and values of m w in the range 0.08–0.15. In our reference case presented above, the values taken by ma and m w are 0.09 and 0.02, respectively, which correspond to rather weak bottom friction. We have then studied three more cases, one at lowest extreme where ma s m w s 0.02 with friction at a very low level, one at the opposite side with ma s m w s 0.25 with high bottom friction, and an intermediate case with ma s m w s 0.09. All such cases have been run with the parameter values listed in Table 1, with the exception of the coefficients C1w and C2w , which assume the values 0.004 and 1.0, respectively. Results are illustrated by means of Figs. 5 and 6, where the velocity of the centre of mass of the slide vs. time, and the instantaneous profile of the mass at various times are respectively diagrammed. Bearing in mind that slide simulations are interrupted as soon as the slide front reaches the end of the slide surface, we understand why the velocity curves of Fig. 5 terminate at different times. Only the slide affected by the highest bottom coefficients Ž ma s m w s 0.25. is progressively slowed down until it stops within the computation domain. All curves have the same qualitative trend, with a quicker increase and a slower descent. Expectedly, cases with smaller friction correspond to higher slide velocities. It is remarkable that the difference between the times, in which the mass peak velocities are attained, is large Žranges from about 50 to 90 s.. Peak velocities vary even more sensibly between ca. 30 and 80 mrs. Looking

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

111

Fig. 5. Velocity of the centre of mass of the slide, computed as the weighed average of the velocities of the slide blocks Ž Õi k ., vs. time for different sets of bottom friction coefficients. The mass slides on surface II.

are at most slightly different at any times. They share the same value of m w s 0.02, whereas ma changes from 0.02 to 0.09. This is suggestive that the effect of ma in the simulations is rather weak and that, correspondingly, the friction parameters that matters most is only m w . This is quite understandable in our case since the subaerial path of the slide is remarkably smaller, and consequently, less influential, than the underwater path. The second set of experiments regard the resistance coefficients C1 and C2 , strongly dependent upon the properties of the ambient fluid. Their values are known to be much smaller for subaerial than for submarine slides. In this paper we have neglected the effect of these resistance terms for the motion outside the sea Ž C1a s C2a s 0. and have considered a total of 5 different sets of values for the submarine coefficients C1w and C2w , which are assumed to range respectively from 0 to 1.5 and from 0.004 to 2.

at the slide profiles, which are illustrated by plotting the instantaneous height of the blocks vs. the position of their centres of mass, it may be appreciated that the initial profile exhibits one absolute Žblock no. 7. and one relative Žblock no. 3. height maximum, which tend to be maintained during the motion. An exception to these is the highest friction case, where the relative maximum converts progressively to a change of profile slope. The difference between the various cases reflects the difference already observed for the velocity graphs. At a given time, high-speed slides precede slower slides with separation increasing as time elapses and with block thickness decreasing, owing to the progressive widening of the slide surface at more distal regions. The difference in the position of the masses in the last shown snapshot ŽT s 80 s. between the extreme cases is very large, separation being almost complete. It can be however remarked that slide profiles are shown plotted vs. time in Fig. 6 to magnify the discrepancies. If they were seen as function of the slide position rather than vs. time, the differences among the various cases would appear not so great. The final important observation that can be made is that the two cases with the smallest friction coefficients give very similar results: velocity curves are almost indistinguishable in Fig. 5 and slide profiles

Fig. 6. Profiles of the mass slide at different times after the slide initiation. The height h i k of block k is assigned to the position of the centre of mass of the block j i k and is plotted with a symbol. The profile is obtained by joining symbols together. Each panel corresponds to a different time and displays four curves relative to the four examined sets of friction coefficients. The slide moves on basal surface II.

112

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

These intervals contain the values that can be found in the literature Žsee e.g., Slingerland and Voight, 1979 and Watts, 1998.. In the course of these numerical runs, all other parameters have the values given in Table 1. The basal surface considered is again slide surface II. The results are illustrated in the same form as the those of the previous sensitivity tests, that is, by means of two Figures: Fig. 7, concerning the velocity of the slide barycentre, and Fig. 8, showing profiles of the slide at various times. The diagram with the curves of the slide speed vs. time enables us to rank straightforwardly the five used sets of the resistance coefficient. The weakest resistance is provided by the couple C1 s C2 s 0.004, while the largest one is produced by the set C1 s C2 s 1.5 Žnotice that in this discussion the coefficients are simply designated by C1 and C2 , with the subscript ‘‘w’’ dropped.. None of the slides comes to rest within the computational domain used for the slide motion, reaching the more distant border at different times according to their speed. Peak velocities differ very much: the highest value, corresponding to the set C1 s C2 s 0.004 is unrealistically large, being in the order of, or being larger than, the value estimated for subaerial slides in volcanic environment ŽVoight et al., 1981; Naranjo and Francis, 1987; Labazuy, 1991.. Profile panels shown in Fig. 8

Fig. 7. Velocity of the centre of mass of the slide, computed as the weighed average of the velocities of the slide blocks Ž Õi k ., vs. time for different sets of water upper Ž C1 . and frontal Ž C2 . resistance coefficients. The slide occurs on surface II.

Fig. 8. Profiles of the mass slide at different times after the slide initiation. Curves are plotted as explained in the caption of Fig. 6. Each panel contains five curves corresponding to the five sets of water resistance coefficients but the last one, due to the fact that at t s100 s, the fastest slide corresponding to the set C1 sC2 s 0.004 has already reached the distal border of the slide surface and the computations stopped. The slide moves on basal surface II.

are similar to those displayed in the corresponding Fig. 6. An observation that can be made here is that separation of slides relative to different cases occurs at later times. Notice that all curves are practically still superimposed at 10 s. The reason is that, water resistance grows quadratically with slide speed Žsee Eq. 9 in Appendix A., and at the beginning of the motion it is very small, irrespective of the values of the resistance coefficients. While the results of the numerical runs we have carried out in the previous section showed that using two different slide surfaces I and II leads to quite similar slide evolutions, the sensitivity analysis performed here on the friction and resistance coefficients of the model proves that they can affect greatly the resulting slide speed. On the other hand, it seems that other slide characteristics such as slide length and profile Žif seen as functions of the posi-

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

tion rather than time. do not change significantly. Practically it can be stated that the slide motion is quite similar in all runs but occurs on different time scales.

5. Simulation of the tsunami Tsunami generation by landslides or rockfalls is well known and the process is being given increasing attention since failure of large masses can produce giant sea waves that can devastate coastal areas, as has been once more demonstrated by the last large tsunami that took place in July 1998 in Papua New Guinea. Here a tsunami, very likely generated by a submarine slump triggered by an earthquake, devastated the nearby coasts causing total destruction and thousands of victims ŽDavies, 1998; Kawata et al., 1999; Tappin et al., 1999.. One other good example comes from Italy, where the tsunami with the highest toll of victims Žover 1500. in the known Italian history was caused by the collapse of a mountain sector into the sea triggered by a strong earthquake sequence. This occurred in February 1783 in Calabria, southern Italy, where most inhabitants of the small village of Scilla were swept away by the waves produced by avalanches that detached from Mount Pacı. ´ The event was very minutely reported in numerous coeval documents and in texts and monographs by several scholars Žsee the recent Italian catalogue of tsunamis by Tinti and Maramai, 1996.. Numerical methods to study landslide-triggered tsunamis are generally based on various forms of the shallow-water approximation to the Navier–Stokes equations of fluid mechanics Žsee e.g., Aida, 1975; Kienle et al., 1987; Ma et al., 1991; Harbitz, 1992; Harbitz et al., 1993; Jiang and LeBlond, 1992, 1993; Jones and Mader, 1996; Imamura and Gica, 1996; Johnsgard and Pedersen, 1996; Tappin et al., 1999. that are used both in the tsunami generation as well as in the tsunami propagation phase. These are generally considered adequate, though lack of field data makes it problematic to carry out sound validation tests. One of the few examples of a satisfactory field data set is the known case of the tsunamigenic Skagway Harbor landslide, which occurred in November 1994 in Alaska ŽKulikov et al., 1996.. It is remarkable that recent simulations performed by

113

Fine et al. Ž1999. succeeded in reproducing the available experimental data, such as the record of the tide-gauge placed close to the source in the harbour and eyewitness accounts. They used a model, based on shallow-water approximation both for the slide, modelled as a viscous fluid, and for the sea waves. It is a modified version of a previous model introduced by Jiang and LeBlond Ž1994., who claimed erroneously and misleadingly Žsee their paper title. that it was a 3D model, although it is essentially 2D, being based on vertically integrated equations for the slide and on purely classically shallow-water equations for the tsunami. The computer simulations run here are based on the shallow-water model solving the following equations: Et z s Et h y = P Ž h q z . u Et u s yg=z y Ž u P =u . y Cf u Ž h q z .

Ž 1. y1

< u<

Ž 2.

where the z Ž x, y,t . and uŽ x, y,t . are respectively the unknown sea surface elevation and the unknown horizontal velocity vector at the point Ž x, y . and at the time t, and hŽ x, y . is the basin depth with respect to the sea surface at rest. Furthermore, g is the vertical component of the gravity acceleration, and Cf is the bottom friction coefficient that is taken to be equal 0.005 in the numerical experiments. Numerical simulations deal necessarily with a limited domain and complementary boundary conditions are needed to complete the above system of equations. Different conditions are imposed on coastal and on open sea boundaries. While open boundaries are assumed to permit free propagation of waves and to be perfectly transparent, coastal boundaries are taken to cause full wave reflection, which mathematically is expressed by: u P n s 2 Ž c1 y c 0 . , u P n s 0,

coastal boundaries

c1 s Ž g Ž h q z . . c 0 s Ž gh .

open boundaries

1r2

1r2

Ž 3a . Ž 3b . Ž 3c . Ž 3d .

where c1 and c 0 are the local phase velocities for non-linear and linear waves respectively.

114

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

Notice that in shallow-water approximation the vertical fluid velocity is assumed to be negligible and the horizontal velocity u is independent from the vertical coordinate, since it is the average taken from the sea bottom to the surface. The idea that in the process of generating water waves by slides vertical velocities may be neglected sounds unnatural and odd, and deserves some comments. Traditionally, hydraulic engineers dealing with interaction between flow currents and bodies and with waves produced by underwater body movements have shown preference to use models accounting for vertical velocities. These models are usually 2D Žone vertical plus one horizontal dimension, which is aligned in the direction of the slide motion. and cover a variety of cases from pure inviscid potential theory to more complex equations accounting for nonlinearity, wave dispersion, etc. ŽNoda, 1970, 1971; Heinrich, 1991, 1992; Assier Rzadkiewicz et al., 1997; Watts, 1998.. The conditions of applicability of shallow-water approximation to underwater body motion were first explored by Hammack Ž1973. in an important pioneer work where theoretical analysis was supported by hydraulic tank experiments. Considering the case of the vertical displacement of a sector of an ideally horizontal sea floor, he was able to prove that shallow-water theory and potential theory provide the same results if the following main conditions are satisfied: Ž1. the size of the displaced sector is large compared to the water depth, and Ž2. the bed motion is impulsive, i.e., occurs over a time scale much shorter than the typical time scale of water waves free propagation. These conditions were recognised to be very likely fulfilled by dislocations produced by submarine earthquakes, and sometimes by underwater landslides. Practically, the vertical water velocities that are accounted for by potential models correspond to quick Žinstantaneous. vertical mass movements produced by the sea bed dislocation, and have the ultimate effect of displacing the sea surface approximately by the same amount as the sea bed. The consequence is that the kinetic energy associated with the vertical water motion is converted to potential energy available for the evolution of the water waves, that can be in turn satisfactorily described by shallow-water theory. Extension of Hammack’s work has been attempted by following studies such as Iwasaki’s Ž1997., where suitable

corrections to the shallow-water theory were proposed to cover small-scale quick underwater landslides, which is a very common case. This point will be elucidated better in the following. Here it is worth recalling that shallow-water approximation has proven to be accurate even against laboratory tests concerning wave generation by moving walls and submerged wedges ŽVilleneuve and Savage, 1993., which is a further explanation of the reason why it is so widely adopted in modelling waves excitation by landslides. Of course there are extreme cases where it cannot be applied, such as the classical box-drop problem treated by Noda Ž1970. and Wiegel Ž1970., that is, the laboratory simulation of vertical rockfalls into the ocean, but these cases cannot hamper its applicability to wave-exciting masses sliding along underwater inclines, which is the main topic of this section. In the set of shallow-water equations ŽEqs. 1 and 2., the term Et h represents the forcing due to the underwater landslide and in the present model it is

Fig. 9. Grid with triangular elements used in the FE simulations of tsunami waves. It consists of 2906 nodes and 5526 triangles.

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

related to the landslide height, say h lŽ x, y,t ., by the following expression: Et h s f Et h 1

Ž 4a . y1

f s w cosh R x Ž 4b . where R is the ratio of the local water depth hŽ x, y . to the horizontal characteristic scale of the slide, depending on the slide dynamics. It was first introduced in modelling water waves generated by tectonic dislocations as a term dependent on the progressive rupture of the seismic fault by Hwang and Divoky Ž1970. who studied the 1964 Alaskan tsunami. They used the instantaneous vertical displacement of the ocean floor as h lŽ x, y,t . and assumed the coefficient f s 1. The reducing factor f is introduced here in view of the aforementioned correction needed by shallow-water theory to deal with small-case disturbances. Notice that when the slide is large compared to the water depth, the ratio R nearly

115

equals 0, and accordingly the reducing factor f is very close to unity, with the consequence that the exciting term is almost identical to the rate of change of the slide height. On the other hand, for very small-scale slides, the factor f is quite large and no appreciable forcing can take place. The reducing factor has a clear physical meaning: disturbances on the sea floor produce effects on the sea surface but with amplitude depending on the thickness of the water mass, that is, the sea depth, since they tend to decay with the distance from the ocean bed Žsee Kajiura, 1963; Iwasaki, 1997.. Therefore f plays the role of a transfer function, which accounts for vertical water movements, and serves to introduce the necessary physical information on typically 3D processes into the 2D shallow-water model. The finite-element ŽFE. method has been used to solve the above governing equations, by means of a code expressly developed to deal with tsunami simu-

Fig. 10. Ža. Slide surface for case I, showing the assumed seaward continuation of the slide area. Its eastern and western limits are plotted with crosses, while the intermediate line is the assumed slide path. Groups of the tsunami grid points aligned on four different cross-sections of the slide surface, are also given. Labels of isobaths are in meters. Žb. The same as Ža. but for case II. Notice that slide surfaces I and II have the same east boundary, but surface I is significantly narrower.

116

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

lations, and has been extensively applied to study historical as well as recent tsunamis Žsee Tinti and Piatanesi, 1996; Piatanesi et al., 1996; Piatanesi and Tinti, 1998a; Tinti et al., 1999., and to study theoretical key cases ŽPiatanesi and Tinti, 1998b.. The technique is based on a double-step time discretisation, accurate to the second order in the time step size, D t, and on a space discretisation with the domain being covered by a mesh of triangular elements, and with the discrete equations being obtained through the weighed residuals method and through linear shape functions ŽTinti and Gavagni, 1994, 1995; Tinti et al., 1994.. The grid for the numerical simulations was constructed by using an original mesh-building algorithm ŽTinti and Bortolucci, 1999a., and is shown in Fig. 9. It is a mesh of unequal triangles, with triangle size being balanced against a local scale that is proportional to the local phase velocity of the wave c 0 . Since c 0 is an increasing function of water depth, this implies that triangles in deep waters are larger and triangles in shallow coastal waters are smaller, thus providing better grid resolution to describe coastal details. Notice that the possibility of using grids with different local resolution is a natural feature and at the same time a great advantage of the FE methods over finite-difference methods. The latter use regular generally square-cell grids and to change resolution are bound to weld together different meshes with different cell size. The grid is formed by 2906 nodes and 5526 triangles, with a resolution on the Stromboli coast ranging 50–100 m. The forcing term Et h has been evaluated by means of Eq. 4 on the tsunami grid nodes by applying suitable space–time interpolation to convert its values from the landslide computational mesh to the tsunami mesh of Fig. 9. To understand better the characteristics of forcing, it is convenient to consider some typical plots of Et h vs. time. For both reference simulation cases corresponding to parameters given in Table 1, and to slides riding on slide surfaces I and II, we have selected some points within the slide surface located on four distinct cross-sections: three grid points for each cross-section gathered in four distinct groups, with groups ordered according to increasing distance from the island Žsee Figs. 10a and 10b.. Curves corresponding to the two cases have the same general trend, and

thus only curves referring to case I, namely to the narrower slide surface, are illustrated in Fig. 11. Considering these diagrams, we see that: Ž1. points on the same cross-section experience very similar forcing; Ž2. forcing amplitude decreases with distance from the island, being smaller when the slide reaches deep waters, according to the effect of the reducing factor f ; Ž3. forcing initial time increases along with the distance from the coast; Ž4. forcing duration is an increasing function of the distance; and Ž5. forcing is initially positive, in correspondence with the passage of the frontal blocks of the slide, then negative, when the rear slide blocks are passing through this cross-section. The above considerations are valid even for the simulation corresponding to the slide surface II, which are therefore not shown in this paper. The tsunami model was run for 800 s for both reference cases. Time histories computed on some characteristic points on the coast of Stromboli and on points off-shore around the island ŽFig. 1. are graphed in Fig. 12 for the entire computation period. The results of the two simulations are very similar for all selected points. This means that the difference in the extent of the forcing area is not sufficient to produce different water waves in most of the basin, apart as

Fig. 11. Diagrams of Et h vs. time evaluated for groups, each of the three grid points, on four distinct cross-sections of slide surface I.

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

117

Fig. 12. Time histories of water elevation in various points of the tsunami grid for both cases I and II. Coastal points are denoted by the corresponding geographic name, while symbols P1, P2, etc. are used for offshore nodes. The position of all points may be seen in Fig. 1.

expected, from the area that is involved by source II, but not by source I. Water elevation fields are plotted only for case I in Fig. 13, with eight snapshots given every 50 s up to 400 s. The main inferences from Figs. 12 and 13 are expressed in the following. The landslide generates progressively a principal pattern of water waves that in the mature stage consists of three main elements: a leading positive front, a following deep trough almost double in amplitude, and a final positive wave. This system is clearly recognisable in its full form in the snapshot at t s 100 s, while at t s 50 s only the frontal wave with the following trough had time to form. Furthermore, if we consider the time histories of the offshore points P1–P4, which are located within the slide trajectory ŽP2 and P3, see Fig. 1., or not too far from it to the east ŽP1. and to the south–west ŽP4., we can see that this triple signature is quite evident and is followed by a long train of secondary oscillations. Notice that

118

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

Fig. 13. Fields of water elevation computed for case I every 50 s with origin time being the landslide initiation instant. Positive elevations are contoured with solid lines, negative elevations with dashed lines. Labels are in cm. The island of Stromboli is dark grey. The instantaneous underwater position of the landslide is also shown in pale grey Ž t s 50–150 s.. Landslide computations are artificially stopped as soon as the slide reaches the seaward end of slide surface I, which occurs around 160 s. After this time, there is no forcing anymore but only free tsunami propagation.

the system is created by the landslide, but after generation, it propagates with its own velocity governed by water depth. In deep waters it is much faster than the landslide, and precedes it. The wavelength of the tsunami has a scale comparable with the landslide and with the coastal morphology of the island, and this factor dominates the propagation of waves near-shore. In fact, part of the wave energy is trapped near the coast due to very strong wave refraction caused by the island slope. Wave trapping around ocean islands and sea-mounts was first theoretically investigated by Longuet-Higgins Ž1967., and then explored both analytically and numerically by several authors ŽLautenbacher, 1970; Tinti and Vannini, 1994, 1995; Piatanesi and Tinti, 1998b. who showed it to be a very efficient process of maintaining energy for a long time around these small-scale topographic features. In our case, wave trapping is able to transform completely the characteristic pat-

tern that the waves have offshore. Along the shore, the positive leading front advances slowly. At t s 100 s, it has already travelled along the northern coast of Stromboli from Piscita` to Ficogrande and is about to inundate Punta della Lena Žsee Fig. 1 for names of localities and Fig. 12 for the corresponding time histories.. In the following, we refer to this front as the one propagating clockwise in the plots. On the opposite coast, the front is between Ginostra and Secche di Lazzaro, progressing anticlockwise toward Punta Lena. This leading front is minor and is followed by much larger water oscillations. At around t s 200 s both fronts have reached the lee of the island and the entire island coast is affected afterwards by tsunami waves. From the fields corresponding to times greater than 200 s, it is evident that the main wave system has left the basin and the residual energy is trapped in proximity to the island in the form of waves that possess larger amplitude near-

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

shore. Their period is in the order of 100–150 s as deduced from the plotted time-histories. Their wavelength is about one third of the total length of the island coast, as may be roughly estimated from the number of lobes Žabout three positive and three negative lobes. in the field pattern of snapshots corresponding to t s 300 s and onwards. Considering the time histories of the coastal points, we can distinguish them into three different groups. The first is formed by the points that are close to the source, on the island side involved by the landslide. Here the typical tsunami signature consists of a large depression followed by a large water elevation and a series of remarkably smaller oscillations Žsee plots of Punta Frontone, Sciara del Fuoco, and Punta Chiappe, ordered north to south on the cost.. In the second group we class points on the northern and southern Stromboli coast ŽPiscita, ` Scalo dei Balordi, Ficogrande to the north, and Ginostra, Secche di Lazzaro to the south.. Time histories here are characterised by a large depression followed by a very narrow large positive peak, with a very rapid transition from negative to positive elevation values. Following oscillations are smaller, but not too much. In particular, the presence of a wide saddle Ždepth less than 100 m. in the submerged NE prolongation of Stromboli Island, supporting a shallow-water depositional area near-coast, in front of Scalo dei Balordi and Ficogrande, explains the higher elevation here Ž35–40 m. of the positive peak as due to the lower water depth. The third group comprises the points on the lee of the island, with the exception of Malpasseddu, that is, from north to south, Punta della Lena, Punta dell’Olmo and Punta Lena, where there is no large difference between the amplitude of the first and the following oscillations, though a general trend of amplitude decay can be observed. The time history of Malpasseddu is very peculiar because the amplitude of the positive peak is extraordinarily high at about 50 m. This extreme value can be physically explained by invoking the superposition of the wave fronts travelling around the island clockwise and anticlockwise and meeting here with constructive interference. The final observation is that from all time histories it appears that the waves excited by the landslide attack the coast of Stromboli with very large amplitude and within a very short time: 20–30 m in most places, 10 m in few others, and more than

119

50 m, as already noted, in Malpasseddu which is inundated by the positive large pulse about 5 minutes after the landslide initial time.

6. Sensitivity analysis for the tsunami The sensitivity analysis performed on the bottom friction and water resistance parameters of the landslide model cannot be considered complete if the effects on the ensuing tsunami are not investigated. Let us start with the bottom friction coefficients ma and m w . The landslides computed in previous tests and illustrated in Figs. 5 and 6 have been used as input of the tsunami model to simulate water waves. Fig. 14 displays the calculated tide gauge records in all localities shown in Fig. 1, and permits an easy comparison with the results corresponding to the reference case shown in Fig. 12 and commented in the above section. It may be seen that all cases exhibit similar pattern of oscillations but there are important differences that have to be pointed out. The most relevant disparity concerns wave amplitude. Slides with higher speed Žcases with ma s m w s 0.02, and ma s 0.09, m w s 0.02. excite larger water waves. This suggests that in general the wave amplitude goes along with the slide speed. However, this is true only in certain circumstances. As shown in idealised simple cases by Harbitz and Elverhøi Ž1999. and by Tinti and Bortolucci Ž1999b., a decisive factor for water waves generation is the Froude number Fr, defined in this context as the ratio between the slide speed and the local water wave celerity, depending on local water depth. The closer Fr is to the critical value of 1, the larger is the tsunami potential of the slide. Therefore, for subcritical flows Ž Fr - 1., an increase in slide speed corresponds to an increase in wave amplitude, while the opposite is true in supercritical regimes Ž Fr ) 1.. Typically, for masses sliding along a slope, Fr initially grows from zero to a maximum value Žthat may or may not be larger than 1 depending on slide dynamics., and then comes back progressively to zero in deep water. Therefore, the initial and final regimes are certainly largely subcritical, while at intermediate times a critical condition may be approached. In all our examples, since the slope is very steep, the slide reaches deep ocean very soon, and

120

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

Fig. 14. Time histories of water elevation at various points of the tsunami grid for four cases corresponding to different sets of the bottom friction coefficients for the slide motion. See caption of Fig. 12 for other details.

tends to be subcritical even in correspondence with its peak velocities. Looking at the graphs of Fig. 14, it may be seen that in the generation region, both near-shore ŽSciara del Fuoco, Punta Frontone and Punta Chiappe. and offshore Žrecords of P1–P4 within or close to the slide path., the main wave signal produced by the fastest slides may be twice larger than the one caused by the slowest ones. The second deviation regarding the occurrence times is a more obvious consequence of the slide different speeds. Peak values of troughs and crests of the more rapid slides are attained at earlier times. Notice however that arrival times of the signal, that is, the times when the record first departs from zero, do not change very much from case to case, since the initial disturbance of the water surface occurs exactly at the slide initiation time in all tests and its propagation is mostly governed by sea bathymetry rather than slide dynamics. A further remark concerns records exhibit-

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

121

Fig. 15. Time histories of water elevation in various points of the tsunami grid for three cases corresponding to different sets of the water resistance coefficients for the slide motion. See caption of Fig. 12 for other details.

ing very large and narrow peaks, such as first of all, Ficogrande and Scalo dei Balordi, and also Ginostra, Punta Lena. Here the discrepancies in the main oscillations between the cases of larger and smaller slide bottom friction are the most relevant. This may be explained particularly by the presence of very shallow waters near-shore. The incoming front tends to become steeper and steeper, turning to an almost vertical wall of water Žan image that is often mentioned by eyewitnesses of large tsunamis; see e.g., accounts on the recent 1998 Papua New Guinea catastrophic event reported by Davies, 1998.. Wave steeping is testified in tide gauge records by the strong asymmetry of the wave oscillations with almost instantaneous water level rise, followed by a comparatively much slower water withdrawal. This process is inherently nonlinear and is extremely magnified by large approaching waves. Finally we observe that wave amplitudes prove to be very sensi-

122

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

tive to the values of the bottom friction coefficients, since they change significantly in the various cases. It is however worth emphasising that even in the most conservative case, the waves observed at the coast are extremely high, in the order of several meters up to more than 20 m in worst case ŽMalpasseddu.. The simulations of tsunamis carried out for various sets of water resistance coefficients give results that are shown in Fig. 15. Since the set C1 s C2 s 0.004 yields slide peak velocities that are unrealistic Žcf. Fig. 7., they have been not taken into account for tsunami calculations. We have also omitted to plot the results corresponding to the set C1 s 0.008, C2 s 2, since they have been displayed as case II of Fig. 12. The main considerations made for the experiments on the bottom friction coefficients apply also here. The main difference between the two sets of tests concerns an aspect that was already underlined in commenting the slide speed plots of Figs. 5 and 7. We noticed that it takes time before the water resistance terms can exert a sensible action on the slide, and therefore they cannot affect significantly its initial motion. On the contrary, bottom friction is active since the very beginning of sliding. This reflects on the first manifestation of the tsunami at the various locations. In all panels of Fig. 14, the curves associated with different sets of m x depart exactly at the corresponding arrival times, while in Fig. 15 the records tend to be superimposed for a longer time, covering most of the first water rise. Analogously to the previous analysis, the main conclusion of the sensitivity study on the water resistance coefficients is that they affect remarkably the amplitude of tsunami waves. However, even in the case with the maximum resistance, the engendered waves reach absolutely conspicuous amplitude at the coast. In the foregoing, much attention has been given to the tsunamis calculated along the coast of Stromboli, but it is of interest also to have a more general view of the wave field within the computational domain. To this purpose, Figs. 16 and 17 display the total energy of the tsunami in the upper diagram, and the instantaneous maximum and minimum values of the sea water elevations, in the lower plot, calculated in the course of the experiments concerning the coefficients m x and C1 , C2 , respectively. The total energy of the tsunami in the framework of the shallow water

Fig. 16. Analysis of cases with different values of bottom friction coefficients. Ža. Tsunami energy vs. time. Around 100 s curves reach their peak value. The following decrease is due to waves leaving the basin covered by the computational mesh. Žb. Water elevation maxima and minima computed in the grid vs. time. Notice that values plotted at two consecutive times may refer to different nodes of the grid.

approximation is defined as the sum of its potential energy, which is proportional to the integral of z 2 , and of its kinetic energy, which is proportional to the integral of the product h times u P u, where both integrals extending over the entire domain of computation Žsee Tinti and Bortolucci, 1999b.. All energy curves exhibit an increasing trend up to about 100 s, followed by a progressive diminution, owing to the fact that the offshore wave system crosses the boundary of the basin used for the computations. Notice that at 200 s most of the energy has already left the basin. From our previous analysis, we know that small amount of remaining energy concentrates in shallow water around the island, and several coastal places are attacked by the main wave system at later times. Observe that the difference between

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

123

7. Conclusions and implications for hazard assessment

Fig. 17. Analysis of cases with different values of water resistance coefficients. Ža. Tsunami energy vs. time. Žb. Water elevation maxima and minima computed in the grid vs. time. Other details in Fig. 16

the most and the least energetic wave field is about one order of magnitude. The bottom diagrams of Figs. 16 and 17 are built by plotting the maximum and minimum values of z computed in the field at any given time. Thus, two consecutive points in the same curve may refer to two different nodes, possibly even very distant apart, of the numerical grid. By construction, they convey information on the whole wave field. It may be seen that large positive waves are present somewhere in the domain in the time interval 100–300 s, and that they decay significantly at later times. Observe that the main signal reaches Malpasseddu, the place of Stromboli where the highest wave is computed, just about the time of 300 s. From these graphs it may be also seen that the most rapid slides produce positive crests that are by far larger than the negative troughs, while crest and trough amplitudes are comparable when the slide is relatively slow.

The simulations of a tsunami induced by a massive landslide on the western flank of Stromboli have been performed according to a scenario envisaged on the basis of studies concerning the most recent episode of sector collapse, which took place less than 5000 years ago and formed the SdF scar. The mass involved in the said event was huge, with volume estimates ranging from 0.97 to 1.81 km3. However, at least two previous sector collapses of similar entity Žvolumes of the order of 1–2 km3 . occurred in the last 13 ka, which means one event every 4.3 ka on the average. Moreover, the poor stability of the present-day Sciara del Fuoco, whose subaerial part Žmostly made of loose clastic material. has anomalously grown over an exceptionally steep slope without adequate buttressing reinforcement at the foot, represents a major potential hazard ŽBarberi et al., 1993.. In this view, computer experiments have explored the consequences, in the present-day bathymetric framework, of sliding of a mass with volume close to the smaller of the above estimations, but they show that very high amplitude waves are excited in the sea even in this case. Two distinct slide surfaces have been assumed; one wider, the other narrower, both oriented toward the NE and N for bathymetric constraints. The resulting tsunami seems to be insensitive to the considered difference in the slide surface extent, at least on coastal locations. A thorough sensitive analysis has been performed on the coefficients of the slide model that influence slide energy dissipation, that is, the bottom friction and the water resistance coefficients. The motion of the slide and the consequent size of the tsunami have been proven to depend significantly on the values assumed by these parameters, while other features, such as tsunami general wave pattern, are almost insensitive. Nevertheless, even the most conservative cases have produced very important water waves. The most interesting results of the computations are that wave energy propagates mainly toward NE and N, in the form of a circular wave system consisting of two positive crests and a deep trough in between, and that considerable energy remains confined close to Stromboli in the form of very large

124

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

trapped waves, travelling round the island in both directions with slow loss of energy. The amplitude of the waves depends strongly on the set of parameters used to compute the slide motion, but in all cases the waves are very energetic. These waves attack all the coast of Stromboli with amplitude changing from place to place. The computed amplitude is about 10–30 m in the source region on the NW coast, 10–40 m in front of the N coast Žnear Ficogrande, where part of the Stromboli village lies close to the sea. and as much as 20–50 m on the opposite SE coastal segment close to Malpasseddu. These waves reach even the most remote part of Stromboli within 200 s of landslide initiation with largest effects being produced within 300 s. The recognition that very large waves can quickly attack all coasts of Stromboli as the consequence of the collapse, poses problems that cannot be neglected in view of evaluation of hazard and of mitigation of risk on the island, and demands for further investigations to be carried out. The main questions to face on the tsunami concern both the past and the future. There is some chance that evidence of tsunami occurrence may be found from possible deposits on the island: the record of erosional and depositional processes related to tsunami waves run-up should have been preserved, buried by volcanic ashes, other ejecta, or lava. Recent research has shown that paleotsunamis can be successfully identified by searching for such coastal geological signatures. Marine sand layers found on eastern Scotland and western Norway, and attributed to the tsunami generated by a large underwater Holocene slide in the Norway Sea Žsee Dawson et al., 1988; Long et al., 1989; Bondevik et al., 1997., supports the hope that similar deposits might be found at Stromboli if adequate research programmes are purposely planned and undertaken. As for the future, if we assume that another gravitational collapse occurs at Stromboli Žfor example triggered by seismicity or by some stronger explosion in the volcano activity., we could ask how hazardous the generated waves could be at the island Žwhich, during the summer season is quite populated. as well as at other islands of the Aeolian Archipelago, or at the coast of the Italian mainland. The answer to this requires that further landslide as well as tsunami simulations be run to investigate far-field large-scale propagation of water waves and coastal wave pro-

cesses such as wave draw-down and wave run-up, and wave flooding and wave breaking. This paper represents only a contribution to research that is far from being concluded. Future simulations should yield advanced water-impact computations using more sophisticated hydrodynamical models and more accurate databases for the coastal bathymetry of Stromboli. 8. Notation ai k Ai k bi k ci k C1a C2a C1w C2w dik ei k fi k g hik Ii k lik mk pi k ri k Ti k uik Õi k Õi k Vk wi k xik a ik g

acceleration term due to gravity and buoyancy block basal area acceleration term due to bottom friction acceleration term due to frontal and upper surface stresses upper resistance coefficient for subaerial slide frontal resistance coefficient for subaerial slide upper resistance coefficient for submarine slide frontal resistance coefficient for submarine slide effective interaction acceleration term interaction coefficient effective interaction coefficient gravity acceleration block height deformability interval block length block mass ratio of the centripetal acceleration to the gravity acceleration pore pressure coefficient bottom friction force pre-interaction barycentre velocity post-interaction barycentre velocity velocity of the block front block volume block width position of the block front inclination of the slide path shape coefficient of the block–block interaction law

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

da dw l ma mw jik r rX rU s s1 D ti

friction angle for subaerial slide friction angle for submarine slide intensity coefficient of the block– block interaction law bottom friction coefficient for subaerial slide bottom friction coefficient for submarine slide barycentre position density of the landslide density of the ambient fluid effective density ratio: r Xrr dynamic deformability coefficient of the block–block interaction law static deformability coefficient of the block–block interaction law time step at time t i

125

The effective interaction acceleration term may be computed as: d i k s Ž Õi k y u i k . rDt i

Ž 3.

The velocity, Õi k , of the frontal surface of the block, is computed by linear interpolation from the block barycentre velocities: Õi k s Õi , ky1 q Ž Õi k y Õi , ky1 . Ž x iy1, k y j iy1, ky1 . = Ž j iy1, k y j iy1, ky1 .

y1

Ž 4.

The position on the slide path of the block frontal surface and of the block barycentre are computed by means of the expressions: x i k s x iy1, k q 0.5 Ž Õi k q Õiy1, k . D t i

Ž 5.

j i k s 0.5 Ž x i k q x i , ky1 .

Ž 6.

The acceleration term, a i k , due to gravity and buoyancy forces, is: a i k s Ž 1 y r U . g sin a i k

Acknowledgements This work was funded through a research grant of the Italian «Gruppo Nazionale di Vulcanologia» ŽGNVrCNR. and through a contract of the European Communities ŽINTAS-RFBR-95-1000.. We are indebted to G. Macedonio and A. Sodi ŽCNR, Pisa. for software elaborations of the Sciara del Fuoco collapse scar.

The equations to compute the dynamic of variables of the block k at the time t i are given in this appendix. The first step of the computation loop at the time t i is the calculation of the block barycentre pre-interaction velocity u i k :

Ž 1.

The change from pre- to post-interaction velocity Õi k is assumed to be expressed by the following law: Õi , kq1 y Õi k s e iy1, k Ž u i , kq1 y u i k . ,

Ž 2a .

which applies to all blocks of the landslide under the constraint of total momentum conservation, i.e.: Ým k Õi k s Ým k u i k

The bottom friction term is given by: bi k s Ti k Ž 1 y r U . Ž r A i k h i k .

y1

Ž 8a .

where the friction force, Ti k , has the expression Žsee Hungr, 1995.: Ti k s r gA i k h i k m Ž cos a i k q pi k . Ž 1 y ri k . tan d x

Ž 8b . and the bottom friction coefficient m x is related to the friction angle d x by the usual relationship:

Appendix A

u i k s Õiy1, k q Ž a iy1, k y biy1, k y c iy1 , k . D t i

Ž 7.

Ž 2b .

m x s tan d x

Ž 8c .

The subscript x is used to denote the ambient where the slide motion takes place: it is ‘‘a’’ in the subaerial part of the motion and ‘‘w’’ when the blocks are underwater. The resistance term due to the stresses on the frontal and upper block surfaces is: c i k s 0.5 r U Ž C1 x Õi2krh i k q C2 x n i2, kq1 D h i k wi krVk .

Ž 9. where D h i k s maxŽ0, h i k y h i, kq1 . and the subscript x is used for the resistance coefficients to allow them to assume different values outside Ž x '‘‘a’’. and inside the water Ž x '‘‘w’’.. The interaction coeffi-

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

126

cient e i k is calculated by means of a four-parameter law: g

e i k s f i k Ž 1 y sy1 <1 y D j i krD j iy1, k < . , D j i k g Iiy1, k e i k s 0,

Ž 10a.

D j i k f Iiy1, k

Ž 10b.

where D j i k s j i, kq1 y j i k is the separation between the barycentres of blocks k q 1 and k with respective positions j i, kq1 and j i k , and Iiy1, k is the dynamic deformability interval corresponding to the time t iy1 defined by: Iiy1, k ' Ž 1 y s . D j iy1, k , Ž 1 q s . D j iy1, k

Ž 11 .

The effective interaction coefficient f i k is given by: fi k s Ž 1 y l. ,

if D j i k G D j 0 , k

Ž 12a.

f i k s Ž 1 y l . Ž 1 y sy1 1 Ž 1 y D j i k rD j 0, k . . , if Ž 1 y s 1 . D j 0, k - D j i k - D j 0, k f i k s 0,

if D j i k F Ž 1 y s 1 . D j 0, k

Ž 12b. Ž 12c.

where D j 0, k s j 0, kq1 y j 0, k is the distance between the barycentres of blocks k q 1 and of block k at the initial time t s 0.

References Aida, I., 1975. Numerical experiments of the tsunami associated with the collapse of Mt. Mayuyama in 1792. Zisin 2 Ž28., 449–460, Žin Japanese with abstract in English.. Assier Rzadkiewicz, S., Mariotti, C., Heinrich, Ph., 1997. Numerical simulation of submarine landslides and their hydraulic effects. J. Waterw., Port, Coastal, Ocean Eng. Div., Am. Soc. Civ. Eng. 123, 149–157. Barberi, F., Rosi, M., Sodi, A., 1993. Volcanic hazard assessment at Stromboli based on review of historical data. Acta Vulcanologica 3, 173–187. Bondevik, S., Svendsen, J.I., Mangerud, J., 1997. Tsunami sedimentary facies deposited by the Storegga tsunami in shallow marine basins and coastal lakes, western Norway. Sedimentology 44, 1115–1131. Campbell, C.S., 1990. Rapid granular flows. Annu. Rev. Fluid Mech. 22, 57–92. Campbell, C.S., Clearly, P.W., Hopkins, M., 1995. Large-scale landslide simulations: global deformation, velocities and basal friction. J. Geophys. Res. 100, 8267–8283. Davies, H.L., 1998. The Sissano tsunami. Published by the author, pp. 34.

Dawson, A.G., Long, D., Smith, D.E., 1988. The Storegga slides: evidence from eastern Scotland of a possible tsunami. Mar. Geol. 82, 271–276. Fine, I.V., Rabinovich, A.B., Kulikov, E., Thomson, R.E., Bornhold, B.D., 1999. Numerical modelling of landslide-generated tsunamis with application to the Skagway Harbor tsunami of November 3, 1994. In: International Conference on Tsunamis, Paris, France, May 26–28, 1998. Commissariat a` l’Energie Atomique ŽCEA., pp. 212–223. Francalanci, L., Manetti, P., Peccerillo, A., 1989. Volcanological and magmatological evolution of Stromboli volcano ŽAeolian Islands.: the roles of fractional crystallization, magma mixing, crustal contamination and source heterogeneity. Bull. Volcanol. 51, 355–378. Gabbianelli, G., Romagnoli, C., Rossi, P.L., Calanchi, N., 1993. Marine geology of the Panarea–Stromboli area ŽAeolian Archipelago, Southeastern Tyrrhenian Sea.. Acta Vulcanologica 3, 11–20. Gillot, P.Y., Keller, J., 1993. Radiochronological dating of Stromboli. Acta Vulcanologica 3, 69–77. Hammack, J.L., 1973. A note on tsunamis: their generation and propagation in an ocean of uniform depth. J. Fluid Mech. 60, 769–799. Harbitz, C.B., 1992. Model simulation of tsunamis generated by Storegga slides. Mar. Geol. 104, 1–21. Harbitz, C.B., Elverhøi, A., 1999. On tsunami characteristics and submarine dynamics. In: International Conference on Tsunamis, Paris, France, May 26–28, 1998. Commissariat a` l’Energie Atomique ŽCEA., pp. 257–266. Harbitz, C.B., Pedersen, G., Gjevik, B., 1993. Numerical simulations of large water waves due to landslides. J. Waterw., Port, Coastal, Ocean Eng. Div., Am. Soc. Civ. Eng. 118, 249–266. Heinrich, Ph., 1991. Nonlinear numerical model of landslide-generated water waves. Int. J. Eng. Fluid Mech. 4, 403–416. Heinrich, Ph., 1992. Nonlinear water waves generated by submarine and aerial landslides. J. Waterw., Port, Coastal, Ocean Eng. 118, 249–266. Hopkins, M.A., 1987. In: Particle Simulation 1 Dept. Civ. Environ. Eng., Clarkson University, Potsdam, NY, Rep. No. 87-7. Hornig-Kjarsgaard, I., Keller, J., Koberski, U., Stadlbauer, E., Francalanci, L., Lenhart, R., 1993. Geology, stratigraphy and volcanological evolution of the island of Stromboli, Aeolian arc, Italy. Acta Vulcanologica 3, 21–68. Hungr, O., 1995. A model for the runout analysis of rapid flow slides, debris flows, and avalanches. Can. Geotech. J. 32, 610–623. Hutchinson, J.N., 1986. A sliding-consolidation model for flow slides. Can. Geotech. J. 23, 663–677. Hutter, K., Koch, T., Pluss, ¨ C., Savage, S.B., 1995. The dynamics of avalanches of granular materials from initiation to runout. Part II: Experiments. Acta Mech. 109, 127–165. Hwang, S.L., Divoky, D.J., 1970. Tsunami generation. J. Geophys. Res. 75, 6802–6817. Imamura, F., Gica, E.C., 1996. Numerical model for tsunami generation due to subaqueous landslide along a coast. A case of the 1992 flores tsunami — Indonesia. Science Tsunami Hazards 14, 13–28.

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128 Iverson, R.M., 1997. The physics of debris flows. Rev. Geophys. 35, 245–296. Iwasaki, S.I., 1997. The wave forms and directivity of a tsunami generated by an earthquake and a landslide. Science Tsunami Hazards 15, 23–40. Jiang, L., LeBlond, P.H., 1992. The coupling of a submarine slide and the surface waves which it generates. J.Geophys. Res. 97, 12731–12744. Jiang, L., LeBlond, P.H., 1993. Numerical modeling of an underwater Bingham plastic mudslide and the waves which it generates. J.Geophys. Res. 98, 10303–10317. Jiang, L., LeBlond, P.H., 1994. Three-dimensional modeling of tsunami generation due to submarine mudslide. J. Phys. Oceanogr. 24, 559–572. Johnsgard, H., Pedersen, G., 1996. Slide generated waves in near-shore regions. A Lagrangian description. Phys. Chem. Earth 21, 45–49. Jones, A.T., Mader, C.L., 1996. Wave erosion on the southeastern coast of Australia: tsunami propagation modelling. Aust. J. Earth Sci. 43, 479–483. Kajiura, K., 1963. The leading wave of a tsunami. Bull. Earthquake Res. Inst. 41, 535–571. Kawata,Y., Benson, B.C., Borrero, J.C., Davies, H.L., De Lange,W.P., Imamura,F., Letz, H., Nott, J., Synolakis, C.E., 1999. Tsunami in Papua New Guinea was as intense as first thought, EOS, Transactions, AGU, 80, pp. 101, 104–105. Kienle, J., Kowalik, Z., Murty, T.S., 1987. Tsunamis generated by eruptions from St. Augustine volcano, Alaska. Science 236, 1442–1447. Kokelaar, P., Romagnoli, C., 1995. Sector collapse, sedimentation and clast-population evolution at an active island-arc volcano: Stromboli, Italy. Bull. Volcanol. 57, 240–262. Korner, H.J., 1976. Reichwette und Geschwindigkeit von ¨ Bergsturzen und Fleishshneelawinen. Rock Mech. 8, 225–256. Kulikov, E., Rabinovich, A.B., Thomson, R.E., Bornhold, B.D., 1996. The landslide tsunami of November 3, 1994, Skagway Harbor, Alaska. J. Geophys. Res. 101, 6609–6615. Labazuy, Ph., 1991. Instabilites au cours de l’evolution d’un ´ edifice volcanique, en domaine intraplaque oceanique: Le ´ .. PhD Thesis, Piton de la Fournaise ŽˆIle de La Reunion ´ University Blaise Pascal, Clermont-Ferrand II. Lautenbacher, C.C., 1970. Gravity wave refraction by islands. J. Fluid Mech. 41, 655–672. Long, D., Smith, D.E., Dawson, D.E., 1989. A Holocene tsunami deposit in eastern Scotland. J. Quat. Sci. 4, 61–66. Longuet-Higgins, M.S., 1967. On the trapping of wave energy round islands. J. Fluid Mech. 29, 781–821. Ma, K.F., Satake, K., Kanamori, H., 1991. The origin of the tsunami excited by the 1989 Loma Prieta earthquake — faulting or slumping. Geophys. Res. Lett. 18, 637–640. McEwen, A.S., Malin, M.C., 1989. Dynamics of Mount St. Helens’ 1980 pyroclastic flows, rockslide-avalanche, lahars and blast. J. Volcanol. Geotherm. Res. 37, 205–235. Mohrig, D., Whipple, K.X., Hondzo, M., Ellis, C., Parker, G., 1998. Hydroplaning of subaqueous debris flows. GSA Bulletin 110 Ž3., 387–394.

127

Naranjo, J.A., Francis, P., 1987. High velocity debris avalanche at lastarria volcano in the north Chilean Andes. Bull. Volcanol. 49, 509–514. Noda, E.K., 1970. Water waves generated by landslides. J. Waterw., Harbors Coastal Eng. Div., Am. Soc. Civ. Eng. 96, 835–858. Noda, E.K., 1971. Water waves generated by a local surface disturbance. J. Geophys. Res. 76, 7389–7400. Pasquare, ` G., Francalanci, L., Garduno, V.H., Tibaldi, A., 1993. Structure and geologic evolution of the Stromboli volcano, Aeolian Islands, Italy. Acta Vulcanologica 3, 79–89. Perla, R., Cheng, T.T., McClung, D.M., 1980. A two-parameter model of snow avalanche motion. J. Glaciol. 26, 197–207. Piatanesi, A., Tinti, S., 1998a. A revision of the 1693 eastern Sicily earthquake and tsunami. J. Geophys. Res. 103, 2749– 2758. Piatanesi, A., Tinti, S., 1998b. Finite-element numerical simulations of tsunamis generated by earthquakes near circular islands. Bull. Seismol. Soc. Am. 88, 609–620. Piatanesi, A., Tinti, S., Gavagni, I., 1996. The slip deistribution of the 1992 Nicaragua earthquake from tsunami data. Geophys. Res. Lett. 23, 37–40. Pierson, T., 1995. Flow characteristics of large eruption-triggered debris flows at snow-clad volcanoes: contraints for debris flow models. J. Volcanol. Geotherm. Res. 66, 283–294. Romagnoli, C., 1990. Caratterizzazione morfostrutturale e vulcanologica sottomarine delle Isole Eolie ŽTirreno meridionale.: ipotesi di stadi evolutivi dei complessi. Unpublished Ph.D. Thesis, University of Bologna, 141 pp. Romagnoli, C., Kokelaar, P., Rossi, P.L., Sodi, A., 1993. The submarine extension of Sciara del Fuoco feature ŽStromboli isl..: morphologic characterization. Acta Vulcanologica 3, 91– 98. Rosi, M., 1980. The island of Stromboli. Rend. Soc. Ital. Mineral. Petrol. 36, 345–368. Savage, S.B., Hutter, K., 1991. The dynamics of avalanches of granular material from initiation to runout. Part 1. Analysis. Acta Mech. 86, 201–223. Slingerland, R.L., Voight, B., 1979. Occurrences, properties and predictive models of landslide-generated water waves. In: Voight, B. ŽEd.., Rockslides and Avalanches, 2. Engineering Sites. Developments in Geotechnical Engineering 14B Elsevier, pp. 317–397. Sousa, J., Voight, B., 1995. Multiple pulsed debris avalanche emplacement at Mount St. Helens in 1980: evidence from numerical continuum flow simulations. J. Volcanol. Geotherm. Res. 66, 227–250. Takahasi, T., 1981. Debris flow. Annu. Rev. Fluid Mech. 13, 57–77. Tappin, D.R., Matsumoto, T., Watts, P., Satake, K., Mc Murtry, G.M., Matsuyama, M., Lafoy, Y., Tsuji, Y., Kanamatsu, T., Lus, W., Iwabuchi, Y., Yeh, H., Matsumotu, Y., Nakamura, M., Mahoi, M., Hill, P., Crook, K., Anton, L., Walsh, J.P., 1999, Sediment slump likely caused 1998 Papua New Guinea tsunami, EOS, Transactions, AGU, 80, pp. 329, 334, 340 Tinti, S., Bortolucci, E., 1999a. Strategies of optimal grid genera-

128

S. Tinti et al.r Journal of Volcanology and Geothermal Research 96 (2000) 103–128

tion for finite-element tsunami models. In: International Conference on Tsunamis, Paris, France, May 26–28, 1998. Commissariat a` l’Energie Atomique ŽCEA., pp. 269–294. Tinti, S., Bortolucci, E., 1999b. Energy of water waves induced by submarine landslides. Pure Appl. Geophys., Žin press.. Tinti, S., Gavagni, I., 1994. A method for reducing the propagation noise in FE modelling of tsunamis. Science Tsunami Hazards 12, 77–92. Tinti, S., Gavagni, I., 1995. A smoothing algorithm to enhance finite-element tsunami modelling: an application to the 5 February 1783 Calabrian case, Italy. Natural Hazards 12, 161–197. Tinti, S., Maramai, A., 1996. Catalogue of tsunamis generated in Italy and in Cote ˆ d’Azur: a step towards a unified catalogue of tsunamis in Europe. Ann. Geofis. 39, 1253–1299, ŽErrata Corrige, Ann Geofis, 40, 781.. Tinti, S., Piatanesi, A., 1996. Numerical simulations of the tsunami induced by the 1627 earthquake affecting Gargano, southern Italy. J. Geodynamics 21, 141–160. Tinti, S., Vannini, C., 1994. Theoretical investigation on tsunami induced by seismic faults near ocean islands. Marine Geodesy 17, 193–212. Tinti, S., Vannini, C., 1995. Tsunami trapping near circular islands. Pure Appl. Geophys. 144, 595–619. Tinti, S., Gavagni, I., Piatanesi, A., 1994. A finite-element numerical approach for modelling tsunamis. Ann. Geofis. 37, 1009– 1026. Tinti, S., Bortolucci, E., Vannini, C., 1997. A block-based theoretical model suited to gravitational sliding. Natural Hazards 16, 1–28.

Tinti, S., Bortolucci, E., Piatanesi, A., 1998. Un modello Lagrangiano per lo studio della dinamica delle frane. In: Atti Convegno Internazionale sulla Prevenzione delle Catastrofi Idrogeologiche, Il Contributo della Ricerca Scientifica, Alba, Italy, 5-7 Novembre 1996. Žin press; in Italian with abstract in English.. Tinti, S., Bortolucci, E., Armigliato, A., 1999. Numerical simulation of the landslide-induced tsunami of 1988 in Vulcano island, Italy. Bull. Volcanol. 61, 121–137. Villeneuve, M., Savage, S.B., 1993. Nonlinear, dispersive, shallow-water waves developed by a moving bed. J. Hydraul. Res. 31, 249–266. Voight, B., Elsworth, D., 1997. Failure of volcano slopes. Geotechnique 47 Ž1., 1–31. ´ Voight, B., Sousa, J., 1994. Lessons from Ontake-san: a comparative analysis of debris avalanche dynamics. Eng. Geol. 38, 261–297. Voight, B., Glicken, H., Janda, R.J., Douglass, P.M., 1981. Catastrophic rockslide avalanche of May 18, USGS Prof. Paper, 1250. In: The 1980 eruptions of Mount St. Helens, Washington. pp. 347–378. Watts, P., 1998. Wavemaker curves for tsunamis generated by underwater landslides. J. Waterw. Port Coastal Ocean Eng. 124, 127–137. Wiegel, R.L., 1970. Water waves generated by landslides in reservoirs. J. Waterw., Harbors Coastal Eng. Div., Am. Soc. Civ. Eng. 96, 307–333. Zanchi, A., Francalanci, L., 1989. Analisi geologico strutturale dell’isola di Stromboli: alcune considerazioni preliminari. Boll. GNV 5, 1027–1044.