KDV THEORY AND THE CHILEAN TSUNAMI OF 1960

We investigate the Chilean tsunami of 1960 to determine the role of KdV dynamics. On the basis of the scales involved, and making use of recent advances, we put on a rigorous foundation the fact that KdV dynamics were not influential in this catastrophic event.


1.
Introduction. Tsunami are without a doubt among the most infamous and least understood natural disasters today. Often referred to in the popular literature by the misnomer "tidal wave", tsunami are generated by large displacements in the sea level, often via seismic activity. Most tsunami -a term from the Japanese for "harbor wave" -are caused by vertical movement along a break in the earth's crust. Other causes can include volcanic collapse, subsidence, as well as landslides. Contrary to popular imagination, a tsunami need be neither large nor destructive -classification is based on origin of the wave or wave period rather than on size. Though between 1861 and 1948 there were more than 15,000 earthquakes recorded, there were only 124 tsunami [3]. Indeed, off the west coast of South America, 1,098 earthquakes have led to only 20 recorded tsunami [3].
As waves of such great scale, generated by complex movements of the earth, and with such devastating consequences for populations surrounding the world's oceans, accurate modeling of tsunami is of utmost importance. One question which has been raised repeatedly is whether the behavior of tsunami at sea can be described by the Korteweg-de Vries equation (see the reviews [27], [31], and [30]). We will pursue this question for one of the greatest tsunami of recorded history -generated by a series of earthquakes in southern Chile on May 22, 1960 -as it propagated from Chile to Hawai'i. These earthquakes, among them the largest ever recorded, resulted from a rupture about 1000 km long and 150 km wide along the fault between the Nazca and South American plates, at a focal depth of 33 km. The principal shock occurring on May 22 at 19:11 GCT registered at 9.5 on the moment magnitude scale, and led to changes in land elevation ranging from 6 m of uplift to 2 m of subsidence -which has been modeled to correspond to an average dislocation of 20 m along the fault, with peaks of more than 30 m [2]. This subsidence extended as far as 29 km inland, resulting in some 10 km 2 of forest around the Río Maullín being submerged by the tides and consequently defoliated [4].
Not only was the principal earthquake at 39.5 • S, 74.5 • W especially powerful, it generated tsunami with an average run-up of 12.2 m and a maximal run-up on the adjacent Chilean coast of 25 m. Over the course of the next day, a number of tsunami wreaked havoc upon the Pacific, taking the lives of more than 2000 people and causing millions of dollars in damages. The initial wave traveled between 670 and 740 km/h, with a wavelength of between 500 -800 km and a height in the open ocean of only 40 cm [3], [20]. Borrowing an example from [30], sitting in a boat in the Pacific, the tsunami wave would take between 45 min to an hour to pass one by while raising the boat by less than one centimeter per minute -hardly noticeable on the open sea. Nevertheless, the tsunami reached amplitudes of 7 m in Kamchatka and 10.7 m in Hilo, Hawai'i, 1 where it caused widespread destruction after traveling 10,000 km in just under 15 hours.

