A Lagrangian oil spill transport model for the Red Sea

Ocean Engineering 217 (2020) 107953 Available online 28 August 2020 0029-8018/© 2020 Elsevier Ltd. All rights reserved. Contents lists available at ScienceDirect Ocean Engineering journal homepage: www.elsevier.com/locate/oceaneng A Lagrangian oil spill transport model for the Red Sea R. Periáñez Dpt Física Aplicada I, ETSIA, Universidad de Sevilla, Ctra Utrera km 1, 41013 Sevilla, Spain ARTICLE INFO Keywords: Red Sea Tide Baroclinic circulation Oil spill Transport Oil evaporation ABSTRACT An oil spill model was developed for the Red Sea. It uses detailed data on water circulation: instantaneous tidal currents, tidal residuals (both from a tidal model specifically applied to the Red Sea) and baroclinic circulation derived from HYCOM ocean model. The model can also incorporate forecasts of local (in the spill area) winds. The transport model is Lagrangian and includes advection/diffusion plus specific processes for oil, as buoyancy, decomposition and evaporation. Results of the tidal model were compared with observations in the sea. The transport model was first applied to passive particles to compare the significance in advective processes of the components of the current (tidal, tidal residual, baroclinic currents) over several areas and seasons and then it was applied to some oil spills: the model output was compared with results from previous simulations and it was also applied to hypothetical oil spills. 1. Introduction Lagrangian transport models are useful tools which have been extensively used in oil spill modelling in different marine areas Korotenko et al.,2004;Periáñez and Pascual-Granged,2008;Proctor et al.,1994b,a;Jordi et al.,2006 among many others). Efforts were also made after the Deepwater Horizon oil spill in the Gulf of Mexico in 2010 (Mariano et al.,2011) and the Prestige accident in NW Spain in 2002 (Carracedo et al.,2006;Sotillo et al.,2008). The Red Sea is located between the Indian Ocean and the Mediterranean Sea, connected to the last through the Canal of Suez. It is connected with the Gulf of Aden (Indian Ocean) through Bab el- Mandeb (BeM) Strait. Because of the dry climate there is an excess of evaporation over precipitation and there is not river runoff. Dense water flows out the Red Sea through BeM, and this subsurface outflow is compensated by a surface inflow from the Indian Ocean (Pickard and Emery,1982). This exchange pattern is modulated by seasonal (monsoon) variations (Guo et al.,2018). The winter flow regime is apparent between October and May, and it is characterized by a twolayer inverse estuarine flow. This pattern is replaced in summer by a three-layer flow characterized by a shallow outflow of the Red Sea surface waters, an intrusion of relatively fresh and cold Gulf of Aden intermediate waters, and finally a deep hypersaline outflow of the Red Sea waters (Sofianos and Johns,2007). Observations made by these authors reveal a very intense and complicated circulation pattern, with a permanent cyclonic gyre in the northern part and a number of them in the central and southern basins. Tides (Pugh,1987) are closely represented by a standing wave with a single central node in the case of semidiurnal tides. Diurnal tides are similar, but the node is displaced to the south. E-mail address: [email protected]. Due to its strategic location between the Mediterranean Sea and the Indian Ocean, the Red Sea suffers intense shipping activities, including of course the traffic of oil tankers. Thus, it is an area potentially exposed to oil spills due to accidental releases which could damage the environment. Therefore, it is useful to have knowledge on the transport processes which occur within the Red Sea as well as predictive tools which could support decision-making in case of an accidental oil spill. However, few oil spill models for the Red Sea are described in literature. Nasr and Smith (2006) present a local spill model for the Gulf of Suez, but they do not provide details on how water circulation (the essential transport driving mechanism) is obtained: only that a water current map was constructed integrating data from different sources, including observations by sea captains. Ahmed et al. (2012a) describe a local model as well, applied to a hypothetical spill in Sudanese coast. Again, surface water currents are provided to the transport model in a ‘‘map form that identifies the spatial variation in the water velocities’’ (sic), but details are not given. The model was run for prevailing (steady) wind conditions of summer and winter. Ahmed et al. (2012b) present an Eulerian model for the same area, which is used to estimate a diffusion coefficient. Again, it is only mentioned that currents were obtained from field measurements. In this paper, an oil spill transport model which covers the whole Red Sea is presented. The most significant differences with respect to previous works are: 1. Model domain: it covers the whole sea instead of being a local model. https://doi.org/10.1016/j.oceaneng.2020.107953 Received 23 April 2020; Received in revised form 11 August 2020; Accepted 12 August 2020 Ocean Engineering 217 (2020) 107953 2 R. Periáñez 2. Water circulation: the model uses detailed data for water circulation. It includes five tidal constituents and their corresponding residuals derived from a high resolution tidal model, previously tested and validated in other regions. Tidal analysis is used to calculate the tidal current at any time and position in the Red Sea. In addition, baroclinic circulation (daily values) is obtained from the widely tested and applied HYCOM ocean model. 3. Winds: forecasts of local winds (in the spill area) may be included in the model to deal with changing wind conditions both in speed and direction, caused for instance by a local storm not described within HYCOM. Of course, a model is just a description of very complex environmental processes. If the purpose of a given model consists of its use as a tool to support decision-making after an accidental oil spill in the sea, it has to be relatively simple to run (required input data should be as few as possible) and running times should be as short as possible while offering reasonable predictive skills. These conditions imply that approximations and simplifications are required; thus some processes have to be parameterized or neglected, retaining the most significant ones. The different simplifications used in the present model are described along the paper. Nevertheless, it goes further than considering the oil spill as a number of passive particles experiencing advection and diffusion as in Bourgault et al. (2014) and Díaz et al. (2008) for instance. The model is described in Section 2, where the hydrodynamic and spill transport components are presented separately. Results are described in Section 3; results from the tidal model are initially presented and compared with observations, as well as some examples of HYCOM baroclinic circulation (Section 3.1). Results from the transport model are given in Section 3.2. First, the model was applied to passive particles to assess the significance of the components of the water currents (tidal current, tide residual current, baroclinic circulation) in transport processes in different regions of the Red Sea and different seasons. Then some examples of applications to oil spills are shown. The model output was first compared with results of some other earlier simulations and then to hypothetical oil spills to illustrate model capabilities. 2. Methods 2.1. Tidal modelling A two dimensional depth-averaged model was used to simulate tides in the Red Sea. The model provides the sea surface elevations and water currents at each point in the model domain and for each time step. Elevations and currents are treated through standard tidal analysis (Pugh,1987, Chapter 4) and tidal constants are then calculated and stored. Five constituents were considered: three semidiurnal (𝑀2,𝑆2 and 𝑁2) and two diurnal (𝐾1and 𝑂1). Thus, the tidal model equations (see for instance Periáñez,2012) are solved for each constituent and tidal analysis is also carried out for each constituent separately. A residual transport cannot be produced by the pure harmonic currents given by the tidal analysis, thus tidal residuals were calculated as well, in the form of Eulerian residual transport. The procedure may be seen in detail in Periáñez (2012). Hydrodynamic equations (which essentially express conservation of mass and momentum) are solved using explicit finite difference schemes, with a second order accuracy scheme for nonlinear terms. Boundary conditions consist of specifying water surface elevations, from measured tidal constants, along the open boundaries of the domain. Measurements are given, for instance, in Madah et al. (2015). The model has been successfully tested in the past for several regions at different spatial scales (Periáñez,2007,2009,2012;Periáñez et al.,2013;Periáñez and Abril,2014). 2.2. Baroclinic circulation HYCOM (Hybrid Coordinate Ocean Model, Bleck,2001) output was used to describe baroclinic circulation in the Red Sea. It is a primitive equation general circulation model with 40 vertical layers whose thickness increase from surface to the sea bed (2 m near the surface to 1000 m in the deepest ocean). Horizontal resolution is 0.08◦ both in latitude and longitude. A wide number of HYCOM model applications over all the world may be seen in the model web page (https://www.hycom.org/). Daily currents for the area of interest were downloaded from HYCOM data server. 2.3. Particle-tracking oil-spill model A pollutant release into the sea is simulated in a particle-tracking (Lagrangian) model by means of a number of particles, each one equivalent to a number of units, whose trajectories are calculated along time during the simulated period. The spill model considers physical transport (advection due to water currents and mixing due to turbulence) plus some specific processes for oil, as buoyancy, evaporation and decomposition. The transport of the oil spill as a surface film and the effects of surface tension are not included in this model, although their theoretical description may be seen in Maderich et al. (2012). Advection is computed solving the following equation for each particle: 𝛥𝑥 =𝑢 𝛥𝑡 +𝜕𝐾ℎ 𝜕𝑥 𝛥𝑡 (1) 𝛥𝑦 =𝑣 𝛥𝑡 +𝜕𝐾ℎ 𝜕𝑦 𝛥𝑡 (2) where 𝛥𝑥 and 𝛥𝑦 are the changes in particle position (𝑥, 𝑦); and 𝑢 and 𝑣are water velocity components in the west-east and south-north directions, respectively, at the particle position and depth; and for the moment when the calculation is done, since currents are changing in time. Derivatives of the horizontal diffusion coefficient (𝐾ℎ) above prevent the artificial accumulation of particles in regions of low diffusivity (Proehl et al.,2005;Lynch et al.,2015). Constant diffusion coefficients are used in the present work, thus these term are not relevant here. The use of constant diffusivities is just a simplification to speed up calculations, although of course more complex descriptions could be implemented, as the Smagorinsky scheme (Cushman-Roisin and Beckers,2011). Actually, Lagrangian transport models were applied to the Pacific Ocean (Periáñez et al.,2019), some with constant and some with Smagorinsky diffusivities, providing comparable results. It must be commented that a first order accuracy equation is used to describe advection. However, Elliott and Clarke (1998) did not find improvements in results when a second order accuracy scheme was used to simulate the movement of surface drifters, with respect to the first order equation. Moreover, in marine transport problems, the effects of turbulence will mask any small errors in the advection scheme. Vertical currents in the sea are small; dimensional analysis (Dyke, 2001) shows that they are typically in the order of 10−5to 10−4 m∕s. Thus vertical advection is masked by vertical mixing due to turbulence. As a consequence, it is a common approach to neglect vertical advection in marine transport modelling. Water currents, 𝑢(𝑥, 𝑦, 𝑧, 𝑡)and 𝑣(𝑥, 𝑦, 𝑧, 𝑡), are the addition of tidal current, tidal residual, baroclinic current and local wind induced current. Of course, this is an approximation which neglects the non-linear interactions which may occur between baroclinic currents and tides, as well as non-linear interactions between the different tide constituents. The baroclinic current for each day of the year is obtained from the downloaded HYCOM data, as mentioned above. Tidal sea surface elevation for the corresponding instant of time 𝑡at a given location is obtained from the calculated tidal amplitudes and phases as: 𝑍(𝑡) = 𝐻0+ 5 ∑ 𝑖=1 ℎ𝑖𝑓𝑖cos(𝑤𝑖𝑡−𝑔𝑖+𝑉𝑖)(3) Ocean Engineering 217 (2020) 107953 3 R. Periáñez where 𝐻0is the location datum, 𝑤𝑖is frequency of constituent 𝑖(five are included as mentioned before), ℎ𝑖and 𝑔𝑖are amplitude and phase (adapted phase, i.e., for the local time meridian) for the corresponding location (calculated from the tidal analysis), 𝑓𝑖is nodal factor and 𝑉𝑖the equilibrium argument of the constituent at Greenwich. Nodal factors and equilibrium arguments for year 2020 are used. Thus, the beginning of this year is 𝑡= 0. Of course, values for any other year may be used. The same treatment is given to tidal currents, as usual in this type of models (Proctor et al.,1994b;Elliott et al.,2001;Periáñez and Pascual- Granged,2008). A standard current profile is used to generate a vertical structure in tidal currents (depth averaged currents are provided by the tidal model), since these currents decrease from sea surface to the bottom due to friction with the seabed. Details may be seen in Pugh (1987) and Periáñez and Pascual-Granged (2008) and are not repeated here. The wind-induced current decreases, below a depth 𝑧1(wind-driven surface layer thickness), logarithmically to zero at a depth 𝑧2. The mathematical form of this profile may be seen in Pugh (1987), for instance. As in Proctor et al. (1994b), it was used 𝑧1= 0.1m and a typical value of 20 m was selected for 𝑧2. The surface current is a percentage of the wind speed, usually in the range 2%–3%. In order to account for Stokes drift, the wind-driven layer speed was increased to 3.5% of the wind speed (Proctor et al.,1994b). Of course, this is just a parameterization of a very complex non-linear process which is often applied in this type of models. Although the surface current should be deflected to the right (northern hemisphere), this angle is small enough to be safely neglected (Proctor et al.,1994b). The maximum size of the horizontal step given by the particle due to turbulence, 𝐷ℎ, is (Proctor et al.,1994b;Hunter,1987;Periáñez and Elliott,2002): 𝐷ℎ=√12𝐾ℎ𝛥𝑡 (4) in the direction 𝜃= 2𝜋𝑅𝐴𝑁, where 𝑅𝐴𝑁 is an uniform random number between 0 and 1. This equation gives the maximum size of the step. In practice, it is multiplied by another independent random number to obtain the real size at a given time and for a given particle. This is required to ensure a simulation of a Fickian diffusion process (Proctor et al.,1994b;Hunter,1987). Similarly, the size of the vertical step is (Proctor et al.,1994b; Hunter,1987;Periáñez and Elliott,2002): 𝐷𝑣=√2𝐾𝑣𝛥𝑡 (5) given either towards the sea surface or the sea bottom. A new independent uniform random number 𝑅𝐴𝑁 between 0 and 1 is generated and used to decide the jump direction through the following equation: 𝛥𝑧 = 2 ⋅(𝑖𝑛𝑡(2 ⋅𝑅𝐴𝑁)−0.5)𝐷𝑣(6) In these equations 𝐾ℎand 𝐾𝑣are the horizontal and vertical diffusion coefficients respectively and 𝛥𝑡 is time step. Constant typical values are used for the diffusion coefficients: 𝐾𝑣= 1.0 × 10−5 m2/s and 𝐾ℎ= 10 m2/s (Elliott et al.,2001;Mariano et al.,2011;Periáñez et al., 2019). Time step was fixed as 120 s. There is not a stability condition in Lagrangian models equivalent to the Courant–Friedrich–Levy condition (Kowalik and Murty,1993) of Eulerian models, but it is required to ensure that each particle does not move through a distance that exceeds the grid spacing during each time step. This was satisfied with the selected time step. Linear interpolation in time and space is used to determine the currents at each particle position at the desired time. Note that time interpolation is required only for baroclinic currents, since tidal currents at a given 𝑡are determined from an analytical equation in the form of Eq. (3). Specific processes for oil are described in what follows. Oil droplets have a size distribution so that larger ones tend to remain in the water surface and move in the direction of wind. Smaller ones, in contrast, tend to mix downwards. The diameter of each oil droplet in the simulation is assigned randomly between a minimum and maximum Fig. 1. Evaporation curves for Arabian Heavy and Light oils and regular Diesel at 25 ◦C temperature (from equations in Fingas,2015). The dashed line is the evaporation curve for Diesel as given by the stochastic method described in Section 2.3. diameter, 𝑑𝑚𝑖𝑛 −𝑑𝑚𝑎𝑥. Then, the procedure described in Korotenko et al. (2004) and Periáñez and Pascual-Granged (2008) is used to estimate the vertical velocity of oil droplets due to buoyancy forces. Of course, oil droplets have a size distribution. Selecting their sizes between 𝑑𝑚𝑖𝑛 and 𝑑𝑚𝑎𝑥 is just a simplification which minimizes the required input data to run the model. These numbers should represent the quality of the modelled oil, but if they are unknown quoted values may be used. For instance, Proctor et al. (1994b) suggest 𝑑𝑚𝑖𝑛 = 60 μm and 𝑑𝑚𝑎𝑥 = 120 μm. The model includes the effect of oil evaporation from the surface and decomposition within the water column (because of biodegradation). These processes are treated in a stochastic form based on decomposition and evaporation e-folding times (details may be seen for instance in Proctor et al.,1994b). Only particles which are in a surface layer of 25 cm thickness may evaporate (Proctor et al.,1994b), while oil at any depth can experience decomposition. The age of each released particle (time since it was released) is calculated. Particles cannot decay after an specified age. This simulates the rapid evaporation which occurs after the release and ensures that some particles remain in the simulation as old and weathered oil (Proctor et al.,1994b). The e-folding time for evaporation can be deduced from the generic equations proposed for a wide number of oil types by Fingas (2015). These equations provide the fraction of evaporated oil (%𝐸𝑉 ) as a function of time and temperature. As examples (Fingas,2015, table 7.2), for Arabian Heavy oil the equation is : %𝐸𝑉 = (2.71 + 0.045 ⋅𝑇) ln(𝑡)(7) and for regular Diesel it is: %𝐸𝑉 = (0.31 + 0.018 ⋅𝑇)√𝑡, (8) where 𝑡is time in minutes and 𝑇is water temperature in ◦C. The evaporation e-folding times may be estimated from this type of equations. Evaporation curves are represented in Fig. 1 for a temperature of 25 ◦C as an example for Diesel and Arabian Light and Heavy. In the case of Diesel, evaporation e-folding time is about 4 days, but it is much longer in the case of Arabians. Decomposition times are usually longer than evaporation ones (Proctor et al.,1994b). Once an e-folding time 𝑇𝑒is known, a probability is defined as (Proctor et al.,1994b): 𝑝= 1 − 𝑒−𝛥𝑡∕𝑇𝑒,(9) and an independent uniform random number 𝑅𝐴𝑁 is generated. If 𝑅𝐴𝑁 ≤𝑝then the particle is removed from the computation. This procedure gives a realistic representation of the evaporation process. As Ocean Engineering 217 (2020) 107953 4 R. Periáñez an example, the evaporation of Diesel (e-folding time 4 days) calculated with this method is represented in Fig. 1 as the dashed red line, which can be compared with the solid red line representing Eq. (8). If a particle reaches the shoreline during a simulation it is considered as ‘‘beached’’ and it does not move anymore. The shore characteristics, as material composition or its slope, are not included in the model as a simplification. Otherwise, a huge amount of information (not always easily available) should be provided as input data. Nevertheless, these simplifications may lead to a overestimation of beaching. The model provides the final positions of particles and the mass of oil per unit of surface over the domain. 2.4. Computational scheme The tidal model was run in advance to calculate and store tidal constituent amplitudes, phases and residuals in files which are read by the transport code during the spill simulations. The tidal model domain extends from 10◦to 30.5◦N in latitude and from 32◦to 45◦E in latitude. The same area is used to download HYCOM circulation and it is also used in the oil-spill model. Bathymetry is obtained from NOAA ETOPO1 data, using its on-line tool to design a grid1at a resolution of one minute of arc in both longitude and latitude. The tidal model was started from rest and run over 20 tidal cycles; which is enough to obtain a steadily oscillating solution. It was also run separately for each constituent. Tide elevations were specified along the only open boundary of the domain (south-east) as mentioned before. The component of the current which is parallel to the open boundary was set to zero as usual in tidal computations. As mentioned before, standard tidal analysis was used to calculate the tidal constants for each point in the domain and constituent. A number of input files are read by the particle-tracking transport code. These are the following: 1. Baroclinic currents previously downloaded from HYCOM website. 2. Tidal amplitudes and phases for the five included constituents, as well as Eulerian residual transports for each of them. 3. Tidal currents and residuals may be individually switched on and off (for instance to evaluate the significance of a given constituent or to make comparisons between them). These switches are provided in a specific file as zeros or ones for each constituent tidal current and residual. 4. Equilibrium arguments and nodal factors for each constituent and for the year of simulation. 5. Bathymetry data with 1 min of arc resolution from ETOPO1 dataset. This seems a reasonable compromise between resolution and computational requirements, although bathymetric data with higher resolution exist. 6. Release data: date, time, position, depth of release and release magnitude and duration. The duration of the simulation in days is also specified. 7. Oil properties: density, droplet size range and e-folding times. 8. Local (in the area of the accident) wind data if required. This information may be obtained from local observations or weather forecasts. Of course the wind is not uniform in speed and direction in such a large domain as the Red Sea. However, the local wind in the release area is considered uniform. Local wind data are provided in a file as a number of different ‘‘wind episodes’’ (any number can be used), each one characterized by a wind speed, direction and start and end of the episode in hours after the beginning of the spill. This way timeevolving wind conditions may be used (which can be obtained from weather forecasts). It must be noted that winds are already included 1https://maps.ngdc.noaa.gov/viewers/wcs-client/ Fig. 2. General scheme of the model working procedure. in the HYCOM calculations of circulation. However, with the present definition of ‘‘wind episodes’’, we have the opportunity of describing oil transport (in a short temporal scale, thus the spill does not travel very far and the wind may be considered uniform) in case that the accident occurs during a local storm, for instance, which is not described in HYCOM. But specifying this local wind is not a requirement to run the model. Once all files are read, the transport simulation starts. The concentration of oil is calculated from the number of particles within each grid cell (of ETOPO1 bathymetry) as usual (Periáñez and Elliott,2002). A general scheme of the model working procedure, summarizing required inputs, is presented in Fig. 2. 3. Results Results of the tidal model will be described before, and compared with observations and previous models. Next, advection processes are analysed and some examples of applications of the spill model are presented. 3.1. Hydrodynamics Tidal charts (corange and cotidal maps) were obtained for the five considered constituents. We only present an example for a semidiurnal tide (𝑀2,Fig. 3) and another (Fig. 4) for the diurnal constituents 𝐾1. These maps are in close agreement with the previous calculations in Madah et al. (2015), which were carried out with the commercial software Delft3D; Gharbi et al. (2018), with COHERENS model, and Guo et al. (2018): an amphidromic system is apparent in the central Red Sea (around 20◦latitude) in the case of the semidiurnal constituents. This amphidrome is displaced to the south (around 16◦latitude) in the diurnal tides. In general our model reproduces the behaviour of tides in the region. A comparison between the observed and calculated tide amplitudes for the 5 considered constituents may be seen in Fig. 5. The calculated amplitude of the tidal current is shown in Fig. 6 for the 𝑀2constituent as an example. The largest amplitudes are found in the area of BeM Strait; but also in areas near the east shore, where currents increase their magnitudes in the straits between the many small islands. Significant currents are observed in the Gulf of Suez as well. The small tidal currents calculated in the central Red Sea are in agreement with the observations in Churchill et al. (2014), which indicate values below 0.04 m/s. Two examples of baroclinic circulation, as downloaded from HYCOM model website, may be seen in Fig. 7, corresponding to the end of January and July. Thus, they are representing winter and summer conditions. These conditions are clearly different: the flow is mainly Ocean Engineering 217 (2020) 107953 5 R. Periáñez Fig. 3. Calculated chart for the 𝑀2tide. Fig. 4. Calculated chart for the 𝐾1tide. Fig. 5. Comparisons between calculated and observed tide amplitudes for several locations in the Red Sea and for semidiurnal (left) and diurnal (right) tides. Ocean Engineering 217 (2020) 107953 6 R. Periáñez Fig. 6. Calculated amplitude of the 𝑀2tidal current (m/s). directed to the north along the east side of the sea in winter. This circulation essentially disappears in summer, when the flow reverses and is directed to the south in BeM Strait. This effect, due to monsoon reversal, is known since a long time (Patzert,1974). 3.2. Transport processes Some numerical experiments were initially carried out to evaluate the relative significance of tidal currents, tidal residuals and baroclinic circulation in the transport of contaminants in the Red Sea. An instantaneous release of a passive tracer (only suffers advection and diffusion) was simulated in four different areas of the sea. Six simulations were carried out for each site: one included the three sources of currents (tidal currents, tidal residuals and baroclinic circulation), in another one tidal currents were excluded and the third one excluded both tidal currents and tidal residuals. Finally, these simulations were carried out for a typical situation in winter (February 1st, 2020) and another in summer (August 1st, 2020). The positions of particles 30 days after the release were obtained. Results of the six simulations for a release occurring in the region near BeM Strait are presented in Fig. 8. First, it may be seen that in summer particles travel a short distance towards the south, but in winter there is significant transport to the north. There are essentially not differences between the three hydrodynamic configurations in summer, thus baroclinic circulation is dominant. If tides are removed in the winter simulation, then northwards transport is reduced. This indicates that tidal mixing is relevant in this area during winter. Results of simulations for a release in the central Red Sea may be seen in Fig. 9. Again, particles travel a short distance towards the south in summer, but in winter there is significant transport to the north. The different hydrodynamic configurations do not lead to noticeable differences in results. Thus, baroclinic circulation is clearly dominant. Similar conclusions can be obtained in the case of a release occurring near the east shore (Fig. 10). However, in this case, differences between summer and winter conditions are not as significant as in previous cases. Only a higher transport towards the central sea is observed during summer. In contrast, transport is mainly directed to the north in winter. The last experiment results, for a release in the northern Red Sea, are presented in Fig. 11. Again there are significant differences between winter and summer conditions. There is a weak circulation in winter, which leads to a short travel distance of particles towards the south. However, in summer there is a significant transport towards the south, but then particles enter in a gyre which direct them to the north again. Some particles reach the northeast shore of the Red Sea. As general conclusions, it may be said that tidal residuals constitute a weak forcing agent in transport processes within the Red Sea. Tides are significant only in the region near BeM Strait, where actually tidal Fig. 7. Baroclinic currents as calculated by HYCOM model for the end of January and July. Only one of each 16 vectors are drawn for more clarity. Ocean Engineering 217 (2020) 107953 7 R. Periáñez Fig. 8. Position of particles 30 days after an instantaneous release in the BeM Strait area in summer (top) and winter (bottom) for several hydrodynamic configurations. The black dot is the release point. currents are stronger (Fig. 6) and thus enhance initial mixing of the spill. Clearly, baroclinic circulation is the dominant forcing of transport within most of the sea, both for typical winter and summer conditions. The oil transport model was compared with the results of the simulation presented in Nasr and Smith (2006). It consisted of a hypothetical surface release of Arabian Light crude occurring in the Gulf of Suez in a point located in the shipping path line. The accident occurred on July 29th and only one day was simulated in Nasr and Smith (2006), consequently the same time was used. Wind data was provided in Nasr and Smith (2006) as well, obtained from the Egyptian Meteorological Authority. Wind was changing along the simulation time, both in speed and direction, thus this simulation additionally allowed us to test our method to specify local wind conditions. Light crude oil density was fixed as 830 kg/m3(Hollebone,2015). Evaporation and decomposition e-folding times are long in comparison with the simulated time and thus these processes were neglected in the present simulation. Results of our simulation may be seen in Fig. 12 (top), which are in good agreement with those of Nasr and Smith (2006) — their figure 3. The spill travels about 20 km in southeast direction and the Egyptian shore is contaminated along several km (red dots in Fig. 12, which indicate beached particles). The same simulation was repeated excluding the local wind and considering only baroclinic currents. These results are presented in Fig. 12 (bottom), which show a totally different transport pattern. Thus, although atmospheric forcing is included within HYCOM ocean model, the specification of local wind conditions may be required in some cases, specially in coastal locations where the resolution of the ocean model used is not high enough. Ahmed et al. (2012a) simulated the transport of oil in Sudanese coastal waters after some hypothetical releases in winter and summer. They used mean prevailing wind directions in each season (thus constant wind was assumed) and did not provide details on the release date and how water currents were specified; only that ‘‘they are provided in map forms that identifies the spatial variations in the water velocities’’ (sic). We have compared our model with results for their winter experiment (previous simulation in Fig. 12 was done for summer). We also assumed that the release occurred on December 30th (since a date was not provided as mentioned before). Simulation time was half a day, thus evaporation and decomposition were neglected, and Sudanese oil with density 854 kg/m3was spilled. Results of our simulation may be seen in Fig. 13, which shows that about 10 km of shore was contaminated with oil. Our results generally are in agreement with Ahmed et al. (2012a). An oil spill occurred on October 11, 2019, in the point with coordinates 38.33◦E, 21.13◦N after the explosion of the oil tanker Sabiti. The resulting spill was simulated, and results were compared with a satellite picture taken two days after the incident (see Fig. 14). The comparison is quite limited because the spill most likely occurred as a continuous release, but there is not information about timing and release rates and thus we simulated a single event (instantaneous release). Nevertheless, the model indicates that there is some transport to the north of the spill point, although most of the oil travels to the south. This branch is curved to the east due to the presence of an eddy. The satellite image shows that the transport is directed more in northwest and southeast directions rather than north and south, but the overall picture provided by the model seems realistic taking into account the limitations of Ocean Engineering 217 (2020) 107953 8 R. Periáñez Fig. 9. Position of particles 30 days after an instantaneous release in the central Red Sea area in summer (top) and winter (bottom) for several hydrodynamic configurations. The black dot is the release point. Fig. 10. Position of particles 30 days after an instantaneous release in the east coastal region of the central Red Sea in summer (top) and winter (bottom) for several hydrodynamic configurations. The black dot is the release point. Ocean Engineering 217 (2020) 107953 9 R. Periáñez Fig. 11. Position of particles 30 days after an instantaneous release in the northern area of the Red Sea in summer (top) and winter (bottom) for several hydrodynamic configurations. The black dot is the release point. the simulation. Only baroclinic circulation was used and some generic parameters for oil (density 940 kg/m3, particle size range 100−300 μm, evaporation and decomposition neglected due to the short simulation time) were used. A hypothetical oil spill was additionally simulated to show the functioning of the model at longer temporal scales than simulations shown before. A hypothetical heavy oil was used for the simulations, with evaporation and decomposition e-folding times fixed as 150 days and 400 days respectively, since the last process typically has a scale of years (Mariano et al.,2011). Density of heavy oil is in the range of 880 to 1000 kg/m3(Hollebone,2015), and an intermediate value of 940 kg/m3was used here. Oil droplet size is in the range of 100 to 300 μm. This selection allows studying particles which will mainly stay close to the surface (Mariano et al.,2011) and consequently will be more efficiently transported by currents. The spill occurred during spring in this case (April 15th) in the central Red Sea, at the surface, and tides were not included since their effects in this area are not significant as shown before. However, tidal residuals were included. Results of a 20 day long simulation can be seen in Fig. 15. The spill is directed to the north by baroclinic currents and during this transit some particles are beached in the eastern shore, approximately at 20◦N latitude. An eddy turns particles back to the south and some oil reaches the islands located close to the western shore at about 19◦N. A balance indicates that 14.5% of oil was beached, 70.0% remains in water and 15.5% has decayed or evaporated after 20 days. If evaporation and decomposition were neglected, then the fraction of oil remaining in water is 84.6%; and 15.4% is beached. Given the used droplet size range (see above), all the oil remains in the surface, within a 5 cm thick layer. This thickness increases to 44 cm if droplets are in the range 60–120 μm. Even smaller droplets (30–60 μm) lead to a penetration of 1.9 m. The same simulation was repeated but with the release occurring in fall (October 15th). Results, which may be seen in Fig. 16, are significantly different to those of Fig. 15. Again the oil is initially directed to the north, but the change of direction to the south occurs at a shorter distance from the release point than before, indicating that the eddy is displaced to the south in this season. Also, the coast is not contaminated. As another example, a map showing the final oil amount per unit surface (kg/m2) in this simulation, resulting from the particle distribution shown in Fig. 16, is presented in Fig. 17. These concentrations were calculated over a surface water layer operationally defined as 1 m thick. It was assumed that 106kg of heavy oil were released (1064 m3, since oil density is 940 kg/m3as mentioned before). 4. Conclusions A Lagrangian model which simulates the transport of oil spills in the Red Sea was developed. Water circulation was obtained from two sources: tides and tidal residuals were computed with a two dimensional barotropic model and baroclinic circulation was obtained from HYCOM ocean model, downloaded from its web site. The transport model includes advection by currents and local, episodic, winds, three-dimensional diffusion and specific processes for oil, as buoyancy, evaporation, decomposition and beaching. The transport model was initially applied to passive particles (only suffering advection and diffusion) to evaluate the relative significance