TSUNAMI SHELF RESONANCE ALONG THE NORTHERN CHILE COAST

This work presents the analysis of resonance in two of the main cities along the northern coast of Chile, where a large tsunamigenic potential remains despite recent earthquakes. By combining a modal analysis solving the theory of free and forced oscillations, with the analysis of background spectra derived from in situ measurements, the spatial and temporal structure of the modes are recovered. Comparison with spectra from three tsunamis of different characteristics shows that the modes found have been excited by past events. Moreover, a difference in the overall response between locations is found. Arica is more sensitive to the characteristics of the tsunami source, whereas Iquique shows a smaller dependency. Results are further contrasted with other methodologies with good agreement. These findings are relevant in characterizing the tsunami hazard in the area, and the methodology can be further extended to other areas along the Chilean coast.


INTRODUCTION
The controlling role of local bathymetry and topography over long waves and tsunamis near the coast has been suggested since early the 20th century (see Rabinovich, 1997, for a review).This controlling role manifests itself as similar spectral characteristics between different tsunamis at the same location, but varying behavior and response for the same event at different locations.Recently, resonant modes, both stationary and progressive, have been invoked to account for the late arrival of devastating tsunami waves observed at several Chilean coastal locations.For example, the third wave was the most devastating in the port of Talcahuano during the tsunami generated by the Maule 2010 earthquake (M w 8.8) (Yamazaki and Cheung, 2011) whereas a wave arriving 90 min after the main earthquake was the largest at Arica during the tsunami generated by the 2014 Pisagua earthquake (M w 8.2) (Catalán et al., 2015).The ubiquity of these bathymetric control effects in areas prone to be affected by tsunami, makes the identification of the resulting resonant modes and their spatial characteristics to be of great significance to characterize and quantify tsunami hazard at coastal communities around the world, be it for hazard mitigation, preparedness or emergency response.This is especially relevant for areas like Northern Chile, where in spite of recent earthquake events, the tsunamigenic earthquake potential remains large (Hayes et al., 2014).
However, identification of the resonant properties requires not only identification of the peak resonant periods, but also it would be desirable to know their spatial characteristics.On the other hand, in analyzing past events it is important to separate source from bathymetric effects.The latter can be also separated in shelf and local scale (i.e., bay) effects.Different methodologies exist to analyze resonance in bays and shelves, which can offer varying degrees of detail.Possible classification criteria of these methodologies include the input data required, the analysis methods used or whether they constitute hindcast or forecast.In what follows, we group them depending on the input data.
The first choice corresponds to cases when actual tsunami data is used in the analysis.Typically, tide gauges or other tsunameters records are used, even in the absence of bathymetric data.This is perhaps the more common approach and it has been used extensively (e.g.Rabinovich, 1997;Vich and Monserrat, 2009;Zaytsev et al., 2016, and references therein).Its main advantage is that no prior assumptions on the tsunami propagation physics are required.However, it can be used only for past events, and typically suffers from poor spatial resolution due to the discrete nature of the measurements.Regarding analysis methods of these type of records, the most common procedure is to use spectral analysis to identify spectral peaks, but more recently the wavelet analysis has also been used to characterize the time evolution of the energy distribution (e.g.Borrero and Greer, 2013).
In order to have a spatial characterization of the tsunami, numerical models can be used to generate synthetic tsunami waveforms over large portions of the ocean, thus bathymetric data are also required.These synthetic input data can also be analyzed by spectral methods to identify the resonant conditions.The results can be coupled and compared with actual tsunami data, thus effectively enhancing the spatial coverage.For example, Munger and Cheung (2008) studied the tsunami generated by the Kuril Islands earthquake (M w 8.3) and found large scale modes extending through the Hawaiian archipelago, and local modes on the island's bays, which were in good agreement with measured run up and inundation patterns.Roeber et al. (2010) analyzed the Samoan Islands tsunami, generated by a M w 8.1 earthquake along the Tonga Kermadec Trench.Their study provides both the amplitudes of the resonant modes and their phase angle, which allows differentiation between modes composed by standing and/or progressive waves.Similar analyses were applied to the more recent Chilean tsunamis of Maule, 2010 andPisagua, 2014, to the analysis of the large scale (i.e.shelf) effects.For the former, Yamazaki and Cheung (2011) identified modes composed of both standing and progressive waves extending roughly through 500 kilometers of coast.They found that the inner shelf (0 to 200 m depth) is the key factor that determines modal characteristics.Catalán et al. (2015) obtained similar results for the 2014 Pisagua tsunami.Specifically, edge waves were identified by means of the modal phase angles and constrasted to analytical models, yielding good agreement.This emphasizes the relation between inner shelf properties to the spatial extension of the modes and resulting resonant frequencies.
Whether real or synthetic data are used, one difficulty lies in separating source from bathymetric events, which in turn, may prevent its generalization to future events unless multiple tsunamis have been recorded at the same location.As an alternative to this, Rabinovich (1997) proposed the use of other long wave data present in tide records to define a background spectra, from which the components dominated by the geomorphological configuration could be identified.This requires sea level data at times when no tsunamis are present, and it provides an event-specific local estimate, lacking spatial coverage.Tsunami source effects can be then identified by comparing these background spectra to tsunami spectra.The method has been subsequently used by Rabinovich et al. (2006); Vich and Monserrat (2009); Abe (2011) and recently by Zaytsev et al. (2016), whom also coupled background spectra from multiple locations to average out local effects and isolate the tsunami source spectra.
Another alternative to establish local effects is to use numerical models to estimate the propagation of synthetic waveforms, for instance, Gaussian humps, over the spatial domain of interest (Tolkova and Power, 2011;Aránguiz, 2015).Once the propagated time series are estimated, Empirical Orthogonal Function (EOF) analysis is used to determine the eigenvalues which represent the resonant periods of the domain.This method allows for isolating local (for instance, at bay scale) from shelf bathymetric effects, and their spatial structure regardless of the tsunami source.However, it may require use of several sources of different shape and size, to minimize synthetic source dependencies.
Finally, there are methods that do not require any tsunami input data or waveform data, but rely solely on the bathymetry.The first approximation considers using hydrodynamic theory for free and forced oscillations to canonical morphologies (see Rabinovich, 2009, for a review).This is the case of early analyses of the bay of Concepcion in Chile (for instance Farreras, 1978).However, actual morphologies tend to depart significantly from these idealized cases.Nevertheless, the underlying basic theory of free and forced oscillations can be extended to more complex geometries by studying the eigenvalue problem associated with the linear long wave equations.Input forcing information is not required to obtain the shelf/bay characteristics, enabling generalization of the results and identification of modes induced by the characteristics of the bathymetry.To solve this modal analysis technique, either finite differences (e.g.Sobey, 2006) or finite element methods (FEM) (e.g.Horrillo et al., 2008) can be used.Further, Horrillo et al. (2008) aims to correlate modes obtained as solution to the eigenvalue problem with spectra from tide gages.One difficulty was that the set of boundary conditions used allowed estimation of standing waves only (see below for details).Bellotti et al. (2012b) expanded the approach by proposing boundary conditions that enabled identification of modes composed by standing and progressive waves, which were later compared to far field tsunami spectra in tide gages within the computational domains.When applied to real cases in Poverty Bay (New Zealand) and Adak Island (Alaska), the proposed scheme yielded modes that were in good agreement with spectra from far field tsunamis (Bellotti et al., 2012a).
Finally, for large portions of the coast the obtained modes will depend on the extension of the domain analyzed, since it is possible that by enlarging the domain various resonators (i.e.several bays) are included.This is relevant, since it requires at least two spatial domains in the analysis, one to study shelf-spread modes, corresponding to a large regional resonator; and another to study local modes of bays and harbors.While results are encouraging, one important issue is the proper selection of the modes with physical significance, among all the eigenvalues obtained.Moreover, it is not clear which of the physical modes will be the most important during a tsunami.Bellotti et al. (2012a) proposes that the technique should be used in conjunction with another approach to asses the physical validity of the obtained modes.
The Chilean coast is an important site to study its resonant characteristics with these methodologies.Not only resonant conditions have been identified to be of importance during recent tsunami events, but the seismic and tsunamigenic potential remains very large, for instance, in northern Chile (Hayes et al., 2014).Consequently, identifying the possible modes that could be excited during a tsunami is relevant to asses which coastal communities could be affected.Consequently, the focus of this work is to identify the resonant conditions in the northern Chile coast, both in terms of frequency and spatial structure, as an important step to improve tsunami hazard assessment in the country.
This work is structured as follows.Section 2 offers a brief review of the spectral method of Rabinovich (1997), the free oscillation problem (Bellotti et al., 2012b) and the criteria to contrast results from both methodologies.Data sources and application to areas of interest are described in Section 3. Next, a discussion on the applicability of the results to forecast tsunami resonance is presented.Conclusions of this study are found in section .