2.
Modeling. The Chilean tsunami of 1960 had wavelengths in excess of 500 km and amplitudes of less than one meter [3] while propagating over the Pacific Ocean, which, though the deepest of the world's oceans, has an average depth of only 4.3 km. These scales lend themselves to modeling with shallow-water long-wave theory, i.e. water depth is small compared to wavelength and depth is large compared to amplitude. We note also that the depth of open ocean across which the 1960 tsunami traveled is relatively uniform, and given that the rupture length exceeded the wavelength of the resulting tsunami, it is reasonable to assume the waves as twodimensional; this is borne out (at least between Chile and Hawai'i) by consulting travel time charts (see [3]).
2.1. The governing equations. Following e.g. [22] we introduce the governing equations for water waves, and then show how these apply specifically to the tsunami of 1960. Since we consider two-dimensional waves, let x be the direction of wave propagation, and y the vertical. Let u denote velocity in x-direction, v that in y-direction. The density ρ, we assume to be constant, for simplicity setting ρ ≡ 1. The pressure P (x, y, t) as well as u(x, y, t) and v(x, y, t) we assume to be suitably differentiable functions. Finally g denotes the acceleration of gravity. The equation of mass conservation arises from the stipulation that the density of water be constant, giving u x + v y = 0, (1) where subscripts denote partial differentiation with respect to x resp. y. Neglecting viscosity, we use the Euler equations where ∂y is the so-called material derivative. We will further consider only irrotational flow -meaning that we set the vorticity, which measures a local spin of the fluid, equal to zero: u y − v x = 0 (3) Although in real water flows, vorticity is rarely absent, it is usually small enough that it does not play a major role in water wave dynamics unless we wish to account for the presence of underlying non-uniform currents [25,24,13,7]. As it is, for the rigorous results we rely on below, irrotational flow is a prerequisite.
We write y = h 0 + H(x, t) for the free surface, where h 0 is an average depth for the water under consideration, and denote the flat bed by y = 0. Note that it is possible to include an analysis of these dynamics with a non-flat bed as in [12,21], provided the variations in bottom topography are limited appropriately. In order to make the derivations that follow as transparent as possible, we will restrict ourselves to the case of a flat bed.
To decouple the motion of the water from that of the air, we introduce the dynamic boundary condition P = P atm on the free surface y = h 0 + H(x, t). (4) Further we assume that particles comprising either of the two fluid surfaces (the free surface and the bed) must stay there -physically this means that the bed is impenetrable and that no particles leave the body of water, implying therefore that velocities along these surfaces have no normal component. This is called the kinematic condition: These equations together comprise the governing equations for water waves.

2.2.
Non-dimensionalization and scaling. Roughly speaking, one might say that nonlinear phenomena arise through the interaction of physical parameters on scales of differing magnitude. As such, given that we are working with physical variables, in order to compare magnitudes meaningfully, the first step is to get rid of their units. Experience has borne out the fact that nonlinear problems are often fruitfully tackled by approximation -by introducing special scales to an otherwise too expansive problem and then considering regimes corresponding to certain values of these scales. Simpler approximate equations in such regimes permit an in-depth analysis of waves enjoying special attributes (such as the existence of solitons [8,23] or stability properties [14].) 2.2.1. Non-dimensionalization. In keeping with this, we introduce h 0 as the typical depth of the water and λ the typical wavelength. These two scales provide the basis for a nondimensional version of the governing equations. The characteristic speed for long gravity waves is taken to be √ gh 0 , and together with the wavelength λ this gives us a time scale for horizontal propagation of the wave, λ/ √ gh 0 , as well as horizontal speed. Care must be taken with the vertical speed v in order to be consistent with (1). These considerations lead us to the following non-dimensional variables, which we denote with the usual variable names: Accordingly, we transform the pressure P into a perturbation of the hydrostatic pressure as follows where p is a new non-dimensional pressure variable. Lastly we set where η is the nondimensional surface profile and a is a typical amplitude.

RAPHAEL STUHLMEIER
The components of the Euler equation (2) under these transformations become where δ = h 0 /λ is the long wavelength or shallowness parameter. Owing to our abuse of notation, the equation of mass conservation (1) remains unchanged. The irrotationality condition (3) becomes The second characteristic parameter enters via the transformation of the boundary conditions in accordance with (7) - (10): where the = a/h 0 is the so-called amplitude parameter.

