Multidata Study to Evaluate the Impact of Submarine Outfall in a Beach Sedimentary Dynamic: The Case of Samil Beach (Galicia, Spain)

Abstract

The Ria de Vigo (NW Iberian Peninsula) is one of the most impacted coastal areas of Galicia, due to demographic and industrial pressure. One of the main consequences of this pressure is the need to extend the current wastewater treatment plant of the city of Vigo (295,000 inhabitants). This extension includes a new submerged pipeline construction to discharge the treated water in the central channel of the Ria. The new planned pipeline must cross Samil Beach, the most important urban beach of the city. Based on a multitool strategy, this work characterizes the interactions between the new pipeline route alternatives and the sediment dynamics of Samil Beach. This approximation improves the reliability of the results in the subtidal area of the beach, where studies are scarce due to the complexity of the data acquisition. The present study is based on high resolution bathymetry data, seabed physical characterization, a granulometric study of the superficial sediment, and a numerical simulation of the tide, wave climate, and sediment transport in low and high energy conditions using open source Delft3D software. The results showed that the area of interest is a low energy area, which is significantly shielded from wave attack, where fine sand predominates. However, the field data indicated an interaction (accretion-erosion) in the submerged obstacles between 0 and 12 m deep. The model revealed that there is significant sediment movement above a 7.4 m isobath, and that the pipeline would not alter the general transport dynamics of the beach, but would interact in the shallowest section. The main conclusion of this work states that the future structure would not alter the global sediment dynamics of the beach. In addition, in order to guarantee the safety of the new pipeline, it should emerge above an 8 m isobath. The multiapproach methodology presented can be applied to other studies of the interaction between coastal structures and the environment.

1. Introduction

At present, shores worldwide undergo very high pressure due to their inherent touristic interest and industrial development favored by access to the sea, and the effect of climate change [1,2,3,4]. A deficient management of coastal activities has led to many environmental problems, such as serious morphologic coastline alterations, pollution, and the loss of important habitats for coastal ecosystems and the local economy. Therefore, it is essential to have a good knowledge of the coastal environment in order to guarantee an appropriate management to preserve both environment and economic activities [5,6,7].

