Characterization of a submarine landslide in Santa María del Oro lake using a numerical model

==CHARACTERIZATION OF A SUBMARINE LANDSLIDE IN SANTA MARÍA DEL ORO LAKE USING A NUMERICAL MODEL==

'''Diego A. Pantoja, Anatoliy Filonov, Noel Gutiérrez 

Departamento de Física  

Universidad de Guadalajara, CUCEI 

Direccion postal 

Guadalajara, Mex. 44430 

e-mail: diego.pantoja@academicos.udg.mx, 

afilonov@prodigy.net.mx, 

nogutier@ucsd.edu  '''

==Abstract==

During summer of 2016 an event was register in the volcanic lake of Santa María del Oro, a perturbation that triggers a considerable rise of amplitude on the lake system for at least two days. Due to the scarce observations in the lake, yet it is not definitive what event took place, a storm, a landslide or other. Based on numerical modelling a series of experiments were carry out in order to determine the kind of forcing that created those anomalies. This study is focused in characterizing a submarine landslide and the effects of the wind valley breeze. The results indicate that the lake acts as a dynamical membrane resonating to frequencies at 3.2 min (external seiche) and in a band of the 2-8 hrs (internal seiches) with the 2.64 hrs the most intense. After every forcing the water  starts to oscillate uniformly and in a clockwise sense until the energy is dissipated by friction at the internal seiches frequencies.

''KeyWords. ''Internal and external Seiches, Crater Lake, Delft3D model.

==INTRODUCTION==