METHODOLOGY
Background Spectra Rabinovich (1997) proposed a methodology to discern between spectral peaks related to either source or site effects, by studying spectra on tide gauges located at the coast.The tide gauge spectra are assumed to be composed by the linear sum of contributions by tsunami (t), background oscillations (b) and instrumental error ( ).
where the subscript c denotes spectra at the coast.The last component is usually regarded as being orders of magnitude smaller than the other two, and is thus discarded.Moreover, if a measurement without a tsunami signal is analyzed, then which means that spectral peaks thus obtained could be associated with local resonant periods.When the tsunami is present, spectral peaks would include both local and source effects.
The main hypothesis is that tsunami and background coastal observations are the result of a transformation of the deep water long waves where Y, P and Q are transfer functions related to the propagation over the ocean, shelf and bays respectively.Since the background spectra in deep water is known to not have any distinctive energy peaks, it is proposed that the characteristics (peaks) of coastal spectra are caused entirely by the transfer functions, which would act as frequency filters on the tsunami source signal (Rabinovich, 1997).However, only the combined effect of these functions can be observed at the coast, since their individual characteristics are unknown.

Modal Analysis
Since the previous method lacks spatial resolution and can be sensitive to the location of the sensor, results using this approach are contrasted with those obtained by a modal analysis method.It consists in the estimation of the normal modes of the system, by solving the eigenvalues associated to the homogeneous Linear Shallow Water Equations (LSWE), given by where η is the free surface displacement, g is the acceleration due to gravity and h the water depth.A typical solution for this problem is given by using separation of variables where α and β represent the spatial and temporal dependencies, respectively.
In the method, selection of appropriate boundary conditions is crucial.At the shoreline, usually a no flux condition of the form ∂η ∂n shoreline = 0, ( 6) is used, where n corresponds to the normal vector pointing outside the domain.
For the boundary open to the ocean, a zero displacement condition (η = 0) is usually employed (e.g.Horrillo et al., 2008).However, the disadvantage of this approach is that by fixing the sea surface deformation, it defines the outer boundary as a nodal line.Therefore, the maximum possible wavelength becomes fixed, and results consists of stationary modes only, which depend directly on the domain size.Bellotti et al. (2012b) introduced the use of an approximated radiation condition for the open ocean boundary in the eigenvalue problem formulation, provided that this boundary is defined as a semi circle of radius R whose center is located in the source of the scattered/radiating waves.This condition corresponds to where C w = gh is the wave celerity; and R is the characteristic length (radius) of the domain.Eq. 7 prescribes the solution at the boundary so that waves can pass trough it.This approach does not restrict wavelengths and can produce modes composed by both stationary and progressive waves.It also allows estimation of the damping factor of each mode given by radiation through the offshore boundary, providing a better representation of the physical phenomena.Combining Eqs. 4 and 5, with the boundary conditions given by Eqs. 6 and 7, it is possible to find a weak formulation of the LSWE.To this end, the Galerkin method (He and Fu, 2001) is applied to the spatial function, so that α(x, y) can be approximated by {X}, representing the spatial distribution of the amplitude of the solution.Upon discretization of the domain, {X} is a vector of N elements, corresponding to the number of nodes.After some algebra, it is possible to recast Eq. 4 as Assuming a canonical time decay temporal component β(t) = e st , the following linear-time independent system is obtained Eq. 9 is similar to structural dynamic systems in that mass [M], damping [C] and stiffness [K] matrices are obtained (He and Fu, 2001;Bellotti et al., 2012b).This expression can be solved as a quadratic eigenvalue problem, whose solutions are composed by a set of complex eigenvalues s and their respective complex eigenvectors {X}, representing oscillations that decay in time.(See Bellotti et al., 2012b, for details on the formulation of the problem.).
Selection of the domain size is a relevant aspect to consider if modes with different scales are expected within the studied zone, making it necessary to explore more than one domain in order to identify relevant results (Bellotti et al., 2012a).For example, Bellotti et al. (2012a) found that the first six long-period modes were relevant, although it is noteworthy that this could have been facilitated by a relatively simple bay shape, that closely resembled an enclosed semi-circle.In a more general situation, the coast will be composed of multiple bays of different scales, thus eigenmodes may include the behavior of multiple bays at once, hampering identification of local characteristics.As it will be seen below, this effect is reduced by using a multi-scale approach.
Modal analysis can provide information about the spatial shape of the resonant modes, but because results can encompass the effect of several resonators simultaneously, its results must be contrasted with the local resonant modes found by the background spectra method.However, the location of the sensor can play an important role.For instance, if the sensor is located in a zone where no significant amplification occurs for a given oscillation mode, no noticeable peaks would be present in the background spectra.An ideal case would be to have a dense array of sensors, which would enable to both identify all the modes that can occur in a particular zone, and compare spectral amplitudes of measured time series with modal amplitudes, thus making the latter results more reliable.However, this is seldom the case.