In this context, the appropriate design of coastal structures is central, since they are indispensable for socioeconomic development, but may cause serious alterations in the coastal environment. According to the Recommendations for Maritime Works collection (RMW), published by Puertos del Estado (http://www.puertos.es/es-es/ROM; Spanish Ministry of Public Works), these structures must simultaneously satisfy three requirements; security, service and exploitation. The security requirement states that the design of the construction must guarantee the physical integrity of the structure itself. The service requirement states that the design must guarantee the service for which the structure is built, and the exploitation requirement states that it must be exploited without serious repercussions in the socioeconomic framework and environment.

From the geotechnical and environmental point of view, one of the requirements in coastal projects is the study of the interaction produced between the structure and the physical environment. These interferences are produced by physical agents (wind, waves, and currents) due to the presence of the structure, and must be predicted in the design; otherwise, it may lead to serious structural problems, and failings in one or more RMW requirements. Accordingly, it is very difficult to predict the interaction between the structure and the local sedimentary dynamics, a process that can have serious consequences both for the environment and the structure itself [8]. From the environmental point of view, the lack of forecast can cause habitat and water quality losses, and permanent morphological alterations [3]. Economically, this lack can also lead to marine (i.e., fisheries) and touristic resource loss, which has a negative social impact. Finally, it can also produce the total failure of the structure functionality. Therefore, it is of great importance to diminish to minimum levels the uncertainties in the forecast of the interaction between the structure and the environment during the design process.

Several methodological approaches have been proposed for the study of the structure—environment interaction [9,10]. Some studies have used photographic interpretation and Geographic Information System (GIS) methods to evaluate the coastal evolution [11]. Other authors used beach and oceanographic parameters monitoring [12,13]. Some studies are based on numerical modeling to study the morphological and sediment transport impact of the coastal structures on beaches [14,15,16,17].

The water-treatment plant E.D.A.R.-Lagares treats the sewage of the Vigo (295,000 inhabitants) coastal city, and empties into the central channel of the Ria de Vigo at a 35 m depth. This plant was built in 1997, with a treatment capacity of 400,000 inhabitants equivalent (IE). Nevertheless, the increase of 14,000 inhabitants in 11 years (1998–2009), and the entry into the force of the E.U. Water Framework Directive necessitated the extension and modernization of the plant. The new treatment station has an 800,000 IE capacity, and treats the storm runoff waters without directly emptying into the Ria during large rain events. This extension implied the construction of a new pipeline larger than the two current ones, which must cross the Samil Beach until the central channel bottom of the Ria is reached, at a 35 m depth.

The main purpose of this work was to determine and characterize the interactions that can be produced between the new pipeline route alternatives and the sediment dynamics of Samil Beach, providing greater reliability to the project. This research was performed using a multitool strategy in order to improve the results, given the complex nature of the environment, and to support the decision-making in the first steps of the construction project, allowing for the mitigation of beach erosion and the protection of valuable infrastructure.

2. Study Area

The Ria de Vigo is the southernmost coastal inlet of the Rias Baixas. It presents a funnel shape in plan view, and its axis is oriented NE–SW (head-mouth). The Ria is 30 km long and 0.6 km wide in its inner part, and 15 km wide at the outer part. This inlet presents rocky islands at the entrance, Cíes Islands, and a central deep channel that reaches a 50 m depth in the southern mouth.

Samil Beach is located in the southern margin of the outer sector of the Ria de Vigo (Figure 1), in a small bay oriented to the NW. The bay is bounded in the SW by the Toralla Island and the large rocky outcrop of Estai Cape-Canido, and in the east by the Mar Cape. Samil Beach is 1.9 km long, occupying half of the bay coast. This beach presents reflective morphodynamic elements and fine sand (209 µm), typical of the middle Ria sector, based on exposure to the open sea waves and sediment provenance [18]. The seabed in the southern margin of the Ria is dominated by sands, with the presence of gravels [19]. Samil Beach is separated by a rocky outcrop from the Vao Beach, and in their submerged sector there are also outcrops presenting magnitudes in the order of tens of meters.

The currents in the Rias Baixas are controlled by the tide, the fluvial discharge and the shelf wind pattern. The tide presents a semidiurnal and mesotidal regime, with a mean range of 2.2–2.4 m [20,21], a mean velocity of 5–10 cm/s and a maximum velocity of up to 30 cm/s in the shallowest and narrowest sections [22]. The fluvial discharges in the Rias Baixas are seasonal; they are at their maximum in winter and minimum in summer [23]. These discharges are in general of minor importance, given the size of their drainage basins (620–3600 km²) compared to the area of the Rias (106–251 km²). The shelf wind action produces the deep-water upwelling events in the coast under the north winds, prevailing between March and October, and it causes a downwelling of superficial waters under south winds, which prevail for the rest of the year. Such currents are subtidal and show a smaller magnitude than the tidal ones (<6 cm/s).

The wave climate in this region also presents strong seasonal behaviour. In summer, the significant wave height (Hs) is 1–1.5 m, with a peak period (Tp) of 8–10 s, whilst the winter shows higher energy conditions, with an Hs of 1.5–3 m and a Tp of 12–14 s [24]. In extreme conditions, the Hs can reach values of 8 m [25,26]. The predominant incoming direction is WNW–NNW, presenting a 76.3% annual frequency. The SW waves are not very frequent in the region, at 8.0%, but they are usually linked to high energy conditions [27]. Inside the Rias, the wave conditions are modified significantly depending on the incoming wave direction due to the presence of the Cíes Islands and the Ria axis orientation. Therefore, the offshore waves undergo lower changes and easily reach the inner zones of the Ria when their direction matches the orientations of the mouths, i.e., the NW and SW waves. In addition, SW wave events coincide with the natural orientation of the Rias reaching the innermost areas. In addition, such wave events are generally related to winter storms, and their effects in the bottom may be important.

Characteristics and Alternative Routes of the Submarine Pipeline

The new pipeline consists of a high-density plastic pipe with an inner diameter of 1.8 m. In the first section of the route, the pipeline is underground, and it emerges gradually on the seabed from a 3.3 m isobath. From this point, the pipeline is protected by a concrete structure with a trapezoidal section, which increases its size gradually with the depth until it reaches a 20 m base width, a 4 m crest width and a height of 5 m from the bottom.

Two possible routes were designed according to the bay morphology (Figure 1). In the southern alternative, the pipeline departed from the southern sector of the bay, specifically from the rocky outcrop located between the Samil and Vao Beaches, with a direction of 300° (with respect to the north clockwise). In the northern alternative, the pipeline departed from the central sector of the bay, close to the current submerged pipelines, with a direction of 340° (with respect to the north clockwise).

3. Methodology

3.1. Field Work

The bathymetry of the subtidal area of Samil Beach was acquired in July 2011 up to a 38 m depth with a Reson SeaBat 8125 multibeam echosounder (Teledyne, Denmark) mounted in the INNDAGA vessel (GEOMA-University of Vigo). This echosounder emits 240 beams, 40 times per second, with an angular resolution of 0.5° in the transversal direction and 1° in the advance direction. The obtained raw field data was corrected to remove errors of a different nature: (a) the GPS positioning error, corrected by the DGPS system (Differential GPS); (b) the error due to the vessel’s movement (heave, pitch and roll), corrected by RTK (Real Time Kynetics) movement and a course sensor; (c) the error due to the tide-produced water level, corrected by pressure sensor field measurements; and (d) the error due to the sound velocity vaRiation produced by the changes in the water column density, corrected by SVP (Sound Velocity Profiler) measurements. The final resolution of the bathymetry was of 1 m in the horizontal plane and of a few centimetres in the vertical plane.

The detailed characterization of the superficial features of the present-day pipelines and the surrounding seabed up to a 30 m depth was performed with a Klein System 3900 (Klein Associates, USA) Side Scan Sonar (SSS) in December 2011. The GPS positioning during the acquisition was corrected with DGPS (Differential GPS). The acquisition and processing were carried out with SonarPro (Klein Associates) software.

Based on previous information of surface sedimentary distribution in the Ria de Vigo [19], the sediment samples were taken (in 16 points) following a regular grid of approximately 400 m using Van Veen grab from the INNDAGA vessel. The positioning was obtained with a GPS system.

3.2. Laboratory Work

The grain size of the sediment samples was obtained by dry sieving for the coarse fraction (>63 µm) every phi. The fine fraction (<63 µm) was separated by wet sieving, and was analyzed using a sedigraph (Micromeritics Sedigraph 5100, EEUU). The sediment grain size distribution was calculated using Gradistat v.7 software for Microsoft Excel [28]. This software analyses the statistically results of the granulometry using the method of moments, and classifies the samples by applying the graphic method [29,30]. The main parameter used in this work is the median (d50).

3.3. Database

The submerged relief of the Ria was characterized previously with bathymetric data of the nautical charts published by the Spanish Navy Hydrographical Institute (IHM). This data covers the Ria de Vigo above a 100 m isobath in the continental shelf with a highly variable spatial resolution between 10 and 200 m.

The wave characteristics in open waters were obtained from the Puertos del Estado database (http://www.puertos.es/es-es/oceanografia/Paginas/portus.aspx). The coastal buoy of Cabo Silleiro (42.10° N, 8.93° W, Figure 1), located 3 km from the southern mouth of the Ria de Vigo, and the node no. 1044069 of the WANA numerical wave model 6 km from the northern mouth (42.25° N, 9.00° W, Figure 1) were selected. The buoy data supplied information about the seasonal, annual and extreme wave regime in terms of significant wave height and peak period, whilst the WANA numerical model supplied modelled data about the significant wave height, peak period and incoming direction of the different wave conditions.

The tidal data was obtained from the tidal gauge located in the harbour of Vigo (42.24° N, 8.73° W, Figure 1). The data used in this work was the water level, with a temporal resolution of one minute, and the harmonic constituents calculated from the historic time series of the gauge.

3.4. Numerical Simulations

Numerical modelling was used to determine the sediment transport under different wave conditions and bathymetry, which was modified in order to include the different route alternatives of the new pipeline. The numerical simulation of tidal currents, waves and sediment transport was performed using the open source version of Delft3D (Deltares, The Netherlands). The software consists of various modules: hydrodynamics (FLOW), waves (WAVE) and sediment transport (MOR). Delft3D computes the simulations by the finite differences method, and the configuration used in this study was in two horizontal dimensions (2DH).

In the hydrodynamic simulations, the program solves the movement equations of the flow, known as the momentum, continuity and conservative transport equations. The wave module (WAVE) is the Simulating Waves Nearshore (SWAN, Delft University of Technology, The Netherlands) [31]. This program solves the wave energy balance equation and takes into account the bottom friction, the current, the wind effect, the wave-wave interaction, the refraction and the reflection. The FLOW and WAVE modules share the information between them; in this way, the wave module takes into account the results of the flow module and vice versa. The morphodynamics module (MOR) uses the results from the previously mentioned modules to compute the sediment transport, both bed and suspended loads, updating in each time step the bottom sediment level and the bathymetry resulting from these processes. The calculations are based mainly on the critical shear stress of the sediment, which depends on the sediment median, d50. This data was obtained from the field campaign samples. The input parameters necessary to simulate the current and the water level were the main harmonic constituents of the tide available in the database of the tidal gauge of the Vigo harbor.

The model was established with a domain covering the entire Ria de Vigo, which is 22.8 km wide and 37.3 km long. A curvilinear grid was applied on this domain that was adapted to the geometry of the Ria, making the orientation of the cells coincide with the orientation of the south and north Ria mouths, respectively, and with the axis of the estuary. The grid cell size varied from 220 m in the open sea to 8 m at Samil Beach, depending on the Ria’s geometry. The grid quality is critical for accurate simulations; therefore, the grid was made sufficiently dense in the area of interest, with orthogonal values of less than 6% in close boundaries and 2% in open boundaries.

The initial bathymetry was based on three different sources: a) the nautical chart of the Ria de Vigo of Instituto Hidrográfico de la Marina, with a resolution varying from 800 to 10 m, inversely proportional to the depth; b) a topographic survey on the emerging areas of Samil and Vao Beaches; and c) a bathymetric survey on the subtidal area of the Samil and Vao Beaches using a multibeam. Spatial resolution of these sources was higher than the resolution of the model grid, so a grid cell sea bottom was obtained by the “nearest point” option.