2.2.2.
Scaling. At this point, we note that and δ between them determine the type of water wave problem under consideration. Looking at the surface boundary conditions, we see that v and p are both proportional to , the wave amplitude. This is sensible, since as → 0 the vertical velocity v → 0 and of course the pressure perturbation p → 0; the free surface is perfectly flat. Taking advantage of this, we define a set of scaled variables where u is scaled similarly for consistency (note that these formal considerations are supported by rigorous results -see [15,5,32]). This leads to the transformation of the system of equations (11) with D Dt = ∂ ∂t + (u ∂ ∂x + v ∂ ∂y ), which we can write explicitly as The equation of mass conservation (1) again remains unchanged, as does the nondimensional irrotationality condition (12). The boundary conditions (13) -(15) become 2.2.3. Approximate Equations. We now have a series of equations which depend on two parameters and δ which measure contributions of amplitude, respectively wavelength, to the problem under consideration. The most common approximations made are → 0 for fixed δ and δ → 0 for fixed . These are known as the linearized problem and shallow-water (or long wave) problem respectively. As noted above, in the first approach, the hitherto unknown free surface becomes the surface y = 1, and in a first approximation, we have a linear problem with dispersive effects. In the latter, a glance at (18) shows that the pressure becomes independent of y; dispersive effects are neglected. For more details see e.g. [22].
While these approximations have been used extensively (and sometimes with some abandon) in the history of water-wave problems, the question remains to what extent any formal asymptotic model can give us relevant results for water waves. Two questions arise: does an asymptotic model provide a good approximation of a solution to the Euler equations, and is the time scale of the model applicable.
Our goal of understanding the dynamics of the Korteweg-de Vries equation (henceforth abbreviated KdV) and its possible application to the Chilean tsunami of 1960 must therefore proceed cautiously. Roughly, we understand that the KdV describes a balance between nonlinearity and dispersion, which means that we will need to retain both parameters and δ to some order in the above equations.
We note that the parameters occur in our equations to order and δ 2 -it turns out that these are precisely the orders that must be retained. It is possible to identify the regime within which a certain model applies, such as via the Ursell number U = aλ 2 /h 3 0 = /δ 2 [33]. An Ursell number U = O(1) corresponds to the classical approach of choosing the parameters such that δ 2 = O( ) as → 0. It may also be pointed out that for any δ as → 0, there should always exist a region of space and time where this balance occurs [10], [11] -although no rigorous results exist to corroborate this. However, we will stick to the classical point of view, for which rigorous underpinnings are available.