Long standing waves oscillations in enclosed water bodies, or ''seiches'' [1], are forced movements provoked mainly by atmospheric (storms) or topographic (seismic motion) perturbations. These oscillations that can span from a few minutes to several hours strongly depend on the shape of the basin and less in the intensity of the disturbed processes. The modal structure of the membrane-like movements losses its energy by dissipating mechanisms until the system came at rest if the external forcing ceases, [1]. In volcanic lakes those oscillations are typically caused by wind burst or other situation like landslides due to the conic-shape of the basin.  In august 18 of 2006 a device deployment in the Lake Santa Maria del Oro recorded a sudden perturbation on the background ambient configuration for a pair of days. In order to shed some light in what could be happened, the numerical approach was used to simulated the typical behaviour of the lake due to an external rapid forcing.     The lake of Santa Maria del Oro, (21<math>^\circ </math> 22' N and 104<math>^\circ </math> 34' W, Nayarit, México), is a volcanic lake located 750 m over sea level, with nearly a circular shape of <math>\sim </math> 2 km of diameter and a maximum depth of 65 m near the center. The lateral slopes are pronounce, except at the west coast where they end smoothly at a bay, Fig. 1. The atmospheric forcing of the lake comprises a valley breeze circulation due to the mountain chain of 200-300 m high that runs east-west in the northern and southern section of the region. At daylight the wind blows to the west and at night time the wind blows to east [2], Fig.1. The summer stratification  of the lake range from 29 <math>^\circ </math>C on the surface to 22 <math>^\circ </math>C at the bottom. The thermocline (where there is an abrupt jump of temperature) lies around 15 m deep [2], Fig. 1.    <div id='img-1'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|[[Image:draft_Pantoja_832101244-Fig1.png|600px|Bathymetry of Santa María del Oro Lake [m] (upper), temperature profile for summer (right) and typical wind valley breeze (bottom).]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 1:''' Bathymetry of Santa María del Oro Lake [m] (upper), temperature profile for summer (right) and typical wind valley breeze (bottom).
|}
According to [3], a possible landslide took place near the bay at the western section of the lake during august of 2006. The observed perturbations co-oscillated with the background system with a periods between 2 and 6 hrs.

==MODEL==

Delft3D is an open-source numerical model developed by WL/Delft Hydraulics and the Delft University of Technology [4]. It includes implementations of several mathematical models for different physical phenomena (currents, transport, wave propagation, morphological developments, etc.). In the present study, it is used to predict the circulation in a lake system.  The Delft3D model solves the Navier-Stokes equations for an incompressible fluid, under the shallow water and the Boussinesq assumptions using a finite difference scheme (see [5]). The model includes the depth-averaged horizontal momentum equations, (shown here in Cartesian coordinates for the sake of clarity):

<span id="eq-1"></span>
<span id="eq-2"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> \frac{\partial u}{\partial t} + \mathbf{u}\cdot \nabla u + g\frac{\partial \eta }{\partial x} - fv + \frac{u\Vert \mathbf{u}\Vert }{C^2(d+\eta )}  - \frac{F_x}{\rho (d+\eta )} - \nu \nabla ^2 u =0,</math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
|-
| style="text-align: center;" | <math> \frac{\partial v}{\partial t} + \mathbf{u}\cdot \nabla v + g\frac{\partial \eta }{\partial y} - fu + \frac{v\Vert \mathbf{u}\Vert }{C^2(d+\eta )}  - \frac{F_y}{\rho (d+\eta )} - \nu \nabla ^2 v =0, </math>
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
|}
|}

the depth-averaged continuity equation:

<span id="eq-3"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> \frac{\partial \eta }{\partial t} + \frac{\partial (d+\eta )u}{\partial x} + \frac{\partial (d+\eta )v}{\partial y} = Q(d+\eta ), </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
|}

and the vertical momentum equation, which reduces to the hydrostatic pressure relationship via the Boussinesq approximation:

<span id="eq-4"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> \frac{\partial p}{\partial z} = -\rho g. </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
|}

where <math>u,v</math> are the depth-averaged velocity in the <math>x</math> and <math>y</math> directions, <math>C</math> is the Chézy coefficient, <math>d</math> is the water depth, <math>\eta </math> is the free surface elevation above the reference plane (at <math>z = 0</math>), <math>\mathbf{u}</math> is a two dimensional current vector, whose Euclidean norm is <math>\Vert \cdot \Vert </math> , <math>Q</math> are sinks or sources of water, <math>f</math>  is Coriolis force, <math>F_{x,y}</math> is Reynolds stress, <math>g</math> is the gravity, <math>\nu </math> is the horizontal eddy viscosity, <math>p</math> is the pressure and <math>\rho </math> the water density.  In the horizontal direction, the model uses orthogonal curvilinear grids, which support both cartesian and spherical coordinates. Both coordinate systems permit the management of geometrically complex domains. In the vertical direction, Delft3D uses the <math>\sigma </math> or the z coordinate systems. The z coordinate system is made up of layers with predefined depths. In this particular study, the spherical and z-coordinate system are used without boundary conditions.  All the simulations were ran for only six days within a mesh of 156 x 126 number of cell in the (x,y)-directions, respectively, <math>dx \sim dy \sim 16</math> m, with 20 z-layer in the vertical distributed according to [10 10 10 10 6 6 4 4 4 4 4 4 3 2 2 2 3 4 4 4] %  of the water column, where the first value is assigned to the bottom layer and the last value is assigned to the upper  layer (surface). In every cell the number of layers is variable depending on the water depth. The time step used was <math>\Delta t = 3</math> seg.

==MODEL RESULTS==

In order to investigate the possible causes that appears in those of the observational data, three scenarios were set-up for the early stage of the perturbation. A Landslide-forced simulation (Ls), a Wind-forced simulation (Ws), and a Landslide- and Wind-forced simulation (LWs). The Ls is a configuration that characterizes a landslide which comprise an equivalent flux of water displaced by the sediment. Remembering that the intensity of the forcing is not as important as the geometrical shape of the basin, then the input of energy here was chosen rather arbitrary. The Ws is the typical configuration as the lake develops day by day in response of the valley breeze wind, and the LWs is the combination of both forcing.   The numerical event (landslide) for the Ls was configured to happen at hr 1 in the day 2. After the impulsive forcing the water level was set in motion in a back- and forward movement rotating counter-clockwise, Fig. 2.     <div id='img-2'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|[[Image:draft_Pantoja_832101244-Fig2.png|600px|Snapshots of water level six hours later of the simulated landslide. Without considering the first frame, the period expand nearly three hours of simulation. Color represent water level in m.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 2:''' Snapshots of water level six hours later of the simulated landslide. Without considering the first frame, the period expand nearly three hours of simulation. Color represent water level in m.
|}
The system was moving almost by a day-long and then goes back to the state of rest. During the first three hrs the water level was moved with a higher frequency, (Fig 3, upper panels), and then after the third hr, the system damps by 2 orders of magnitude reducing also its frequency considerably, (Fig 3, lower panels). In the east-west section the motion seems more regular and more intense than in the north-south section. However, in both sections the water level shows an opposite behaviour, while in one side the water level is upward and in the other side is downward, resembling the movement show in Fig.2.      <div id='img-3'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|[[Image:draft_Pantoja_832101244-Fig3.png|600px|Time series of water level at four observation points. In the upper panel is shown the early development of the water level at the east-west points (left) and at the  north-south (right) points, note the x-axis is in minutes. At the lower panel are the corresponding sections through the time evolution during the next 19-hours after three hours the initial perturbation. In this case, the time series is in hours. The thick-line represent the eastern and southern points. Note the change in scales in every plot.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 3:''' Time series of water level at four observation points. In the upper panel is shown the early development of the water level at the east-west points (left) and at the  north-south (right) points, note the x-axis is in minutes. At the lower panel are the corresponding sections through the time evolution during the next 19-hours after three hours the initial perturbation. In this case, the time series is in hours. The thick-line represent the eastern and southern points. Note the change in scales in every plot.
|}
With respect to the currents, they behave similar as the water level. The currents are nearly homogeneous in direction but varying cyclonically and lowering its magnitude as the time progress, Fig 4. It also can be seen that in surface layers the circulation is more regular than in lower layers.   <div id='img-4'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|[[Image:draft_Pantoja_832101244-Fig4.png|600px|Horizontal velocity components at the surface layer (top) and at the sixth layer in the bottom of the lake (bottom). In every row, is shown the time series of the components (left), a representative diagram of the vector velocity (center) and the initial snapshot of the horizontal velocity field (right). The thick line in time series represent the north-south component, the vertical lines the initial time marked as  in the first panel whereas the thick vector represents the initial time.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 4:''' Horizontal velocity components at the surface layer (top) and at the sixth layer in the bottom of the lake (bottom). In every row, is shown the time series of the components (left), a representative diagram of the vector velocity (center) and the initial snapshot of the horizontal velocity field (right). The thick line in time series represent the north-south component, the vertical lines the initial time marked as  in the first panel whereas the thick vector represents the initial time.
|}
In order to characterized the whole variability of the lake, the power spectra was computed for the water level time series. This technique was used successfully in others lakes (see ''e.g.'' [6,7]). In Fig. 5, are shown the spectrum of several points in the lake. As can be seen, the motion is concentrated into the frequency of the 3.2 min, and in the frequency band of 2-8 hrs around 2.64 hrs.   <div id='img-5'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|[[Image:draft_Pantoja_832101244-Fig5.png|600px|Power spectra of the water level at the four points in the lake. The upper lines represent the east-west points, whereas the lower lines the north-south points. The vertical line correspond to the 95% confidence interval. The triangles mark the 2.64 hrs and 3.2 min frequencies.]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 5:''' Power spectra of the water level at the four points in the lake. The upper lines represent the east-west points, whereas the lower lines the north-south points. The vertical line correspond to the 95% confidence interval. The triangles mark the 2.64 hrs and 3.2 min frequencies.
|}
In Fig. 6, are shown the power spectra of the simulation Ws and LWs. The interesting fact is that when the lake is forced only with wind (Ws) again the motion is set in the band of the 2-8 hrs (black lines), whereas in the combined simulation appears that both frequencies can coexist (blue lines).    <div id='img-6'></div>
{| class="floating_imageSCP" style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: 100%;max-width: 100%;"
|-
|[[Image:draft_Pantoja_832101244-Fig6.png|600px|Same as Fig. 5, for Ws (Black lines) and LWs (Blue lines). ]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" | '''Figure 6:''' Same as Fig. 5, for Ws (Black lines) and LWs (Blue lines). 
|}

==SUMMARY AND DISCUSSION==

Lake Santa María del Oro can be considered as a circular basin which according to the 2D structure modal configuration, obeys the following standing oscillation from the seiches motion ([1]):

<span id="eq-5"></span>
{| class="formulaSCP" style="width: 100%; text-align: left;" 
|-
| 
{| style="text-align: left; margin:auto;width: 100%;" 
|-
| style="text-align: center;" | <math> \zeta (x,y,t) =  J_s(kr)[A_s cos(s\theta ) +B_s Sin(s\theta )] cos \omega t, </math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
|}

For a basin with radius r, where <math>x=rcos\theta </math>, <math>y=rsin\theta </math> and <math>\theta </math> the polar angle. <math>J_s(kr)</math> is the Bessel function of order s, <math>A_s</math> and <math>B_s</math> are arbitrary constants and <math>s = 0, 1, 2, 3, . . . </math>, with the roots of <math>J_s(ka)</math> given the normal modes of the system. The period of oscillation of the first external model correspond to <math>T = 2\pi / k \sqrt{gH}</math>, where g is the gravity and H the mean depth. If the lake was considered homogeneous in temperature, the first mode of oscillation (external mode) will remain between 2-3 min, whereas for an stratified lake, the first internal mode remains between 2.1-4.5 hrs ([2,8]).   As shown by the combined configuration (Fig. 6), the resonant internal modes can coexist in the same simulation, and in the same band of frequency (2-8 hrs), so there is not a very clear contribution of which forcing caused the sudden perturbation on the observational data.   Both simulations, wind and submarine landslide, can cause resonant amplification on the lake. Then, with the information at hand from the numerical outputs we can not confirm that effectively was a landslide that perturb the lake system, there is more work to do in order to support the idea of a topographical perturbation.

==Acknowledgements==

This work was supported by the Mexican Secretaría de Educación Publica (SEP) through the Program for Teachers PRODEP of DAP. The AF acknowledge the support of CONACyT.

==REFERENCES==

[1] A. B. Rabinovich, "Seiches and harbor oscillations", Handbook of coastal and ocean engineering, (2009).  [2] D. Serrano, A. Filonov and I. Tereshchenko, "Dynamic response to valley breeze circulation in Santa María del Oro, a volcanic lake in México", Geophysical Research Letters, Vol. 29(131649), pp. 27-1-27-4, (2002).   [3] N. Gutierrez, "In-situ measuments of the coupling between internal waves submarine landslides in a volcanic lake", Tesis de Licenciatura en Física (In spanish), Universidad de Guadalajara, CUCEI, Guadalajara, Mexico, (2017).  [4] G.R. Lesser, "Development and validation of a three-dimensional morphological model", Coast. Eng. Vol. 51, pp. 883–915, (2004).  [5] Deltares, "Delft3D-FLOW, User Manual", The Netherlands, (2013).  [6] A. Filonov, I. Tereshchenko, "Thermal lenses and internal solitones in the shallow lake Chapala, Mexico", Chinese Journal of Oceanology and Limnology, Vol. 17, pp. 308-314, (1999).  [7] A. Filonov, I. Tereshchenko, J. Alcocer, "Dynamic response to mountain breeze circulation in Alchichica, a crater lake in Mexico", Geophysical Research Letters, Vol. 33, L07404, (2006).  [8] D. Serrano, "Procesos termodinamicos en el lago volcanico de Santa Maria del Oro, Nayarit", Tesis de Doctorado, UNAM, (2004).