The model was 2DH, and had three open boundaries of two types: the western (open sea) boundary, located 8.5 km away from the Cíes Islands, at a depth of approximately 100 m, is a water level type, and the northern and southern boundaries are of the Neumann type. The water level type boundary calculates the tide by means of the known harmonic components M2 and S2 (main lunar half-day and main solar half-day, respectively). The Neumann-type boundary calculates the flow from water level gradients. This configuration is appropriate in the case that the model contains transverse boundaries to the coast with an open sea boundary, thus minimizing the generation of false waves within the model. Wave parameters are entered from all three boundaries simultaneously.

The input data for the wave modelling (significant wave height, peak period, and provenance direction) were grouped in different cases. In the selection of the cases, low energy (named summer) and high energy (named winter) conditions were taken into account. In addition, the sensitive incoming directions for the area of interest were also selected. The presence of the Cíes Islands in the entrance produces a shielding effect within the Rias with respect to the incident open sea waves. Therefore, the wave conditions most affecting the area of interest are those that coincide with the natural orientation of the Ria mouths; i.e., SW (225°) and NW (315°). The Hs-Tp correspond to the most frequent combination (low energy case) and to the most frequent considering high energy conditions (assuming high energy when Hs > 2 m and Tp > 14 s). The simulated cases are summarized in