The Korteweg-de Vries equation (KdV).
Starting from the equations (17) -(21) we will proceed to derive KdV (introduced in [26]), following the exposition in [22]. As discussed above, we take the classical approach and consider a special choice of parameters, namely δ 2 = O( ) as → 0. In accordance with this choice, we transform the independent variables As above, we require that the condition of mass conservation be satisfied, and so must transform accordingly. Then the equations of motion (17) -(21) along with mass conservation (1) and irrotationality (12) now appear with δ scaled out in favor of : Now we let → 0 and observe that for a first order approximation (24) implies u t + p x = 0 as well as that p is independent of y. Therefore (27) implies that p = η everywhere and (26) that u is independent of y. It then follows from (25) that which satisfies the boundary condition (29). Furthermore, since v = η t on y = 1, (30) implies η t = −u x or η t + u x = 0. This together with mass conservation (25) means that η fulfills the wave equation We know that the wave equation has right-running as well as left-running solutions (cf. [17]), and will follow the right-running waves, consistent with the introduction of the characteristic variable ξ = x − t. We also introduce a slow time scale in order to treat the far-field region, where we consider the regime for ξ = O(1) and τ = O(1). In view of this, we can rewrite the equations of motion (24) -(29) as follows: We would now like to determine an approximate solution in terms of an asymptotic series Ansatz in (for background cf. [29]) by introducing the series expansions: v(ξ, τ, y, ) ∼ n≥0 n u n (ξ, τ, y) p(ξ, τ, y, ) ∼ n≥0 n p n (ξ, τ, y) Notice, however, that we have a problem in (36) and (37): appears both in the coefficients as well as the arguments, making it impossible to equate powers of epsilon as we would like to do. To this end, we perform a transfer of the boundary conditions from y = 1 + η to y = 1 by expanding p(ξ, τ, y), u(ξ, τ, y) and v(ξ, τ, y) in Taylor series about y = 1 as follows: p(ξ, τ, 1 + η) = p(ξ, τ, 1) + p y (ξ, τ, 1) η + 1 2! p yy (ξ, τ, 1) 2 η 2 + . . .
Substitute this into (36) and apply the series expansions for η and p in (39) and (40) to get Analogously substitute the Taylor series for u, v into (37) and expand to get: At leading order ( 0 ) (33) -(38) then reduces to: w 0 = 0 on y = 0 (49) This system is analogous to that in (24) -(29) above -and the analogous arguments lead to the fact that p 0 = η 0 , u 0 = η 0 + C, C a constant which we may assume to be zero, and v 0 = −yη 0ξ . At order 1 we get: Recall that we know from the first order approximation: Taking this into account, and in view of the boundary conditions in the second approximation it is easy to see Now we would like to eliminate η 1 and get an equation solely in η 0 . Notice that Integrating with respect to y yields which on the free surface y = 1 is equal to −η 1ξ + η 0τ + 2η 0 η 0ξ , whereupon the factor −η 1ξ cancels and we have: the Korteweg-de Vries equation. Much has been written about the dynamics associated with the KdV (see e.g. [16]), but we must return to the question of what solutions of the KdV can tell us about water waves -especially whether the governing equations have a solution on the time-scale of our asymptotic model. Recent results in [1] answer this question in the affirmative.
(59) Simply identifying the regime does not suffice however, to see KdV dynamicsindeed only at certain length scales can these dynamics appear. Recall that τ = t and the non-dimensionalisation (7) performed in 2.2.1. The above then transforms into This gives us a length scale for the KdV balance of We note that this length scale was long thought to be x = O(h 0 / ), (cf. the classical results [18], [19]) as is also espoused in the recent survey [30], but (61) provides the correct scale -see also the discussion in [6]. Thus far we have identified a regime, = O(δ 2 ), but what does this mean in practical terms? Given that we need to check whether this regime holds based on real-world data, we take the approach that O(1) allows for deviation by a factor of ten in either direction, as is usually assumed in the hydrodynamical literature (see e.g. [28]). Thus 10 −1 ≤ δ 2 ≤ 10 is a good realization of the KdV regime. Recall that the definitions = a/h 0 , δ = h 0 /λ mean that the above is where we take h 0 to be 4.3 km. Taking a = 0.4 m (cf. [3]), we see that this yields a realistic range of wavelengths between 140 and 1400 km. So far we are in the right regime to see KdV dynamics, but need to find a bound for the distance in which we expect a balance of nonlinearity and dispersion to appear. Given that we consider the 1960 tsunami only between Chile and Hawai'i, a distance of about 10 4 km, (61) means that we need λ < 10 4 km.
(64) However, measurements place the wavelength of the tsunami of May 22, 1960 between 500 -800 km [3] making it unlikely that KdV dynamics played a role. This is further supported by the fact that the first two tsunami waves reaching Hilo, Hawai'i were smaller than the third, most destructive wave -something which should not occur if KdV dynamics were significant for the leading waves of the tsunami. We have deliberately used the wavelength λ because of the relative ease of measurement and error tolerance compared with measuring amplitude. An argument based on the amplitude can be found in [9] and [6]. The propagation distances involved in the 1960 Chilean tsunami are the largest possible on earth, making it one of the best candidates among teleseismic tsunami for the appearance of a KdV balance. Nevertheless, previous studies of this and other large teleseismic events (cf. the recent investigations of the 2004 tsunami [11], [10] concluding that KdV theory did not play a role in this event) support our view that the Chilean tsunami of 1960, though the most widespread ever to affect the Pacific Basin, was not governed by KdV dynamics.