Methodology
It is thus desirable to combine both approaches systematically, to have an improved representation of the whole coast.Therefore, eigenvalues obtained by the modal analysis are contrasted against periods corresponding to significant energy peaks in background spectra of a particular tide gauge.
To account for resonators of different scales, a two-scale approach is used in the modal analysis.First, a large, regional scale, domain is defined, to study modes extending beyond the continental shelf while retaining local morphological features.Next, a second, smaller size domain is employed to effectively isolate a particular bay, to identify local effects.
Considering that eigenvectors lack physical dimension, it is required to establish a standardized procedure to compare among modes.To this end, eigenvectors are normalized by the mass matrix, and taking their absolute value (Bellotti et al., 2012b) as This enables a qualitative comparison of the spatial distribution of amplitudes among different modes.However, the range of periods to be resolved depends on the domain size (i.e.radius).Moreover, at the shorter periods, the solution becomes more noisy (i.e.many modes with close numerical values are found), making identification difficult.Bellotti et al. (2012b) recommends to consider only a few modes taken from those with the lowest eigenvalues.It is expected that the longest periods would be related with large topobathymetric features within the domain.For instance, the continental shelf and bays, for large and small domains, respectively.
The obtained eigenvalues are then compared with periods associated with peaks on the background spectra.If the value of the periods obtained by these two approaches differ in no more than 10% in period (≈ 1-2 minutes), they are considered equivalent.