Table 1. Simulated wave cases, including the frequency of each incoming direction and the frequency of the Hs-Tp values in % of time.

Cases	Direction (°)	Hs (m)	Tp (s)	F. Direction %	F. Seasonal %
SW winter	225	2.5	14	7.5	2.5
SW summer	225	1.5	10	7.5	10.5
NW winter	315	2.5	14	50.6	2.5
NW summer	315	1.5	10	50.6	10.5

.

The model simulated 10 days for each wave condition, with a morphodynamic acceleration factor of 5; that is, the morphologic changes in each simulation were equivalent to 50 days. This time was enough to observe the sediment transport effects in the beach, as well as the interaction between the sediment dynamics and the planned structure under different wave conditions. The morphologic response of the beach to each wave condition was analyzed by transversal profiles obtained from the resulting bathymetry and erosion-sedimentation maps of the seabed. The model ran in a unique bathymetric grid, in which the resolution was locally increased to represent the submerged pipeline.

3.5. Birkermeier (1985) Equation of Closure Depth

The theoretical closure depth in the present work was calculated by the equation proposed in [32]. This equation relates the closure depth profile to the most energetic waves:

$$h^{*}=1.75H_{12}-57.9\left(\frac{H_{12}^2}{g{T}_{12}^2}\right)$$

(1)

where H₁₂ is the significant wave height exceeding 12 h per year, T₁₂ is the period linked to H₁₂, and g is the acceleration due to gravity. The data to determine H₁₂ and T₁₂ were obtained from ROM 0.3–91 [24]. The H₁₂ was calculated from the mean significant wave height regime. The period linked to the wave height is defined in ROM 0.3–91 for Area III as

$$T_p=\left(4.5~9.2\right)\sqrt{H_s}$$

(2)

In the above equation, the peak period is proportional to the Hs by a coefficient that may vary within a range of values, highlighting that there is no lineal relationship between these two parameters.

4. Results

4.1. Bathymetry

The submerged area of the Samil Beach is characterized for combining sandy zones with rocky outcrops (Figure 2a). Above a 10 m isobath, the subtidal zone is bounded in the north by the Cabo de Mar rocky outcrop (a) and in the south by the outcrop dividing the Samil and Vao Beaches (b). In the central sector, there are also a group of outcrops (c, d and e) in which the two old submerged pipelines are located. The slope in the southern zone is less pronounced (0.027) than in the northern zone (0.036). Given an isobath of 10–20 m, the seabed is formed mainly by sediment, excepting a series of rocky outcrops (f) located at a 15 m depth. In this area, there is a group of marks in the seabed that could be anthropic, possibly tracks left by some fishing gear. The slope in this area is 0.009. This sector is bounded in the deepest zone by the escarpment (20 m) that limits the central channel of the Ria. In the escarpment, two sections with different slopes are observed—one much steeper (0.063) (h) than the other one (0.032) (g). In the latter section (g), 300 m long ripples with wavelengths of 15 m are observed. The SW zone of Samil Bay is characterized by the rocky outcrop that forms Toralla Island. In the shallowest sector (0–10 m), between Vao Beach and the island, a lobe-shaped double sand bar is developed (i).

In the bathymetry, the old pipelines can be easily distinguished in the seabed. Both pipelines are completely buried in the shallowest area (Figure 2a–c). The occidental pipeline (OCP) emerges at 4.3 m, whilst the oriental pipeline (ORP) remains buried up to a 6 m depth, the ORP arriving to a 25 m isobar and the OCP to a 30 m isobar, already in the central channel.

The resolution of the bathymetry allows one to distinguish some morphologic characteristics linked to the obstacles present at the bottom (rocky outcrops and current pipelines) of the shallowest zone (0–10 m) at Samil Beach (Figure 2b,c). The rocky outcrops present erosion structures 5–10 m long and 0.5 m deep on the SE side, and sediment accumulation in the NW side. In the pipelines, the OCP presents erosion at both sides up to a 12 m depth (Figure 2d–f) that is 4–6 m long and 0.5 m deep. The ORP presents a larger erosion furrow, up to 8 m, with a similar depth, but in the east flank only. This furrow also disappears at a 12 m depth (Figure 2d–f). The magnitude of these furrows parallel to the pipelines decreases gradually as the depth increases. The space between the two pipelines is characterized by an asymmetric accumulation that is developed in the W–E direction (Figure 2f).

4.2. Surficial Grain Size Distribution

The subtidal sector of Samil Bay is dominated by sands (Figure 3). The superficial distribution presents a textural transition in the SW–NE direction from very coarse sand (>1000 µm) to very fine sand (200–100 µm). In the NE corner of the studied area, the presence of mud was detected in low percentages (<10%).

