CC BY
47
RUSSIAN JOURNAL OF EARTH SCIENCES, VOL. 12, ES2005, doi:10.2205/2011ES000509, 2011
Hydroacoustic effects in the 2003 Tokachi-oki tsunami source
A. Bolshakova1, S. Inoue2, S. Kolesov1, H. Matsumoto3, M. Nosov1, and T. Ohmachi4 Received 21 November 2011; accepted 28 November 2011; published 10 December 2011.
We process the JAMSTEC ocean-bottom pressure gauges and ocean-bottom seismometers datasets obtained during the 2003 Tokachi-oki tsunamigenic earthquake - the first records which have ever been obtained in a large tsunami source. On these records, we discover the unique phenomenon in tsunami source - hydroacoustic resonance, i.e. manifestation of long-lasting elastic oscillations of water column at the minimal normal frequency (« 0.14 Hz). The concept of a weakly coupled system is applied in 3D numerical simulation of the Tokachi-oki event. First, we simulate earthquake ground motion due to seismic fault rupturing. Then, compressible water column disturbance resulting from the dynamic seismic ground motion is simulated using the velocity of bottom deformation as an input to the water column. Comparison between JAMSTEC in-situ measurements and synthetic signals is carried out. KEYWORDS: Tsunami source; bottom earthquake; water compressibility; dynamic bottom deformation; numerical simulation; ocean-bottom pressure gauges; ocean-bottom seismometers; spectrograms.
Citation: Bolshakova, A., S. Inoue, S. Kolesov, H. Matsumoto, M. Nosov, and T. Ohmachi (2011), Hydroacoustic effects in the 2003 Tokachi-oki tsunami source, Russ. J. Earth. Sci., 12, ES2005, doi:10.2205/2011ES000509.
1 Introduction
Majority of tsunamis springs from strong bottom earthquakes. In simulating a tsunami the water column is usually considered as incompressible medium [e.g. Fujii and Satake, 2008; Kowalik et al., 2005; Titov and Gonzalez, 1997]. The issue of accounting for the compressibility of water in the problem of tsunami generation has been raised repeatedly in the literature [e.g. Chierici et al., 2010; Gisler, 2008; Ka-jiura, 1970; Levin and Nosov, 2008; Li et al., 2009; Miyoshi, 1954; Nosov, 1999, 2000; Nosov et al., 2005, 2007; Ohmachi et al., 2001; Ohmachi and Inoue, 2010; Panza et al., 2000; Sells, 1965; Stiassnie, 2010]. In particular, Nosov [2000] showed theoretically that the compressibility of water is significant only at the stage of tsunami generation by an earthquake, while tsunami propagation and onshore wave run-up can be described as incompressible fluid motion. It was also shown that the main difference in the behaviour of a compressible ocean, as compared to an incompressible model
1Faculty of Physics, M. V. Lomonosov Moscow State University, Moscow, Russia
2 Engineering Seismology Section, Takenaka Research & Development Institute, Tokyo, Japan
3Japan Agency for Marine-Earth Science and Technology, Yokosuka, Japan
4Japan Dam Engineering Center, Tokyo, Japan
Copyright 2011 by the Geophysical Center RAS.
http://elpub.wdcb.ru/journals/rjes/doi/2011ES000509.html
medium, is in the formation of low-frequency elastic oscillations of the water column (~ 0.1 Hz) in the tsunami source area. It is worth noting, frequencies that are traditionally considered in hydroacoustics are much higher (1 — 106 Hz) [Brekhovskikh and Lysanov, 2003].
Till recently, only one of the manifestations of water being compressible in the case of bottom earthquakes, namely the T-phase (or T-waves), has been registered many times and therefore studied relatively well [e.g. Ewing et al., 1950; Okal et al., 2003; Soloviev et al., 1968]. T-waves can propagate efficiently only if their wavelengths fit inside the width of the SOFAR channel, so, in practice their frequency exceeds 2 Hz. Being channelled by the SOFAR, hydroacoustic waves from bottom earthquakes may travel thousands of kilometres before dissipating. This remarkable feature facilitates registration of the T-phase. In contrast to the T-phase, hydroacoustic phenomena of low-frequency (< 1 Hz) can be observed only in the immediate vicinity of the tsunami source.
In past decade, when JAMSTEC (Japan Agency for Marine-Earth Science and Technology) deployed a real-time observatory at the continental slope close to the islands of Japan [Hirata et al. 2002], it became possible to investigate a tsunami formation, including hydroacoustic phenomena, just at its source. The Tokachi-oki earthquake of 2003 turned out to be the first large tsunamigenic seismic event, the epicentre of which was located in the immediate vicinity of the JAMSTEC sensors. Spectral analysis of variations of bottom pressure recorded during the earthquake provided a unique opportunity to reveal manifestations of compress-
ibility of water in tsunami source [Li et al., 2009; Nosov and Kolesov, 2007; Nosov et al.; 2005, 2007; Ohmachi and Inoue, 2010]. Moreover, bottom pressure gauges recorded static pressure changes which were interpreted as vertical uplift of the ocean bottom caused by the event [Mikada et al., 2006; Nosov et al., 2005; Watanabe et al., 2004].
In the present study, we continue the line that had been originated from our previous publications [Inoue et al., 2010; Kolesov et al., 2010; Nosov and Kolesov, 2007; Nosov et al., 2005, 2007; Ohmachi and Inoue, 2010; Ohmachi et al., 2001]. In Section 2, we process the dataset recorded by JAMSTEC sensors during the 2003 Tokachi-oki tsunami-genic earthquake. In contrast to our previous publications, which had been mostly based on the processing of 1 Hz bottom pressure dataset, now we provide a comparative analysis of 10 Hz ocean-bottom pressure gauges dataset and 100 Hz ocean-bottom seismometers dataset as well.
In Section 3, we present results of numerical simulation of the 2003 Tokachi-oki event. First, we briefly describe numerical simulation technique which takes into account compressibility of water column and realistic dynamic bottom deformation resulting from fracturing of the seismic fault. After that, a comparison between the in-situ measured and synthetic bottom pressure variations is carried out. It is worth noting here that our first attempts to reproduce bottom pressure variations during the Tokachi-oki earthquake were not successful enough [Inoue et al., 2010; Nosov and Kolesov, 2007] due to lack of realistic time-spatial law of bottom deformation and vertical resolution in numerical modelling. In this joint study, for the first time we succeeded to reproduce some features of the in-situ measured signals more or less correctly.
Conclusions summarize general results of this study. In the Appendixes, descriptions of simulation techniques and mathematical models are provided. In particular, in the Appendix B the applicability of linear model for simulation of compressible water column in tsunami source is substantiated.
2 In-Situ Observations
JAMSTEC is operating some offshore cabled observatories in the seismogenic zone in Japan, of which one system has been deployed off Hokkaido in 1999 [Hirata et al., 2002]. The offshore cabled observatory off Hokkaido consists of two ocean-bottom pressure gauges (PGs) and three ocean-bottom seismometers (OBSs), and their dataset are sent to the JAMSTEC in real-time. The 2003 Tokachi-oki earthquake of JMA (Japan Meteorological Agency) magnitude (Mj) 8.0 occurred off the south-eastern coast of Tokachi district at 0450:06 JST on 26 September, 2003. In this section, we consider in-situ pressure and seismic waveforms recorded during the 2003 Tokachi-oki earthquake. The datasets under consideration are the first records ever obtained in a tsunami source.
PGs measure high frequency quartz oscillations, and then they are converted to the physical value of hydro-pressure with compensating of thermal effect. In order to improve
signal-to-noise ratio, original 10 Hz dataset is usually converted to 1 Hz dataset and then processed. In the present study, large signal-to-noise ratio during the 2003 Tokachi-oki earthquake makes it possible to process original 10 Hz dataset.
OBSs measure three components of bottom acceleration at sampling frequency of 100 Hz. Before processing, we convert acceleration to East-West (EW), North-South (NS), and Up-Down (UD) directions. In what follows we consider the UD acceleration component registered by OBS1 and OBS3 which are located close to PG1 and PG2, respectively. Relevant distances are of 4.0 and 3.6 km, whereas distance between PG1 and PG2 is of 68 km. Location of the sensors under consideration is shown in Figure 1.
Time-series of bottom pressure and UD acceleration since 0450:00 JST on 26 September 2003 and relevant spectrograms are represented in Figure 2 and Figure 3. In order to make time-varying spectral representation (spectrograms) we use the following technique:
1. The frequency domain is divided into a number of equal (in logarithmic scale) non-overlapping subsections;
2. Band-pass filtering associated with each frequency subsection is applied to the entire dataset;
3. Envelope of the filtered waveform for each frequency subsection gives absolute amplitude which is used to make spectrogram;
4. Finally, each spectrogram is normalized by dividing by its maximum value. The 300 s time series of bottom pressure and UD acceleration were also used to calculate the power spectra which are shown in Figure 4.
Before interpretation of the observed data, let us introduce characteristic frequencies for compressible water column [Levin and Nosov, 2008; Nosov and Kolesov, 2010]. Being a natural quarter-wave resonator, the water column exhibits a discrete set of normal frequencies:
c(1 + 2n)
4H ’
(1)
where c is the sound velocity in water, H is the ocean depth and n = 0, 1, 2, 3.... For typical depth of H = 4 km the minimal normal frequency v0 = c/4H ~ 0.1 Hz.
An additional characteristic frequency is related to gravitational waves: fg = \J~gJH, where g is the acceleration due to gravity. In the conditions of the planet Earth, at any point of the World Ocean, the characteristic frequency fg is always smaller than the minimal normal frequency v0. Frequencies fg and v0 determine three frequency bands: “Gravitational waves”, “Forced oscillations” and “Hydroacoustic waves”. Generally, bottom oscillations of any frequency give rise to forced oscillations of water column in the vicinity of tsunami source. Bottom oscillations of low-frequency (f < fg), additionally to the forced oscillations, generate gravitational waves which are radiated from the source. Bottom oscillations of high-frequency (f > v0), additionally to the forced oscillations, radiate hydroacoustic waves. In the case of bottom oscillations of intermediate-frequency (fg < f < v0),
=
n
Figure 1. Location map and calculation areas. Location of ocean-bottom seismometers (OBSs) and ocean-bottom pressure gauges (PGs). The vertical co-seismic bottom deformation is shown by isolines at 0.1-m intervals. The horizontal co-seismic bottom deformation is shown by black arrows.
the water column in tsunami source is involved in the forced oscillations only, whereas neither gravitational nor hydroacoustic waves are radiated from the source.
Within the frequency band “Forced oscillations”, according to Newton’s second law, the following relation between the UD component of bottom acceleration and variations of bottom pressure p must be fulfilled:
p = pHaL
(2)
duce an additional set of the normal frequencies for such a two-layered system which is determined from the following transcendent equation [Nosov et al., 2007]:
tanf 22ninnK 1 tanf2nYnHs
ps cs
Pc
(3)
where p is the density of water. So, within the frequency band fg < f < vo, if PG and OBS sensors are located close to each other, spectra of bottom pressure and bottom acceleration represented in pressure units with use of formula (2) are expected to coincide.
The acoustical basement of the continental slope southeast of Hokkaido is covered with a rather thick sedimentary layer (up to 2 km). Insignificant difference between the acoustic stiffness of water and sediments gives us ground to expect coupled oscillations of water column and sedimentary layer. The above reasoning makes it possible to intro-
where Hs is the thickness of the sedimentary layer, cs is the velocity of elastic longitudinal waves in the sediment, and ps is the density of the sediment. In the case of Hs = 0, the set of frequencies described by the equations (1) and (3) naturally coincides.
In order to calculate values of characteristic frequencies Vn, Yn, and fg the following dataset was used: depth averaged value of the sound velocity in water for September c = 1475 m s-1 (Generalized Digital Environment Model, GDEM-V 3.0); depth of ocean at the locations of the pressure sensors Hpg1 ~ Hpg2 = 2230 m Japan Oceanographic Data Center (JODC) bathymetry, http://www.jodc.go.jp/; water density p = 1030 kg m-3; sediment features: cs = 1740 m s-1, ps = 1816 kg m-3, and Hs = 1000 m
c
c
Figure 2. The bottom pressure recorded during the 2003 Tokachi-oki earthquake by PGs: original time-series since 0450 JST on 26 September 2003 and relevant time-varying spectral representations.
(http://igppweb.ucsd.edu/ gabi/). Numerical values of the frequencies (in Hz) are as follows: v0 = 0.165, v1 = 0.495,
V2
0.825, Y0 = 0.138, Y1 = 0.375, Y2 = 0.572, fg
0.0663.
Figure 2 presents the change in time of the pressure near the ocean bottom, registered by sensors PG1 and PG2 and relevant spectrograms. Before processing, the data were reduced to the zero level by means of subtracting the linear trend. The range of pressure variations, calculated as pmax — Pmin, amounted to 490 kPa at sensor PG1 and ~ 518 kPa at sensor PG2. As for the spectrograms, continuous signals are observed in the frequency band nearly v0 (resolution of the spectrograms doesn’t allow distinguishing v0 and Y0). This is the manifestation of elastic oscillations of water column produced by the dynamic deformation of bottom in the tsunami source. We shall call this effect the
acoustic resonance. From Figure 2 it is seen that the acoustic resonance essentially dominates in the bottom pressure variations. Moreover, the acoustic resonance continues for a long time, whereas the other frequency components of the signal decay relatively soon. At the same time, manifestations of acoustic resonance at other normal frequencies (n = 1, 2, 3.. .) are not observed.
The UD components of bottom acceleration registered by the seismometers OBS1 and OBS3 are shown in Figure 3. The range of acceleration, calculated as amax — amin, amounted to ~ 2.23 m s-2 at sensor OBS1 and ~ 2.20 m s-2 at sensor OBS3. One of the most striking features of spectrograms presented in Figure 3 consists in the manifestation of the long-lasting acoustic resonance at the same frequency v0 in the bottom acceleration records. This feature pro-
Figure 3. The UD acceleration of bottom recorded during the 2003 Tokachi-oki earthquake by OBSs: original time-series since 0450 JST on 26 September 2003 and relevant time-varying spectral representations.
vides clear evidence in favour of considering the fully coupling model as the most adequate approach in description of dynamical processes in the tsunami source.
Along with the spectrograms, we consider power spectra of bottom pressure and UD acceleration converted to pressure units by means of formula (2). The spectra are shown in Figure 4. Vertical dashed lines stand for the positions of characteristic frequencies vn, Yn, and fg. Spectra of bottom pressure (blue curves) exhibit main maxima at approx. 0.15 Hz. Both maxima are located between theoretically predicted values of Y0, and v0. The positions of the maxima are slightly different along the frequency scale. This is in favour of the assumption that the main maxima are not related to the spectral characteristics of the seismic source, but to the
resonance response of the compressible water column coupled with sediment layers at the minimal normal frequency. As for the other normal frequencies (n = 1, 2, 3 ...), none of them are clearly evident in the spectra.
As for the spectra of UD acceleration depicted in Figure 4 by red curves, the main maxima are attributed to the frequency domain “Hydroacoustic waves”. Thus, acoustic waves must be generated in water column anyway. Spectra of bottom acceleration also exhibit maxima at the minimal normal frequency (0.15 Hz). However, the most striking fact here consist in nearly ideal coincidence of the spectra of pressure and acceleration represented in pressure units with use of formula (2) within the frequency band “Forced oscillations” (between fg and Yo): the red curves perfectly follow
Figure 4. Frequency spectra of bottom pressure variations and UD acceleration of bottom recorded during the 2003 Tokachi-oki earthquake: (a) PG1 and OBS1, (b) PG2 and OBS3. Red and blue curves stand for OBSs and PGs, respectively. The UD acceleration is shown in pressure units according to the formula (2). Vertical dashed lines stand for the positions of characteristic frequencies: gravitational wave frequency fg, discrete set of normal frequencies of elastic oscillations of compressible water column vn, discrete set of normal frequencies of coupling elastic oscillations of compressible water column and sediments layer Yn.
the blue ones. Minor differences between the curves certainly arise due to finite distances between the relevant PGs and OBSs. Note that outside the frequency band “Forced oscillations”, where hydroacoustic or gravitational waves contribute to bottom pressure, the spectra turn out to be essentially different.
3 Numerical Simulation
In the present technique, a total system consisting of the water column and the underlying ground is assumed to be a weakly coupled system. On this assumption, earthquake ground motion due to seismic fault rupturing is first simulated by the boundary element method (BEM) [Kataoka, 1996]. Then, water column disturbance including acoustic and gravitational waves resulting from the seismic ground motion is simulated using the velocity of dynamic bottom deformation as an input to the water column. Details of this method are described in the paper [Ohmachi et al., 2001].
3.1 Dynamic Bottom Deformations
A number of fault models have been proposed for the 2003 Tokachi-oki earthquake [Fukuyama et al., 2009]. Among all these models, we adopt the fault model proposed by Koketsu et al. [2004] which is based on inversion of strong motion and geodetic data.
Ground model is treated as an elastic half-space. Calculation area for simulation of the ground motion is shown in Figure 1. It extends from 141.01°E to 145.84°E in longitude and from 40.51°N to 43.66°N in latitude (400 x 350 km). In BEM simulation (detailed description of the method can be found in Appendix A), the grid size and grid spacing were 81 x 71 nodes and 5 x 5 km, respectively. BEM implies calculations in the frequency domain. The maximum frequency of the simulation was set at 0.5 Hz. The sampling frequency of the dynamic bottom deformation data - input data for hydrodynamic simulations - is of 4 Hz. Dynamic bottom deformation dataset represent three (EW, NS, and UD components) 3D matrixes each of which contains 81 x 71 x 512 values of bottom velocity.
The fault projection is shown in Figure 5. Employed parameters of the fault and ground are listed in Table 1. The fault is considered as a heterogeneous dislocation which includes 120 subfaults. Each subfault has its own values of slip and rake. Rupture starts at the south-eastern edge (yellow star) and propagates with velocity of 3 km s-1. White dashed lines stand for the positions of rupture front. For each subfault we assume the rise time of 2 s [Yagi, 2004]. The total duration time of the fault rupturing is 35.7 s. Dynamic bottom deformation is calculated till the moment of 128 s when all the motions halt. Co-seismic bottom deformation is depicted in Figure 1 by contour lines (vertical component) and by vectors (horizontal component). Subsidence occurs near the coastal region whereas uplift can be seen offshore where the JAMSTEC sensors are located.
3.2 Compressible Water Column
We consider water column as a layer of an ideal compressible homogeneous liquid of variable depth in the field of gravity. To find acoustic-gravitational wave perturbations of water column in the case of dynamic deformation of the basin floor we consider the problem with respect to the velocity potential. Relevant mathematical model is presented in Appendix B. Numerical method had been described in details in paper [Nosov and Kolesov, 2007].
As for the mathematical model, in this section, we only introduce innovations which concern the bottom boundary condition. Instead of simplified condition that had been previously used in [Nosov and Kolesov, 2007] we consider the following general formula:
dF f ,
dn _(u’ n)’
(4)
Depth 23 km
Dip 2 o o
Strike 230°
Length 100 km
Width 120 km
Dislocation max_7.1 m, mean=3.0 m
P-waves velocity 8.0 km s- 1
S-waves velocity 4.6 km s- 1
Rupture velocity 3.0 km s- 1
Density 3.3 g cm- 3
Rise time 2s
Figure 5. The fault projection and slip distribution estimated in [Koketsu et al., 2004] by strong-motion waveforms and GPS data. Yellow star shows the epicentre of the main shock. White dashed lines stand for the positions of rupture front at 10-s intervals.
coustic waves due to their refraction in the bottom. The framework of the present technique (weakly coupled system) makes it impossible to introduce the bottom elasticity during the deformation process. However, after the time moment when the bottom deformation stops, we can switch from the boundary condition (4) to the following semi-transparent boundary condition:
dF _ dF 1 R - 1 dn dt c R +1 ’
(5)
where F is the velocity potential, u = (uEW, uNS , uUD) is the velocity vector of bottom deformation, n is the unit vector normal to the bottom surface. We assume the deformation of the basin floor takes place during a time interval t, upon which the basin bottom stops.
A realistic formulation of the problem under consideration should, naturally, take into account the elasticity properties of the ocean bottom. Finiteness of the bottom elasticity leads, in particular, to a gradual damping of hydroa-
Table 1. Fault and Ground Parameters
where R is amplitude reflection coefficient at the “water-bottom” boundary for a normally incident hydroacoustic wave. The reflection coefficient is estimated from the following formula [Brekhovskikh and Lysanov, 2003]:
R=
pbcp pc
pbcp + pc
(6)
where pb is density of bottom and cp is P-wave speed in bottom. In assessment of the reflection coefficient we ignore sedimentary layers and consider elastic properties of acoustic basement and water. In the acoustic basement the P-wave velocity is of 8 km s-1 and density is of 3.3 g cm-3, whereas for the water column the sound speed is of 1.475 km s-1 and density is of p _ 1.03 g cm-3. Thus, the reflection coefficient takes the value of R ~ 0.89.
Numerical simulation of behaviour of compressible water column was carried out in the computational domain which is shown in Figure 1 as “Calculation area of tsunami”. It extends from 143°E to 145.84°E in longitude and from 40.51°N to 43.2°N in latitude (235 x 299 km). The explicit Finite Difference scheme (z-levelled model, rectangular grid) was employed as numerical method. The 500 m gridded bathymetric data set of JODC was used to build the computational grid. In calculations the following values of space increments
Figure 6. Comparison of synthetic (red lines) and in-situ measured (black lines) variations of bottom pressure. Blue lines stand for the synthetic UD bottom velocity represented in pressure units in accordance with formula (7). Zero time corresponds to the moment of 0450 JST on September 26, 2003.
were used: Ax _ 336 m, Ay _ 332 m, and Az ~ 103 m. The total number of the grid points was 701 x 901 x 74. The time increment was At ~ 0.035 s.
The low-pass cut-off frequency of the input data (dynamic bottom deformation) is of 0.5 Hz. This frequency corresponds to the shortest hydroacoustic wave length of 3000 m. The space increments and the time increment are surely small enough to perform adequate modelling.
3.3 Results of Numerical Simulation and D iscussion
In this Section, results of numerical experiments will be presented. We shall restrict ourselves to considering synthetic time-series of UD component of bottom velocity and bottom pressure calculated at the positions of PGs. In order to reveal the role of the horizontal components of bottom deformation special series of numerical experiments was un-
dertaken. The contribution of the horizontal components of bottom deformation (EW and NS) in the amplitude of synthetic pressure signal turns out to be minor (~ 1 percent). Moreover, the accounting of the horizontal components has no visible influence on the wave forms of the signal. This is why EW and NS components of bottom velocity are excluded from further analysis. However, in what follows we shall describe results of numerical experiments carried out with account of all three components of bottom deformation.
In Figure 6, synthetic and in-situ recorded variations of bottom pressure are shown by red and black curves, respectively. The computed and in-situ measured waveforms exhibit a noticeable difference. However, amplitude of synthetic signal turns out to be in a reasonable agreement with in-situ measured one. Double amplitude of synthetic signal, calculated as pmax — pmin, amounted to ~ 408 kPa at the position of the sensor PG1 and ~ 557 kPa at the position of the sensor PG2. The relevant in-situ observed values are of 490 kPa and 518 kPa. As for the signal arrival time, our nu-
Figure 7. Comparison of frequency spectra of synthetic (red curves) and in-situ measured (black curves) variations of bottom pressure. Blue curves stand for the frequency spectra of synthetic UD bottom velocity represented in pressure units in accordance with formula (7). Vertical dashed lines stand for the positions of characteristic frequencies of compressible water column (see Figure 4).
merical simulation demonstrates an acceptable coincidence between the model and in-situ observed situation: the front of disturbance reaches PG1 approx. 10 s earlier than PG2. During the 2003 Tokachi-oki earthquake, pressure oscillations of large amplitude lasted approx. 100 s, afterwards, long-lasting low-frequency signal of small amplitude - manifestation of hydroacoustic resonance - was observed. In reality the resonance oscillations feed on seismic bottom motions which can exist for a long time due to coupling between water column and elastic half space. Present technique of numerical simulation (weakly coupled system) doesn’t allow reproducing the long lasting low-frequency tail. The switching to the semi-transparent boundary condition (5) after the time moment when the bottom deformation stops (128 s) leads to an exponential decay of the elastic oscillations of water column.
Hydroacoustic resonance in compressible water column exists under condition of presence of two boundaries, i.e. free surface and bottom. In the case of water column being unbounded above, bottom velocity and pressure variations
are linked according to the following well-known formula [e.g. Landau and Lifshitz, 1987]:
P _ Pcuud . (7)
In order to demonstrate the effectiveness of hydroacoustic resonance, UD bottom velocity, represented in pressure units with use of formula (7), is also shown in Figure 6 by blue curve. It is seen that the double amplitude of UD velocity (in terms of pressure) is of 84 and 130 kPa at PG1 and PG2, respectively. Due to hydroacoustic resonance, the double amplitude of pressure variations turns out to be approx. 4 times larger.
The 300 s time-series of bottom pressure (both synthetic and in-situ) and UD velocity were used to calculate the power spectra which are shown in Figure 7. In the frequency band “Gravitational waves” (f < fg) an acceptable agreement between in-situ and synthetic variations of pressure is observed, i.e. gravitational waves (tsunamis) are described by our simulation technique quite properly. The main max-
ima of the synthetic spectra are located close to the theoretical value of the minimal normal frequency v0. However, the maxima of synthetic signals, as compared with in-situ measured ones, are shifted towards high frequencies. This shift arises due to neglecting of the sediment layer during numerical simulation. A noticeable overestimation of spectral level in the frequency band “Hydroacoustic waves” (f > v0) results from acoustic modes generated in a shallow water areas. In contrast to reality, within the framework of the model of absolutely rigid bottom, such modes propagate down the slope without decay due to the refraction in bottom.
The spectra presented in Figure 7 also demonstrate well the effectiveness of hydroacoustic resonance. Spectral level of synthetic bottom pressure is by 2-3 orders of magnitude larger as compared with the UD velocity converted to pressure units with use of formula (7).
Within the frequency band “Forced oscillations” (fg < f < v0) spectra of synthetic bottom pressure (red curves) and UD bottom velocity (blue curves) turn out to be very close. Moreover, the red curves seem to exactly follow the blue ones. Similar feature is observed during comparative analysis of the in-situ pressure and UD acceleration datasets (Figure 4). According to Newton’s second law, all peculiarities of these spectra certainly must fit each others. As for the level of spectra, in present case, a sort of occasional coincidence takes place. It must be stressed that in contrast to the case presented in Figure 4 we consider the UD velocity uUD of bottom deformation, but not the acceleration aUD. Bearing in mind obvious link between acceleration and velocity, one readily obtains from relation (2) the following formula:
p _ pH2nfuuD. (8)
In the case under consideration (f ~ 0.1 Hz, H _ 2230 m), the combination H 2nf takes the value of ~ 1400 ms-1 which is approximately equal to the sound velocity in water. Thus, formula (8) turns out to be an approximate equivalent of the formula (7).
4 Conclusions
JAMSTEC records during the 2003 Tokachi-oki tsunami-genic earthquake provided unique data on behavior of compressible water column in the tsunami source. Time-frequency analysis confirmed the existence of the hydroacoustic resonance - long-lasting oscillations of water column in the tsunami source at the minimal normal frequency (~ 0.14 Hz), whereas other normal frequencies were not distinctly observed in spectra. The hydroacoustic resonance -is detected while processing time-series of not only bottom pressure but also Up-Down bottom acceleration. More precise definition of the hydroacoustic resonance involves consideration of coupled oscillations of water column and sedimentary layer.
Withinthe “Forced oscillations” frequency band (~ 0.07— 0.14 Hz), nearly ideal coincidence of spectra of bottom pressure variations and UD bottom acceleration represented in
pressure units was revealed. Thus, theoretical representations have found the confirmation in reality.
Present technique of numerical simulation is not capable to reproduce bottom pressure time-series in tsunami source exactly. However, some features of the synthetic signals turn out to be quite close to the in-situ measured ones. We succeeded to reproduce more or less correctly amplitude and arrival time of bottom pressure variations. It signifies consistency of the employed approach. As for the dominating frequency, the maxima of synthetic spectra are noticeably shifted toward the high frequency. Exact reproducing of the dominating frequency surely requires an advanced technique of simulation, i.e. multi-layered fully coupling model.
In recent years, bottom pressure observations, both in the near and far field, are actively involved in tsunami early warning systems. The practical importance of this study can be judged from the fact that without considering the compressibility of water it is difficult to interpret signals recorded by bottom pressure gages, especially in the near field.
Acknowledgments. This study was partly supported by JSPS (Japan Society for the Promotion of Science, project 100423) and RFBR (Russian Foundation for Basic Research, project 10-0592102) under the Japan-Russia Research Cooperative Program. We are grateful to JAMSTEC and JODC for providing the data.
Appendix A. Simulation Technique of Dynamic Bottom Deformation
In this section we provide description of technique employed for calculation of the dynamic bottom deformation resulting from fracturing of a seismic fault in an elastic halfspace. The technique is based on the direct boundary element method (BEM) [Kataoka, 1996; Yamanaka, 2010]. According to [Aki and Richards, 1980; Kennett, 2001] elasto-dynamic equation in frequency domain can be written as
(A + ^)uj,ij + ^ui,jj + pw2ui + bi _ 0, (A1)
where ui denotes xi component of displacement, A and p, are Lame parameters, p is mass density, u is circular frequency, and bi is body force. The summation convention is assumed for repeated spatial indices in all equations shown in this section. A weighted-residual expression of equation (A1) in the domain Q is given by
/ {(A + ^)ufc,jfc(xo) + ^uj,fcfc(xo) + pu2uj(xo) +
Jn
bj (xo^ Uij(x’ x0)dQ(x0) _ 0, (A2)
where Uij (x, x0) is the fundamental solution of displacement, which means the displacement at a field point x in xi direction under un impulsive force acting at a source point x0 in xj direction. When Gauss theorem is applied to equation (A2), the equation is rewritten as
I { (A + P')Uik,jk (x’ x0) + ^ Uij,kk (x’ x0) + n
pu2Uij (X’Xo)} uj (xo)dQ(xo) — l Tij(X’Xo)uj(xo)dr(xo) +
J Uij (X’Xo)Tj (xo)dr(xo) +
/ Uij (X’Xo )bj (xo)dQ(xo) _ 0’
of an incident wave. When dislocation of the fault is considered, dislocation vector along the fault is given by
di(xF) _ u+ (xf) — u- (xf)’
(A9)
where xF is a point on the fault. Then, the term i/’i(x) for a point dislocation source is expressed as
V>i(x) _ Tij(i’iF)dj(if).
(A10)
When inhomogeneity of fault plane is considered, the en-
(A3)
tire response of incident waves consist of summation of those n from each point source. Therefore, for the linearly elastic
where Ti and Tij are a surface force and its fundamental system under consideration, the total response can be ob-
solution, respectively, expressed as
Ti _ {Aum,m$ik + Mui,k + uk,i)} nk
tained by summing all the differential effects from each point source.
Fundamental solution of displacement Uij is derived by (A4) solving equation (A6). Then, Uijis expressed as
Tij — {AUim,m$jk + ^(Uij,k + Uik,j)} nk.
(A5)
A fundamental solution for equation (A1) satisfies the following equation:
Uij _ Uj + u j
(A11)
((3 — 4v)Sij + r,ir,j) (A12)
(A + ^)Uik,jk(X’Xo) + M Uij,kk (X’Xo) + pu Uij (X’Xo)
— Fo $ij 5d (x — xo)’
(A6)
U( _ 4^ (“Jij + ^r,i'r,j)’
(A13)
where Fo is unit of force, Sij and SD are the Kronecker delta and the Dirac delta functions. Substituting equation (A6) in equation (A3), one arrive at the Somigliana equation
— y Tij(X’Xo)uj(xo)dr(xo) +
where v denotes Poisson’s ratio, r _ |x — xo| is the distance and a and are written as
Uij (X’ xo)Tj (xo)dr(xo) +
(A7)
E
rn-1 + £an(n — 1)rn-3
n=3
^ Ea»(n — 1)(n — 3)rn
T n=4
(A14)
(A15)
Uij (X’Xo)bj (xo)dQ(xo)
ui(x)’ x e Q
0’ x e Q .
(—*)n(fc? — kn)
(A16)
On the upper boundary (ocean bottom), r, the following boundary equation is satisfied:
where kT and kL are wave numbers for S-wave and P-wave, respectively. Equation (A12) is fundamental solution of static problem, which is called Kelvin’s solution. In the same way, fundamental solution of traction Tij is expressed as
Cij (x)uj (x) + j Tij (X’Xo)uj (xo)dr(xo) —
T- . ____ T(s) + T
Tij _ Tij + Tij
(r)
(A17)
Uij(X’Xo)Tj (xo)dr(xo) _ V>i(x)’
(A8)
where cij is a free term determined from the shape of the boundary, and i/>i on the right-hand side represents effects
T(s) ____
1
8n(1 — v )r2
{(1 — 2v)(5ij r,n +
r,j ni — r,iUj ) + 3r,i'r,j r,n }
(A18)
1
r
a
n
n!
n
r
j _ [(71 + Y3/2)(^ijr,n + r,jni) +
{- (71 + Y2 + 73) + 73} r,j ni+
2(72 — Y3)r,ir,jr,n] ’
where
71 _
da dr ’
72 _
dr ’
73 _
2^
(A19)
(A20)
dv , . V p
— + (V’ V )v _-------------------------+ g
dt p
and the continuity equation
dp
« +d'vp v _0’
The bottom boundary condition (B4) is a mathematical statement of the following principle of classical hydrodynamic [e.g. Landau and Lifshitz, 1987]: in case of non-viscous fluid, if an impermeable boundary is moving, then the normal component of the fluid velocity must be equal to the velocity of the boundary normal to itself.
Let us assume that before the earthquake the water column is in a state of hydrostatic equilibrium. In this state, hydrostatic distributions of fluid pressure, po (z), and density, po(z), are linked by the following equation:
Vpo
po
g.
(B5)
Appendix B. Mathematical Model of Compressible Water Column
The goal of this section is the construction of a mathematical model describing the motion of the compressible water column in the field of gravity, g, in the case of dynamic deformation of the basin bottom during an earthquake. The water column is limited from above by its free surface and from below by absolutely rigid bottom of arbitrary topography. We shall consider the amplitude of motions of the basin bottom a small quantity as compared to the depth, Ini C H. The origin of the Cartesian reference system 0xyz will be located on the free non-perturbed surface with axis 0z directed vertically upward.
Considering the water as an ideal fluid, we shall base the mathematical model on the Euler equation
Let us further consider motions of the water column generated by the earthquake as a deviation from the hydrostatic equilibrium. Assuming the equilibrium flow velocity to be zero, we substitute the pressure and density, expressed as p _ po + p', p _ po + p', into the Euler and the continuity equations (B1) and (B2). Upon performing elementary transformations we obtain the following set of nonlinear equations:
dv
^ + (v’ V )v^
V P
+
p g
po + p' po + p'
dp
dt
+ div(pov + p' v) _ 0.
(B6)
(B7)
(B1)
(B2)
Taking advantage of the direct measurements in the 2003 Tokachi-oki tsunami source, we shall assess the role of nonlinearity in the equations (B6), (B7). The double amplitude of bottom pressure variations amounted to p' _ 500 kPa (Section 2 of the present paper). The respective amplitude of variations of water density is estimated from the following formula [e.g. Landau and Lifshitz, 1987]:
p'
'2
:p/c
(B8)
where v = (u V’ w) is the fluid velocity vector, p is the fluid pressure, p is the fluid density. The equations (B1) and (B2) need to be supplemented by boundary conditions on the free surface
where c _ \J(dp/dp)ad ~ 1500 m s-1 is the sound velocity in water (differentiation is taken with respect to adiabatic change). According to the formula (B8), density variations amounted to p' ~ 0.2 kg m-3. So, in the problem under consideration, bearing in mind the value of water density of po ~ 1000 kg m-3, it is reasonable to assume
z _ £(X’ y’t) : — + u^x + Vdyy — w _ 0’ p _ patn
and on the bottom surface
z _ —H(X’ y) : Vn _ (u’ n) ’
p C po.
(B9)
(B3)
(B4)
Let us further assess the relative impact of the convective and unsteady acceleration terms in the equation (B6). Considering an acoustic wave of length A and period T, we derive the following estimations:
where patm is the atmospheric pressure, vn is the fluid velocity normal to the bottom surface, u is the velocity vector of bottom deformation, n is the unit vector normal to the bottom surface. Since we restrict ourselves to small bottom deformations, |n| C H, depth H as well as vector n are considered as time-independent functions.
V 2 V 2
(v’ V )v ~ "X _ cT
dv V di ~ T’
(B10)
(B11)
r
where V is a characteristic fluid velocity. By comparison of formulae (B10) and (B11) one readily arrive at the following relation:
(v’ V ) v
dv/dt
V
(B12)
The value of characteristic fluid velocity in the 2003 Tokachi-oki tsunami source can be estimated making use of the following formula [e.g. Landau and Lifshitz, 1987]:
V _ p'/poc.
(B13)
The interaction of waves with the coast is described as the full reflection from a virtual vertical wall established along a certain isobath. General formula (B19) represents mathematical expression of the full-reflection boundary condition.
Free pass boundary conditions for hydroacoustic waves are specified on the ocean-crossing outer borders of the computational domain.
The dynamic pressure p' and the fluid surface displacement ^ are expressed via the potential of the flow velocity by the following formulae [e.g. Landau and Lifshitz, 1987]:
p _ —po
dt ’
1 dF
c _ — gat
(B20)
According to the formula (B13), the characteristic velocity amounted to V ~ 0.3 m s-1. This value is much less than the sound speed in water, therefore one can neglect the convective term in the equation (B6).
Ultimately, the set of nonlinear equations (B6), (B7) is reduced to simpler linear one:
dv
dt
V p' po
+
p g
poc2
dp'
—- + po div v _ 0. dt
(B14)
(B15)
Along with the equations, the free surface boundary condition (B3) can be linearized as well. We shall not describe here this well-known procedure.
The second term in the right-hand side of the equation (B14) is responsible for the buoyancy effects. In hydroacoustics this term is usually neglected. Indeed, the terms in the right-hand side of the equation (B14) can be comparable only the case of very long hydroacoustic waves: A ~ c2/g ~ 230 km. In reality, the wave length is limited by ocean depth. Excluding the buoyancy term from equation (B14), we arrive at the following equation:
dv
dt
V p' po
(B16)
Introducing the potential of the flow velocity, F (v _ V F), one can easily reduce equations (B15), (B16) to the wave equation
HF — c2 A F _0.
(B17)
The wave equation is supplemented by linearized boundary conditions on the free surface
d2F dF
z_0: nt*_ —gdz
and on the bottom surface
dF
z _ —H(X’ y) : — _ (u’ n).
(B18)
(B19)
References
Aki, K., P. G. Richards (1980), Quantitative Seismology, Freeman and Co., New York.
Brekhovskikh, L. M., Yu. P. Lysanov (2003), Fundamentals of Ocean Acoustics. Third edition, Springer.
Chierici, F., L. Pignagnoli, D. Embriaco (2010), Modeling of the hydroacoustic signal and tsunami wave generated by seafloor motion including a porous seabed, J. Geophys. Res., 115, doi:10.1029/2009JC005522.
Ewing, W. M., I. Tolstoy, F. Press (1950), Proposed use of the T phase in tsunami warning systems, Bulletin of the Seismological Society of America, 40, 53—58.
Fujii, Y., K. Satake (2008), Tsunami Sources of the November 2006 and January 2007 Great Kuril Earthquakes, Bulletin of the Seismological Society of America, 98, 3, 1559—1571.
Fukuyama, E., R. Ando, C. Hashimoto, S. Aoi, M. Matsu’ura (2009), A Physics-Based Simulation of the 2003 Tokachi-oki, Japan, Earthquake to Predict Strong Ground Motions, Bulletin of the Seismological Society of America, 99, 6, 3150—3171, doi:10.1785/0120080040.
Gisler, G. R. (2008), Tsunami simulations, Annu. Rev. Fluid Mech., 40, 71-90.
Hirata, K., et al. (2002), Real-time geophysical measurements on the deep seafloor using submarine cable in the southern Kurile subduction zone, IEEE J. of Oceanic Eng., 27, 2, 170-181.
Inoue, S., T. Ohmachi, T. Imai (2010), Verification of a dynamic tsunami simulation technique based on observation records in near-field, Joint Conference Proceedings, 7th International Conference on Urban Earthquake Engineering (7CUEE) & 5th International Conference on Earthquake Engineering (5ICEE), March 3-5, 2010, Tokyo, Institute of Technology, Tokyo, Japan. 1665-1669.
Kajiura, K. (1970), Tsunami source, energy and the directivity of wave radiation, Bull. Earthq. Res. Inst. Tokyo Univ., 48, 835-869.
Kataoka, S. (1996), Development of simulation methods for earthquake motion based on three-dimensional fault-ground models, Ph. D. thesis, Tokyo Institute of Technology.
Kennett, B. L. N. (2001), The Seismic Wavefield, Cambridge University Press.
Koketsu, K., K. Hikima, S. Miyazaki, S. Ide (2004), Joint inversion of strong motion and geodetic data for the source process of the 2003 Tokachi-oki, Hokkaido, earthquake, Earth Planets Space, 56, 3, 329-334.
Kolesov, S., A. Bolshakova, S. Inoue, H. Matsumoto, M. Nosov, T. Ohmachi (2010), Numerical simulation of hydroacuoustic effects in tsunami source, Joint Conference Proceedings, 7th International Conference on Urban Earthquake Engineering (7CUEE) & 5th International Conference on Earthquake Engineering (5ICEE), March 3-5, 2010, Tokyo Institute of Technology, Tokyo, Japan, 1687-1692.
c
z=o
Kowalik, Z., W. Knight, T. Logan, P. Whitmore (2005), Numerical modeling of the global tsunami: Indonesian Tsunami of 26 December 2004, Science of Tsunami Hazards, 23, 1, 40-56.
Landau, L. D., E. M. Lifshitz (1987), Fluid Mechanics, V.6 of Course of Theoretical Physics, 2nd English edition. Revised, Pergamon Press, Oxford-New-York-Beijing-Frankfurt-San Paulo-Sydney-Tokyo-Toronto.
Levin, B. W., M. A. Nosov (2008), Physics of Tsunamis, Springer.
Li, W., H. Yeh, K. Hirata, T. Baba (2009), Ocean-bottom pressure variations during the 2003 Tokachi-Oki Earthquake, Nonlinear Wave Dynamics (Ed: P. Lynett), World Scientific Publishing Co., Singapore, 109-126.
Mikada, H., K. Mitsuzawa, H. Matsumoto, T. Watanabe, S. Morita, R. Otsuka, H. Sugioka, T. Baba, E. Araki, K. Suye-hiro (2006), New discoveries in dynamics of an M8 earthquake - Phenomena and their implications from the 2003 Tokachi-oki Earthquake using a long term monitoring cabled observatory, Tectonophysics, 426, 95-105.
Miyoshi, H. (1954), Generation of the tsunami in compressible water (Part I), J. Oceanogr. Soc. Japan, 10, 1-9.
Nosov, M. A., S. V. Kolesov, A. V. Denisova, A. B. Alekseev, B. W. Levin (2007), On the near-bottom pressure variations in the region of the 2003 Tokachi-Oki tsunami source, Oceanology, 47, 1, 26-32.
Nosov, M. A., S. V. Kolesov, A. V. Ostroukhova, A. B. Alekseev, B. W. Levin (2005), Elastic oscillations of the water layer in a tsunami source, Doklady Earth Sciences, 404, 7, 1097-1100.
Nosov, M. A. (2000), Tsunami generation in a compressible ocean by vertical bottom motions, Izvestiya - Atmospheric and Ocean Physics, 36, 5, 661-669.
Nosov, M. A. (1999), Tsunami generation in compressible ocean, Physics and Chemistry of the Earth, Part B: Hydrology, Oceans and Atmosphere, 24, 5, 437-441.
Nosov, M. A., S. V. Kolesov (2007), Elastic oscillations of water column in the 2003 Tokachi-oki tsunami source: in-situ measurements and 3-D numerical modeling, Natural Hazards and Earth System Sciences, 7, 243-249.
Nosov M. A., S. V. Kolesov (2011), Optimal Initial Conditions for Simulation of Seismotectonic Tsunamis, Pure and Applied Geophysics, 168, 6-7, 1223-1237, doi:10.1007/s00024-010-
0226-6.
Ohmachi, T., S. Inoue (2010), Dynamic tsunami generation process observed in the 2003 Tokachi-oki, Japan, earthquake, Advances in Geosciences - Ocean Science, 18, 159-168.
Ohmachi, T., H. Tsukiyama, H. Matsumoto (2001), Simulation of tsunami induced by dynamic displacement of seabed due to seismic faulting, Bulletin of the Seismological Society of America, 91, 6, 1898-1909.
Okal, E. A., P.-J. Alasset, O. Hyvernaud, F. Schindele (2003), The deficient T waves of tsunami earthquakes, Geophysical Journal International, 152, 2, 416-432.
Panza, G. F., F. Romanelli, T. B. Yanovskaya (2000), Synthetic tsunami mareograms for realistic oceanic models, Geophysical Journal International, 141, 2, 498-508.
Sells, C. C. L. (1965), The effect of a sudden change of shape of the bottom of a slightly compressed ocean, Phil. Trans. Roy. Soc. London (A), 258, 495-528.
Soloviev, S. L., P. S. Voronin, S. I. Voronina (1968), Seismic hydroacoustic data on the T wave (review of the literature) (in Russian), The Tsunami Problem, Nauka, Moscow, 142-173.
Stiassnie, M. (2010), Tsunamis and acoustic-gravity waves from underwater earthquakes, J. Eng. Math. 67, 23-32, doi:10.1007/s10665-009-9323-x.
Titov, V. V., F. I. Gonzalez (1997), Implementation and testing of the Method of Splitting Tsunami (MOST) model, NOAA Technical Memorandum ERL PMEL-112.
Watanabe, T., et al. (2004), Offshore monitoring system records recent earthquake off Japan’s northernmost island, Eos, 85, 2, 14-15.
Yagi, Y. (2004), Source rupture process of the 2003 Tokachi-oki earthquake determined by joint inversion of teleseismic body wave and strong ground motion data, Earth Planets Space, 56, 311-316.
Yamanaka, H. (2010), Evaluation of Earthquake and Tsunami Hazards, (in Japanese), Asakura Publishing Co., Ltd.
A. Bolshakova, S. Kolesov, and M. Nosov, Faculty of Physics, M. V. Lomonosov Moscow State University, Russia. ([email protected])
S. Inoue, Engineering Seismology Section, Takenaka Research
& Development Institute, Tokyo, Japan.
H. Matsumoto, Japan Agency for Marine-Earth Science and Technology, Yokosuka, Japan.
T. Ohmachi, Japan Dam Engineering Center, Tokyo, Japan.