APPLICATION TO NORTHERN CHILE Area of Interest and Data Sources
The zone studied corresponds to the Northern Chile coast, owing to its history of tsunami events and the sustained tsunamigenic potential.The focus of the study is set on two locations, namely Arica (18 • 29' S, 70 • 20' W) and Iquique (20 • 13' S, 71 • 10' W), see Figure 1.The former has been hit by several destructive tsunamis in the past.The last significant events were the tsunamis of 1868 and 1877, which destroyed the city (Lomnitz, 2004).Arica is located on the Peninsula Alacran, placed just at the corner where the South American coast changes its orientation, and it is flanked to the north by the bay of Arica, and to the south by Ensenada San Martin.To the south, the coast runs with a N-S orientation, and it can be considered as relatively straight down to Pisagua (19 • 36 S, 70 • 13' W).Here, the inner shelf can be characterized by a nearly uniform slope between 0 and 200 m depth.To the north, the coast runs on a NNW orientation, characteristic of the Peruvian coast.In this section, the coastline can also be considered as a straight line with a uniform slope up to Punta Coles near Ilo (17 • 43' S, 71 • 22' W), in Peru.The presence of surrounding bays and long, straight sections of coast with uniform slope allows for the development of standing waves and edge waves.For example, during the 2014 Pisagua tsunami, the largest wave arrived 90 min after the earthquake, and it was significantly larger than preceding waves.Catalán et al. (2015) identified this as the arrival of edge waves propagating northward between Arica and Pisagua.
Iquique area is considered relevant owing to its economic importance, and because of tsunami activity observed during the 2014 Pisagua tsunami, when waves of relatively constant maximum amplitude were observed during the first 5 hours of the event (about 2 m.).This is also considered indicative of resonant conditions (Catalán et al., 2015).Iquique is located on a relatively straight stretch of coastline, but unlike the section close to Arica, it has a wider shelf and the coastline is interspersed by small peninsulas and bays.For this study, the most relevant section of the coast is that comprised between Punta Pichalo, just south of Pisagua, and Punta Gruesa, south of Iquique (20 • 22' S, 70 • 11' W).Iquique itself is located around the Peninsula Serrano, and it is flanked by the small bay of Iquique to the north and Bahia Chiquinata, a 20 km-long bay to the south.
Bathymetry and topography for the area of interest were integrated from three sources.Regional scale analyses considered 30" GEBCO and 1' ETOPO data, and variable resolution data acquired from nautical charts elaborated by SHOA (Hydrographic and Oceanographic Service of the Navy, by its Spanish acronym).On the other hand, detailed local domains are constructed solely from nautical charts data.All the topo-bathymetric information was transformed to UTM coordinates.The coarser bathymetric grid comprises an area ranging from Southern Peru to Mejillones bay in Chile, interpolated to a resolution of 300 m in both directions, whereas in the detailed grids, resolution was set to 60 m.These grids are then interpolated to triangular elements for the finite element numerical implementation of the modal analysis.
In order to study the resonant characteristics of these two locations under a wide set of conditions and to calculate the background spectra, sea level records from may 2014 to may 2015 were used, typically collected with a sampling rate of 1 min.The year-long record was divided into non overlapping 27-hour segments, which were detrended and detided using the UTIDE package (Codiga, 2011).In order to minimize the effect of atmospheric processes such as storms contaminating the analysis, an arbitrary energy variance threshold was set to 0.04 [m 2 ] and applied to the residual time series.Series with energy below this threshold were kept, yielding 35 events for Arica, and 94 for Iquique.Other threshold values were tested showing no significant influence on the results.It is important to highlight that no tsunami activity was registered on the tide gauges during the selected interval.The time series length was chosen based on to the elapsed time between the 2014 Pisagua Earthquake (M w 8.2) main event and smaller tsunami generated by a significant aftershock (M w 7.7) and was the maximum record length possible for this event.All the analysis used the same record length in order to use the same spectral resolution.
A final data set consisted of the tsunami time series from the Maule 2010, Tohoku 2011 and Pisagua 2014 events, at the tide gauges.These time series were processed in the same way as the tide gauge data, with the exception that a single spectrum, without any variance threshold, was estimated for each one of the 27 hour tsunami segments.

Tsunami and Background Spectra
A FFT with a Hanning window filter with 3 degrees of freedom was applied to each tsunami event and to each 27-hour background record to estimate spectra.The aggregate of the 27-hour tide events were next Bartlett averaged to estimate an average background spectra at each location (Bendat and Piersol, pp:280-283, 2010).The main features of background spectra are illustrated in Figure 2, where the average spectra are shown in red, the 95% confidence interval are shown with dashed lines and individual 27hour tide event spectra are shown in gray.Despite the variability present among individual spectra, some characteristic structure is noticeable which is well recovered by the average spectra.In addition, some clear differences can be observed between Arica and Iquique.Spectral peaks at Arica (top panel of Figure 2) are located within the 55-20 minutes band, and energy decreases outside of this interval.On the other hand, at Iquique main peaks are under 20 minutes, but the overall energy content is relatively constant at all periods.Although visually distinguishable, peaks above T ≈ 33 minutes do not show enough contrast over the baseline energy.Even when the overall behavior of the background spectra differs between the two sites, there are common energy peaks near T ≈ 32 and 26 min.A summary with the resonant periods identified is shown in Table 1.

Modal Analysis
To solve the FEM implementation of the modal analysis, triangular-element meshes were constructed as semi-circles enclosing the coastline.First, a single large scale domain is defined to encompass the whole coastline from Southern Peru to Mejillones bay in Chile, with its center at Pisagua, as shown in Figure 1.Next, two smaller and more detailed domains were defined with centers at Arica and Iquique.Domains were then discretized by means of linear triangular elements, and the eigenvalues were calculated using the QARPACK routine (Betcke and Voss, 2011;Betcke, 2013).To validate the numerical implementation of the method, the modal analysis was applied to the same locations studied by Munger and Cheung (2008), Tolkova andPower (2011), andBellotti et al. (2012a), yielding similar results in spite of the different methodologies used (not shown).

Results
The methodology outlined above was applied to the area of interest.A R=390 km radius domain centered at Pisagua is large enough to encompass both Arica and Iquique, whereas the smaller domains cover local features at Arica and Iquique.Figure 3: Comparison between the identified modes for regional (vertical blue lines) and local (vertical green lines) domains and the Background spectra of Arica (left) and Iquique (right).Model label corresponds to those shown in Figures 4 and 5.
In what follows, B is used to denote regional domain modes; and l a , l i denote local modes for Arica and Iquique, respectively.The final result is shown in Figure 3, where background spectra are shown with the periods identified by the modal analysis as vertical lines.The modal analysis shows a reasonable agreement with relevant peaks in the background spectra, although in some cases, modes coincide with valleys.For example, the second mode (2B) corresponds to a valley in the background spectra of Iquique but is nearby a peak at Arica.To understand this apparent contradiction, the analysis must consider the spatial structure of the modes.Spatial mode amplification maps are shown in Figure 4 for the first six regional modes.These modes typically spread from the coastline to the end of the continental shelf (here assumed to be down to ≈ 200 m depth, similarly to Yamazaki and Cheung (2011) and Catalán et al. (2015)).It can be seen that the cross-shore extension decreases with the period value.In some cases, resonant structures can be found at nearly the whole domain.For example, Mode 1B ( T ≈ 54.55 min) excites energetically both Arica and Iquique, and some indication of a less energy resonant pattern offshore the coast between Iquique and Pisagua is present.Each of these zones extends alongshore for about 80 km and are separated by nodal points that could not be correlated to morphological features.This mode does appear in the background spectrum of Iquique as a low energy mode, but it is not readily visible in the spectrum at Arica.In other cases, modes have comparable periods, but excite different locations.Such is the case of mode 2B (T ≈ 48 min) which shows amplification only at Arica, but not in Iquique.This explains the apparent deficiency found above, and highlights the need to use both modal analysis and background spectra in unison, as well as a multi-scale approach, to identify the role of different resonators.
Mode 3B, T ≈ 42.7 min, corresponds to the second peak in Iquique's background spectrum.Although it does not show clearly in Arica's background spectrum, Arica does show some resonant activity for this mode in some tsunami records (see Discussion for details).In this mode, both Arica and Iquique are confined within resonant zones, while Pisagua becomes again a nodal point.
Both modes 4B and 5B represent a peak at T ≈ 33 minutes in Arica and Iquique's background spectra, which is split into two close modes in the modal analysis, as seen in Figure 3, top panel.The close proximity of these modes and similar spatial structure, suggest that both modes (33.96 and 31.63 [min]) are a result of leakage on the modal analysis by contribution of different resonators, and they could be fused on a single mode.Mode 4B concentrates energy in Iquique, whereas Mode 5B does so at Arica.
Finally, Mode 6B seems to represent both the T≈ 25.6 and the 27.0 min peaks at Iquique and Arica background spectra, respectively.Resonant structures are completely confined within the 200 meters depth zone, and the amplification is similar among the different resonant zones.Although not the focus of this study, it is of note that Pisagua is located at a nodal point for this mode also.This shows that the procedure can also be used to identify zones where resonant activity is consistently small, although this should be addressed with a site-specific analysis.
For both local domains, shorter period modes are found that correlate well with morphological features on each site, as shown in Figure 5.One significant feature is that although the Arica radius is shorter than that of Iquique, the longer mode found has almost twice the period in Arica than in Iquique.This situation will be reviewed in detail in the Discussion.Longer modes at each location spread through the continental shelf and correlate well with the 200 m depth contour.Shorter modes are confined to the main bays within the domains, especially in the case of the Iquique where the coast is interspersed by several peninsulas and wide bays.In contrast, at Arica, the shorter modes show large energy over the relatively uniform stretch of southern Peru coast where the shelf is wider.These modes are also much longer than the longest mode found at Iquique.At Arica, modes 1l a (T ≈ 46.2 min), 2l a (T ≈ 30.6 min), 3l a (T ≈ 25.3 min) and 4l a (T ≈ 22.3 min) have a good agreement with peaks found in the background spectrum.For Iquique, mode 1l i (T ≈ 24.9 min) coincides with the first peak, while the others do not.However, it should be noticed that for mode 2l i the tide gauge is located near a nodal line (Figure 5).Modes 3l i (T ≈ 17.8 min), 4l i (T ≈ 14.9 min) appear to be influenced by local effects in Iquique but also from neighboring bays.This could explain the offset between these periods and those found by the background spectrum.

DISCUSSION
One of the overarching goals of the present work is to identify resonant conditions that could be excited by tsunamis and thereby influence the tsunami hydrodynamics at the coast.Therefore, it is of interest to test the correlation between the modes found in the previous analysis with those present in tsunami records Although each tsunami spectrum has its own characteristics, there is a clear consistency between the results found using both background spectra and modal analysis.It also clear that each location shows a different sensitivity to the tsunami source.At Arica, for instance, all three tsunamis concentrate energy at different frequency bands, but each of them correlates with certain portions of the background spectrum and modal analysis modes..For instance, the Tohoku 2011 spectrum concentrates energy at mid frequencies, centered about T ≈ 33 min.Maule 2010 is present over a wide energy band at long periods, above T = 40 min, and shows a decaying trend toward higher frequencies.Pisagua 2014, on the other hand, shows a more uniform distribution of energy.Despite these differences, the modes found at the regional level show a good degree of correlation with the background spectrum, for instance at T ≈ 48.8 min, with the exception of the Pisagua 2014 spectrum.An improved correlation is found at the local domain, where the mode at T ≈ 22.3 min is found for all three events.This suggests that while the local response is relevant, Arica is sensitive to the characteristics of the source.
In contrast, Iquique shows a more clearly defined spectral structure.Spectral peaks at mid frequencies (T ≈32, 25 and 19, 16 and 13 min) appear in all events.Moreover, the shape of the spectra are relatively similar among all three tsunamis, showing amplification near T ≈ 25 min (mode 5B and 1l i ) and near the peaks found close to modes 3I i and 4I i .While some dependency on the source is still present, it seems that the response at Iquique is more clearly controlled by shelf and local effects.
Additional activities were conducted to test the present methodology.First, results were compared against theoretical solutions for the fundamental period, as given in Rabinovich (2009).The relatively straight coastline at both sites of interests allows applying some simple models.In particular, a rectangular, constant depth basin with reflective borders yields a fundamental period defined as whereas a constant slope profile open to the ocean, yields where in both cases L and h correspond to a characteristic domain length and depth, respectively.Since these expressions are valid for simple symmetric geometries, approximated domains were constructed for both Arica and Iquique, as shown in Figure 7.It is of note that areal coverage and representative depth values determine the resulting period, and their selection could be considered arbitrary.Therefore, to provide a meaningful comparison, modes 1l a and 1l i were used as reference as both modes are well aligned with the 200 m contour.Thus, the larger side on the rectangular domains is chosen to be approximately parallel to both the coastline and this contour.The alongshore extent of the domain was selected to coincide with the alongshore extension of the resonant zone that contains Arica in 1l a , and Iquique in 1l i .Evaluating Eqs.11 and 12, requires determining a characteristic depth.The value h R is taken as the mean depth on the grayed zones in Figure 7, whereas h S represents a closure depth, taken as 200 m at both locations.As shown in Table 2, alongshore lengths (L R ) are similar whereas the shelf width (L S ) is almost 60% larger in Arica.The resulting periods agree well with those determined by the modal analysis (see Table 2).In both locations the alongshore period T R is approximately 50 min, value close to the period in mode 1B.Moreover, period T S in Arica is 44.3 [min] while period 1l a is 46.9 [min], and at Iquique period T S corresponds to 27.8 [min] while period 1l i is equal to 25.6[min].These results suggest that the main features conditioning the value and extension of the modes in the northern Chilean coast are the shelf width defined by the 200 m contour.This is equivalent to establish that a relatively constant shelf slope controls the mode characteristics.As a final comparison, it can be seen that some of the modes found here correspond closely to the ones found by Catalán et al. (2015), where spectra were obtained from numerical modeling data of the 2014 Pisagua tsunami.They focused on what is considered here the regional scale, finding the spatial amplification pattern of several modes.Of those, the longer and shorter mode reported by Catalán et al. (2015) correlate very well with modes 1B and 6B from the present analysis, showing a similar spatial structure and resonant periods (54.5 min versus 57 min for model 1B; and 25.5 min vs 27 min for mode 6B).Although the same bathymetry is used, the method is completely independent from the one used here, and it is considered a further validation of the modes identified here.

CONCLUSIONS
The main goal of the present study was to identify the resonant characteristics of the northern Chilean coast.This is achieved by solving for the modal characteristics of the bathymetry by means of an eigenvalue problem implemented using finite element methods, and comparing the obtained modes with peaks in the background noise spectra from tide gauge data.Afterwards, tsunami spectra enables to asses the response of these sites and the correspondence with the modes found.
Two locations in northern Chile were studied, where two types of modes could be identified.For longer periods, modes can extend through multiple bays in the alongshore direction, while their cross-shore extension is defined by the limit of the continental shelf, defined as the 200 m depth in accordance to previous studies.For these modes, the main conditioner of modal period is the extension of the shelf, or the shelf slope.For shorter periods, resonant zones become associated to coastal bathymetric features such as bays and peninsulas.
It was found that found modes agree well with those occurring in at least one tsunami event and there were modes that occurred for all events.One of the site proved to be more sensitive to the tsunami source, whereas the other site seemed to have a more uniform response.This highlights the controlling effect of morphological features on the tsunami response and also the importance of the source, but the relative importance between these depends on the local morphology.Therefore, it might be relevant to conduct similar assessments along the remainder of the Chilean coast.
The approach used allows identification of resonant zones in the absence of a tsunami event.This is critical to preparedness of coastal communities prone to be affected to tsunamis, in the sense that it would allow identification of hot spots, i.e. places subject to amplification due to several modes.Moreover, cold spots (areas where resonance is absent) can also be found.One shortcoming is that the method does not provide the actual amplification factor, in part owing to the potential dependency on the source characteristics.However, identifying resonant characteristics can have an impact on tsunami hazard protocols.For instance, if it had been known the response of Iquique, the sustained energy levels found during the 2014 event could have been better understood and provide guidance during the emergency.
In future work, the methodology will be extended to the whole Chilean coast, as a mean to improve tsunami hazard assessment by identifying coastal cities and communities that could be subject to tsunami resonance impact.Further and denser measurements are strongly encouraged, even in locations a tide gauge is present already.By having several sensors deployed in a basin, it is possible to gain an insight on the actual structure of background oscillations.

Figure 1 :
Figure 1: Bathymetric domains: Large(left), Arica (center) and Iquique (right).Gray line indicates the open sea boundary for the computational meshes to be used in the modal analysis.

Figure 4 :
Figure 4: Normalized amplitude map for each of the first six modes identified for the regional domain (B).White line indicates the 200 meters depth contour.

Figure 5 :
Figure 5: Normalized amplitude map for local domains at Arica (top) and Iquique (bottom).White line indicates the 200 meters depth curve and green circle indicates the tide gauge location.Note that local domains are zoomed in to improve presentation.

Figure 6 :
Figure 6: Tsunami and backgroud spectra (red) of Arica (left) and Iquique (right).Vertical dashed lines indicate modes identified regional domain (blue) and the local domain (green).

Figure 7 :
Figure 7: Rectangular domains for Arica (left) and Iquique (right) used in the theoretical model calculations.

Table 2 :
Parameters and resulting theoretical periods on rectangular domains.