Toralla Island and the submerged outcrops exert important hydrodynamic control on the sediment distribution in the bay. These geographical features produce a shadowing effect, favouring the finest material accumulation on their lee side (<200 µm).

4.3. Side Scan Sonar (SSS)

The bottom characterization by the side scan sonar allowed for an assessment of the current state of the old pipelines and the evaluation of the interaction between the seabed and the structures, in addition to gathering information about the sedimentological characteristics of the seabed.

The shallowest images, at 3.3–7 m deep (Figure 4a), showed the emerging zone of the occidental pipeline (OCP). It can be observed that there is an accumulation in the west flank, at some point burying the pipeline, and erosion in the east flank with parallel accumulation, as was identified in the multibeam echosounder dataset. In some images, this pipeline shows an unearthed concrete basement. In the images at a deeper depth, between 7 and 12 m, the oriental pipeline (ORP) is identified (Figure 4b), and a similar accumulation-erosion pattern in the surroundings of both structures is observed, but of a smaller magnitude than in the shallowest section. At isobaths at least 12 m deep, the interaction between the pipelines and the seabed disappeared (Figure 4c), and the images at a deeper depth did not show interaction marks anymore.

The main sedimentary characteristics of the seabed in this area were well-developed ripple fields due to the wave action. One of those fields was identified 300 m eastward from the pipelines, at an 11 m depth. Another field was observed in the vicinity of the pipelines at a depth between 15.1 and 15.3 m, occupying a 0.41 km² area with a maximum height of 0.2 m.

4.4. Seasonal Dynamics of Samil Beach

The seasonal dynamics of Samil Beach are controlled by the wave action. The wave energy in this area is conditioned by the amplitude and orientation of the Ria mouths, as well as the beach location and orientation. From the wave database, we identified two main incoming directions that directly affect the area of interest. One is the NW direction, the most frequent direction in this area, which coincides with the natural orientation of the northern mouth. This allows these waves to propagate southward inside the Ria and towards Samil Beach, which is exposed because it faces the W-NW direction. The second main incoming direction is SW, less frequent than the previous case, but linked to storms, which coincides with the orientation of the central channel of the Ria. Such waves reach Samil Bay after a refraction processes due to the southern margin coastline in the outer Ria sector and the presence of the Toralla rocky outcrop.

Under NW winter conditions (315°, Hs₀ = 2.5 m) in the outer sector of the Ria, the significant wave height (Hs) diminishes rapidly until it reaches 1 m (

Table 2. Wave condition variation during the wave propagation inside the Ria.

Open Sea		Studied Area		Propagation Coefficient
H (m)	Direction	H (m)	Direction
Winter Conditions
2.5	315° (NW)	0.5–0.8	310°	20–32%
2.5	225° (SW)	0.8–1.0	300°	32–40%
Summer Conditions
1.5	315° (NW)	0.3–0.5		20–33%
1.5	225° (SW)	0.4–0.6		26–40%

). The waves arrive in the area of interest with a 0.5–0.8 m height, presenting a propagation coefficient of 20–32%. The incoming direction changes from 315°, in the open sea, to approximately 310° in the area of interest. Under similar conditions, but in the SW incoming direction (225°), the waves are attenuated as they travel inside the Ria through the central channel, reaching values of 1–1.5 m in the outer sector of the Ria. In the area of interest, they present 0.8–1 m, showing a propagation coefficient of 32–40%. The change in the incoming direction is greater under these conditions due to the refraction along the southern margin—from 225° in the open sea to 300° in the area of interest.

Under summer conditions (Hs₀ = 1.5 m), the Hs values inside the bay are much smaller, although the distribution pattern inside the Ria and the propagation coefficients are similar. Under NW conditions, the values in the bay are 0.3–0.5 m, showing propagation coefficients of 20–33%. Under SW conditions, in the area of interest, Hs values of 0.4–0.6 m are observed, presenting a propagation coefficient of 26–40%.

The study of the morphological response of the beaches to different wave conditions (Table 1) was based on the Delft3D model’s (MOR module) results. The erosion-accretion results under winter conditions (Figure 5) showed that the sediment transport was produced exclusively along a narrow strip 300 m wide in the subtidal sector. In the two cases studied, the main feature was the transport produced in strips parallel to the coastline, accretion in the intertidal and shallow subtidal sector, and erosion in the deepest subtidal sector.

The resulting sedimentation presented values of 0.3–0.5 m, and the erosion presented values of 0.25–0.5 m. Under summer conditions (not shown), the transport pattern was similar, but of substantially smaller magnitude, with an accretion of 0.3–0.35 m, and with an erosion of 0.15–0.35. The location of the accretion and erosion areas, and their magnitude change, was between the NW (Figure 5a) and SW (Figure 5b) wave conditions. The location varies specifically in the deepest subtidal sector, where accretion is produced under SW conditions, and erosion is produced under NW conditions. Based on the bathymetric results of the model, the seasonal changes from the transversal profiles were analyzed; four profiles in Samil Beach, one profile located in the transition zone between the two beaches, and one more profile in Vao Beach (Figure 1c).

In general, the resulting profiles in the sandy area (P1, P2 and P3, see Figure 1), under winter and summer wave conditions, present a higher slope than the original profile (Figure 6). This favours the existence of a cutting point between the original and resulting profile. Thus, in shallower depths than this point, it produces accretion, whilst at deeper depths, it causes erosion in the profile. Comparing the summer and winter conditions, the summer profiles have higher slopes and, consequently, cut the original profile at a shallower depth (between −1 and −2 m) than the winter profiles (between −2 and −4 m). In summer conditions, for both NW and SW waves, the beach response was similar. In winter conditions, the profile response varies slightly, developing steeper slopes with NW waves. The slope of the profiles shows an increasing trend towards the south. In P1 (Figure 6a), the resulting slope varied slightly between 0.040 and 0.047; in P2 (Figure 6b), it varied between 0.038 and 0.070 in winter and between 0.026 and 0.068 in summer; in P3 (Figure 6c), it varied between 0.059 and 0.070 in winter and was 0.050 in summer.

The sedimentation, as well as the erosion, of the original profile also tends to increase towards the south, from P1 to P3. Thus, the accretion values varied between 0.3 and 0.8 m in P1, between 0.6 and 1 m and between 0.7 and 1.6 m in P2, and they were 1.3 m in P3. Regarding the erosion, the values varied between 0.5 and 1.0 m in P1 and P2, and between 0.25 and 2.1 m in P3. Under summer conditions, more sedimentation was observed in the intertidal sector of P1 (0.6–0.7 m) and P2 (0.6–0.75 m), whilst this is linked to winter conditions in P3 (1.0–1.25 m). The highest erosion in the subtidal sector is associated with winter conditions in the three profiles (0.65–2.1 m). The analyses of seasonal changes in the Samil Beach profile allowed us to determine the deepest point at which sediment movement is observed. The maximum movement depth is −7.1 m in P1 and P3 (Figure 6a,c) and is −7.4 m in P2 (Figure 6b).

The profile P4, located in an area of rocky outcrops (Figure 6d), presents erosion over the entire profile, with values of 0.75–1 m under winter conditions and a value of 0.5 m under summer conditions. The maximum sediment movement depth in this profile is located at a −6 m isobath. In all cases, the final slope is similar to the original one (0.033).

The profiles P5 and P6 (Figure 6e,f) are located in a zone between Samil and Vao Beaches, and in Vao Beach, respectively. The seabed of P5 consists of rocks up to a 5 m depth. The profile P6 presents a similar behaviour to P3, showing a cutting point between the original and resulting profiles that defines the accretion zones at shallow depths and the erosion zone in the deepest part of the profile. The accretion is 0.75 m, whilst the erosion is 0.25 m. The slope of the resulting profile in winter conditions is 0.032, higher than the original, which is 0.026. In this case, the maximum sediment movement observed is at a −4.5 m depth.

4.5. Modelling of the Different Route Alternatives

The simulation of the mutual interaction between the structure and natural sediment dynamics in the beach was performed for the two possible pipeline route alternatives. For each alternative, the previously described four wave cases were simulated; NW summer and winter and SW summer and winter. The results described hereafter refer to winter (storm) cases, since in the summer cases the presence of the new pipelines did not cause any significant morphologic alteration.

4.5.1. Southern Route Alternative

Under SW storm waves (not shown), the area of the bay where the pipeline construction is planned is partially sheltered from the wave action by the presence of Toralla Island. Thus, the waves reach the southern corner of the bay with an incoming direction of 300°, practically parallel to the pipeline route. The Hs showed locally higher values (≈0.65 m) over the pipeline, between 5 m and 10 m deep, than in the surrounding area (≈0.55 m). This local increase was also reflected in the orbital velocity, with 0.30 m/s over 0.27 m/s. Along the pipeline route, sediment movement is detected eastward from the shallowest section. The variation in the wave conditions in the shallowest section of the structure produced an erosion of 0.20 m in its east flank with respect to the results without a pipeline presence. This interaction disappears as the depth increases along its route. The results of the model reflected that the structure was located in a low mobility area, since the shallowest sector of the beach, located above the pipeline emerging point, consisted of a rocky seabed. Under these wave conditions, the presence of the structure did not produce any effect in the natural dynamics of the beach (Figure 5a).

Under NW storm wave conditions, the southern corner of the bay is exposed to the incident waves, and Toralla Island does not have a strong influence, as in the previous case. The Hs in the southern corner of the bay is between 0.5 and 1 m, again in a similar incoming direction to the pipeline route (300°). Similarly, with respect to the previous case, in the shallowest section of the pipeline (<8 m), a variation in the Hs over the structure is produced −0.6 m over it and 0.55 m in the surrounding area. This increase was also produced in the bottom orbital velocity −0.25 m/s over 0.22 m/s. Thus, the alterations produced in the waves only cause sediment transport in the shallowest section of the pipeline. The presence of the pipeline generates a flank erosion furrow that is similar to the previous case but of smaller magnitude, 0.11 m, in the east flank, with respect to the results without a pipeline. The analysis of the erosion/sedimentation processes in Samil Beach revealed that the presence of the pipeline did not alter the natural behaviour of the beach. Comparing the erosion/sedimentation results with the simulation with no structure (Figure 5b), it is also observed that the pipeline did not produce a shielding effect in the sediment transport in this case, the erosion/sedimentation zones and values being similar.

Ultimately, the southern route alternative showed erosion on the east flank in the shallowest sector of the pipeline related to the alteration of the waves due to the pipeline’s presence. However, this alteration was local and did not affect the global sediment dynamics of the beach.

4.5.2. Northern Route Alternative

The northern corner of the bay is more exposed to the waves than the southern corner. However, in key terms, the presence of the structure did not affect the waves in the bay. Under SW storm waves, a local increase in the Hs in the shallowest section was also detected (Figure 7a), with values of 0.7 m over 0.6 m, equally reflected in the bottom orbital velocity (Figure 7b) −0.35 m/s over 0.29 m/s. The incoming direction of the waves that reach the beach is 310°, an angle of 10° with respect to the pipeline. The shallowest section of the pipeline (<8 m) is located in an area of general erosion (Figure 7c). The presence of the pipeline in this sector produced an alteration of 0.15 m in its east flank followed by a 0.11 m accretion with respect to the results without the structure (Figure 7d), whilst at deeper depths, no interactions were observed.

The presence of the planned pipeline under SW storm waves did not alter the general transport pattern of the beach (Figure 7c) with respect to the observed pattern without a structure (Figure 5a). The areas where sedimentation/erosion was produced were the same, with a similar range of values, and no shielding effect in the sediment pattern was detected.

The sediment dynamics under NW winter storm waves changes with respect to the previous case (SW). The waves’ incoming direction in this case coincides with the pipeline’s orientation. The values of the wave parameters, Hs, and the bottom orbital velocity are significantly smaller than they are under SW winter storm waves. These parameters increase over the structure in the shallowest section (<8 m), showing a difference of 0.1 m in the Hs and 0.05 m/s in the bottom orbital velocity, with respect to the surroundings (Figure 8a,b, respectively). The shallowest segment is located in an erosion area, similar with respect to the previous modelled case (Figure 8c) but of smaller magnitude (Figure 8d). In this case, the presence of the pipeline produced minor erosion at its flanks (0.07 m) due to the lower wave energy.

The natural beach dynamics in this case were not altered by the presence of the planned pipeline either, since the waves were barely altered. Comparing the resulting erosion/sedimentation with the results without a pipeline, the sediment transport was the same in both cases.

According to the modelled transport results of southern and northern route alternatives under different storm conditions (SW, NW), two areas can be clearly distinguished in the beach—the shallow zone (between 0 and 8 m deep) and the deep zone (more than 8 m deep). The shallow zone is where the main seasonal sediment movement occurs, producing erosion-sedimentation strips parallel to the shoreline, especially in the norther corner. In the deep area, no sediment movement was detected. The presence of the pipeline in the shallow zone produces a local interaction with the sediment transport, causing erosion zones at both pipeline flanks. Nevertheless, the presence of the pipeline does not alter the global sediment transport pattern in the beach.

5. Discussion

When a hard intervention is performed in a highly dynamic sedimentary system, as in a beach, it must fulfil security, service and exploitation requisites. In this work, the security of the pipeline and Samil Beach was analyzed, allowing for the selection of the most feasible route alternative and identifying future interactions that would occur between the planned route alternatives and the sediment beach dynamics.

5.1. Effect of the Pipeline Route Alternatives on the Beach and Alternative Selection

All data shows that Samil Bay is a low energy zone, in which the waves that arrive are very attenuated with respect to the open sea conditions, between 60% and 75%, due to the Ria de Vigo’s geomorphology. In addition, the presence of Toralla Island exerts a local shielding effect in the southern part of Samil Bay. Thus, the most exposed part of the bay to the wave action coincides with depths of more than 20 m. In the same way, the resulting energy due to the wave shallowing process is concentrated in a narrow strip along the coastline at depths less than 8 m. In addition, the beach sediment shows that the beach is a low energy area where the prevailing mean grain size is lower than 300 µm (medium, fine and very fine sands).

The numerical model Delft3D has been shown to be a useful tool for the assessment of submerged pipeline construction alternatives. The results of the numerical model showed that the pipelines of the two alternatives did not alter the global wave pattern within the bay. Regarding the sediment movement, it was determined that the two suggested routes produced a low affection in the subtidal area of Samil Beach. Only in the shallowest (<8 m) sector of the beach did the presence of the pipelines cause small-scale alterations at the bottom. The field data (bathymetry and SSS) supported the model results, determining that the highest interactions are produced locally at depths lower than 12 m, diminishing as depth increases. Accordingly, the presence of rocky outcrops, some of them large, seems to produce a limited effect on the sediment dynamics since the bathymetry does not reveal strong interaction signs. In the assessment of the best-suited pipeline routes, taking into account the structure security criteria, the results showed that, in the two planned route alternatives, the interaction was similar. Therefore, this criterion was not crucial in the route selection. Nevertheless, the different characteristics of the bottom between the northern and southern routes conditioned constructive aspects at each alternative and, consequently, the work cost. Thus, the alternative selection was based mainly on economic criteria, and the northern route was finally selected.

5.2. Pipeline Emerging Depth in the Northern Route

The constructive design of the northern route alternative consisted of a 620 m long, buried section from the water treatment plant to Samil Beach, and of a second, 3040 m long, submerged section that emerged gradually in the subtidal area of the beach, lying on the seabed for the rest of the route.

For this reason, the definition of the pipeline’s emerging depth was a critical constructive point in the project. The selected depth must also fulfil the security criteria of the structure. The beaches are environments of high temporal variability, where small changes in the wave conditions may lead to important changes in the sediment dynamics and in the morphology of the emerged and submerged sectors. Thus, the maximum sediment movement depth is a crucial parameter in the pipeline design. In engineering terms, this depth is known as the closure depth, a theoretical concept that establishes the depth from which we assume there are no significant morphological changes in the bottom. The definition of this depth allows for the establishment of a security depth for the structure placement.

The closure depth was addressed by a multitool approach in order to increase the reliability of the solution. Thus, existent scientific literature data, theoretical calculations by the Birkemeier equation [32], fieldwork data (multibeam echosounder bathymetry and SSS) and numerical modelling results were compared. All of them were calculated to a zero sea level reference of the Vigo harbor.

Previous studies in Samil Beach [33] established that the closure depth was located between 8 and 12 m in low energy conditions, and up to 17 m in high energy conditions.

The theoretical closure depth in the present work was calculated by Equation (1). For the relation between H_s12 and T_p (Equation (2)), a range of values was obtained in the closure depth calculations; a minimum value was associated with the lowest peak period (T_{12, min}), h*_min = 10.36 m, and a maximum value was associated with the highest period (T_{12, max}), h* _max = 12.47 m. A range of values were obtained in the closure depth calculations; a minimum value associated to the lowest peak period (T_{12, min}), h*_min = 10.36 m, and a maximum value associated to the highest period (T_{12, max}), h* _max = 12.47 m. The study of the bathymetry and the SSS images highlighted the existence of erosion structures at both flanks of the old pipelines. These structures were identified up to approximately a 12 m depth. This can be considered to be the maximum depth at which the sediment can be moved due to dynamical processes. Therefore, this can be interpreted as direct data regarding the beach closure depth.

The use of the Delft3D numerical modelling and the analysis of transversal beach profiles and its temporal evolution established a maximum profile variation at around a 7 m depth. From this point on, there were no changes in the morphology. Though the model has not been validated with measurements, comparison with in situ data acquired by other technical means shows that the model slightly underestimates the morphodynamic beach behavior.

In

Table 3. Closure depth determined from different tools and sources.

Data Type	Closure Depth (m)
Scientific literature	8–12
Birkemeier (1985) equation	11.9 (10.36–12.47)
Multibeam echosounder bathymetry	12
Side Scan Sonar	12
Delft3D model	7.4

, depth values determined from different tools and sources are summarized:

The closure depth data determined from the visual inspection of the bathymetry and SSS images coincides with the main terms of the Birkemeier formula [32], even though this was developed for open sea beaches. This is most likely due to the extreme storm events’ effect on the obstacles. When these short but very energetic events take place, these waves interact with the obstacles located at deeper depths than those that usually produce erosion/sedimentation structures. These structures would have a permanent character, since the usual wave conditions could not modify them due to the depth.

The closure depth given by the scientific literature and numerical modelling is smaller than the theoretical and observed one. The scientific literature presents a range of values whose upper limit (12 m) coincides with the closure depth determined by the Birkemeier equation [32]. The numerical modelling establishes the seasonal variability of the beach profile up to a 7.4 m depth, with a horizontal forward-retreat of the profile that can exceed 100 m. This model is not calibrated in absolute terms, but it is important to highlight that the resulting data is within the range of the beach morphodynamics under a similar hydrodynamic regime.

6. Conclusions

Several data collection techniques and numerical modeling approaches have been integrated to explore the interaction between two pipeline construction alternatives and Samil Beach. The submerged pipeline planned for Samil Beach would not affect the global sediment dynamics of the beach, based on the simulated cases. Most of the pipeline is outside the sector where the main transport processes take place, avoiding a shielding effect in the main transport paths. Thus, the seasonal morphodynamic behaviour of the beach would not be affected.

Nevertheless, the structure would be affected by the beach sediment dynamics, producing mainly erosion in its surrounding seabed until the closure depth is reached. This interaction does not hazard the integrity of the pipeline in the case that it emerges deeper than the closure depth. Taking into account the values obtained for the closure depth in the beach (8–12 m), the pipeline emerging point must be located above an 8 m deep isobath in order to avoid its integrity hazard, which represents the seasonal beach sediment movement. The results also indicate that, if the pipeline emerges on a seabed above a 12 m bath isobath, its integrity would be guaranteed.

Thus, the pipeline-beach interaction criterion for the best suited pipeline alternative selection has been discarded. The multitool strategy has proved to be efficient in assessing the sediment dynamics of the Samil Beach environment and the effects of each pipeline route alternative on the beach. It has also allowed us to identify the consequences of the seasonal sediment dynamics on the pipeline alternatives; something essential for the security of the structures.

This study highlights the advantage of exploring the effect of submerged infrastructures with a multidata approach to fully understand the diversity of potential influences on local beach behavior.