Oil Spill Dispersion Forecasting Models

Abstract

Oil spill models are used worldwide to simulate the evolution of an oil slick that occurs after an accidental ship collision or during oil extraction or other oil tanker activities. The simulation of the transport and fate of an oil slick in the sea, by evaluating the physicochemical processes that take place between oil phase and the water column, is the base for the recognition and assessment of its environmental effects. Numerous oil spill dispersion models exist in the bibliography. The contribution of this chapter is the introduction of a 3D oil slick simulation model developed by the Aristotle University of Thessaloniki, which has been recurrently used in different updated forms and applied in operational mode, since 1991 when it was originally created. The model has been tested in various hypothetical scenarios in North Aegean Sea, Greece, and responded with great success. Findings of the present study highlight the existing experience on the subject and denote the applicability of such models in either tracing the source of a spill or in predicting its path and spread, thus proving their value in real-time crisis management.

1. Introduction

Pollution that comes from a single source, like an oil or chemical spill, is known as point source. Commonly, this type of pollution has bigger impact on the marine environment, but fortunately, it occurs less often than land runoff [1]. Oil spills occur after a collision and/or the sinking of oil tankers, under bad weather conditions, and are extremely hazardous in ports with dense maritime traffic. The environmental impact of oil accidents is immense on both the water ecosystems and the coastal environment, including the urban and economic growth of the affected coastal zones. Strong winds and significant water currents contribute to the creation of accidents and concurrently to the spreading of the spills. Ship collision or malfunction, or simply cleaning and sailing, can contribute to marine pollution by spreading garbage, black or gray water, sludge, water ballast, coatings or even air emissions, in deep sea or at seashore. The most common nautical accidents occur due to sinking or foundering, grounding, structural failure, scuttling, by contact or collision, explosion or fire, or after disappearance or abandonment [2]. The discharge from oil pipelines, oil platforms, and vessels is also causing significant damage to the marine environment and the coastal areas, including the urban and economic growth of the affected coastal zones. The Exxon Valdez (1989), Atlantic Empress (1979), Amoco Cadiz (1978), and Torrey Canyon (1967) (Figure 1) are among the most renowned oil spills in history, that alarmed the scientific community and directed all toward the management of anthropogenic environmental disasters in the marine environment.

Figure 1.

(a) Exxon Valdez (1989); (b) Atlantic Empress (1979); (c) Amoco Cadiz (1978); (d) Torrey Canyon (1967).

Spill modeling has a long history [3]. Since 1960s, numerous oil spill models have been developed to simulate weathering processes and forecast the fate of oil spilled, in terms of providing valuable support to both contingency planners and pollution response teams. These models, developed by various organizations, companies, and researchers, vary in complexity, applicability to location, and ease of use [4]. There are two categories according to the Industry Technical Advisory Committee (ITAC) for oil spill response. One category, known as oil weathering models, estimates how oil properties change over time, but does not predict potential migration of the slick. The second category includes trajectory or deterministic models, stochastic or probability models, hind cast and three-dimensional models. In addition to predicting weathering profiles, these models estimate the evolution of a slick over time. The first-generation models were limited to tracing the movement of rigid bodies under the influence of current and wind advection. Oil spreading was justified by the balance of forces, and although chemical processes played significant role, most of them were ignored in this modeling [5]. Second-generation models have been summarized in [6, 7], where some of the principal forces and major phenomena, determining the behavior and the fate of the chemicals at sea, were analyzed. The recent generation models are reviewed in many papers [8, 9, 10, 11, 12]. Some of the oil spill models that are currently available are: General NOAA Operational Modeling Environment (GNOME), MEDSLIK-II, SeaTrackWeb (STW), Model Oceanique de Transport d’ Hydrocarbures (MOTHY), DieCAST-SSBOM (Shirshov-Stony Brook Oil spill transport Model), COastal Zone OIL spill model (COZOIL), etc.

These operational computer models, pondering major physical processes, efficiently simulate the movement, the spreading and the evolution of the floating chemicals. Strong winds and significant water currents, which are contemplated, are major contributors to both the creation of accidents and faster spreading of the oil. Such models can be used either in hind cast mode capable of tracing the source of a spill, or in forecast mode predicting the path, the horizontal dispersion, and the mass balance, assisting in this way the real-time crisis management.

2. Oil weathering

From the time the oil spreads over the water, under the influence of hydrometeorological conditions (wave, winds, currents, solar radiation, etc.), oil properties (density, viscosity, pour point, etc.), and discharge characteristics (instantaneous, continuous, surface, depth), several processes take place, directly or indirectly related, which disperse the oil and change its properties [13, 14]. These time-dependent processes comprise physical, chemical, and biological ones, such as evaporation, dissolution, emulsification, dispersion, etc. (Figure 2). The combination of these processes is called oil weathering. In other words, weathering is a series of chemical and physical changes that cause spilled oil to break down and become heavier than water. Winds, waves, and currents may result in natural dispersion, breaking a slick into droplets which are then distributed throughout the water. These droplets may also result in the creation of a secondary slick or thin film on the surface of the water.

Figure 2.

Weathering processes that oil undergoes in the sea [15].

In general, weathering processes can be divided into three categories. Processes such as spreading, evaporation, dissolution, dispersion, and emulsification are rapidly occurring (within hours) and have immediate effects on the oil slick, as shown in Figure 3. Sedimentation, biodegradation and photo-oxidation operate more slowly (within months) and comprise the long-term mechanisms for the breakdown of hydrocarbons in the environment. Sedimentation, stranding, and oil-ice interaction are important processes, under distinct environmental conditions. Evaporation is probably one of the most significant processes that affect the surface oil particles, in sea or on coast, and therefore, it will play an important role in the construction of the models presented here.

Figure 3.

Duration of oil slick processes in the sea.

In the following chart, the vertical axis portrays the physicochemical processes that oil is going through, and the time each one of them requires is depicted in the horizontal axis, whereas the line thickness denotes the most critical phase of each stage.

Natural processes may act to reduce the severity of an oil spill or accelerate the decomposition of spilled oil. More analytically, these natural processes are described below [15]:

• Evaporation occurs when the lighter substances within the oil mixture become vapors and leave the surface of the water. This process leaves behind the heavier components of the oil, which may undergo further weathering or may sink to the ocean floor. For example, spills of lighter refined petroleum-based products such as kerosene and gasoline contain a high proportion of flammable components known as light ends. These may evaporate completely within a few hours, thereby reducing the toxic effects to the environment. Heavier oils leave a thicker, more viscous residue, which may have serious physical and chemical impacts on the environment. Wind, waves, and currents increase both evaporation and natural dispersion.

• Emulsification is a process that forms emulsions consisting of a mixture of small droplets of oil and water. Emulsions are formed by wave action, and greatly hamper weathering and cleanup processes. Two types of emulsions exist: water-in-oil and oil-in-water. Water-in-oil emulsions are frequently called “chocolate mousse,” and they are formed when strong currents or wave action causes water to become trapped inside viscous oil. Chocolate mousse emulsions may linger in the environment for months or even years. Oil and water emulsions cause oil to sink and disappear from the surface, which give the false impression that it is gone and the threat to the environment has ended.

• Biodegradation occurs when microorganisms such as bacteria feed on oil. A wide range of microorganisms is required for a significant reduction of the oil. To sustain biodegradation, nutrients such as nitrogen and phosphorus are sometimes added to the water to encourage the microorganisms to grow and reproduce. Biodegradation tends to work best in warm water environments.

• Oxidation occurs when oil contacts the water and oxygen combines with the oil to produce water-soluble compounds. This process affects oil slicks mostly around their edges. Thick slicks may only partially oxidize, forming tar balls. These dense, sticky, black spheres may linger in the environment, and can collect in the sediments of slow-moving streams or lakes or wash up on shorelines long after a spill.

3. Oil spill simulation models

A detailed presentation of the oil slick models developed by the Aristotle University of Thessaloniki, which successfully incorporates most of the aforementioned processes, is following. A concise description of the models’ construction and their operational applications are briefly elaborated.

3.1 SOSM (Sea and Oil Slick) 3D Simulation Model

A model suite, comprising of a coastal hydrodynamic circulation, a wave generation model, and a Lagrangian (tracer) model for the transport and physicochemical evolution of an oil slick, is synthesized by the Aristotle University of Thessaloniki [16]. It follows a model that was developed in the context of an EU DG XI project, in cooperation with the UGMM Group of the Belgian Ministry of Public Health and Environment [17]. Following that, the Oil Spill Model (OSM), based on the PARCEL model, was coupled with the hydrodynamic Princeton Ocean Model (POM), and applied on the Navarino bay, in SW Greece, adjacent to the Ionian Sea, to monitor the pollution impact on the benthic structure and the related fisheries [18]. Parallel to that, the Poseidon pollutants transport model (PPTM) was developed to provide real-time data and forecasts for marine environmental conditions and protect the Greek seas from oil spills, as an operational management tool, which consists of a 3D floating pollutant prediction model coupled with a weather, a hydrodynamic, and a wave model [19]. PPTM was applied in four areas of strategic interest in the Greek seas.

The two basic components of the model, that is, the hydrodynamic and the oil transport, are briefly presented here, followed by the applications under the prevailing circulation conditions in a coastal basin of Greece in the NW Aegean Sea, where the second largest port of the country is located in the town of Thessaloniki, with dense maritime traffic [16].

3.1.1 Model construction

The model is composed of two parts: the hydrodynamic part, comprising of the circulation and the wind wave subparts, and the oil slick transport part. The first is aimed at the description of the 3D nearly horizontal flow velocity vector field due to wind and/or tide. It is written with respect to the free surface ζ(x,y,t) and the depth-varying velocity components u(x,y,z,t), v(x,y,z,t)

$$\frac{\partial \zeta }{\partial t}+\frac{\partial }{\partial x}\underset{h}{\int }\mathit{udz}+\frac{\partial }{\partial y}\underset{h}{\int }\mathit{vdz}=0$$E1

$$\frac{\mathit{Du}}{\mathit{dt}}=-g\frac{\partial \zeta }{\partial x}+{{\rm N}}_h\left(\frac{{\partial }^{2}u}{{\partial x}^2}+\frac{{\partial }^{2}u}{{\partial y}^2}\right)+{{\rm N}}_v\frac{{\partial }^{2}u}{{\partial z}^2}=0$$E2

$$\frac{\mathit{Dv}}{\mathit{dt}}=-g\frac{\partial \zeta }{\partial y}+{{\rm N}}_h\left(\frac{{\partial }^{2}v}{{\partial x}^2}+\frac{{\partial }^{2}v}{{\partial y}^2}\right)+{{\rm N}}_v\frac{{\partial }^{2}v}{{\partial z}^2}=0$$E3

The horizontal eddy viscosity-diffusivity N_h is deduced from the 4/3 law with a characteristic length related to the used Dx, and the vertical eddy viscosity-diffusivity N_v is controlled by the surface friction velocity, the mixing length, the wave height and period. The model is solved numerically by finite differences, on a staggered horizontal Arakawa C grid (Figure 4), and a vertical discretization, following the σ-coordinates transform, using the technique of fractional steps (i.e., resolving explicitly for the horizontal advection and implicitly for the vertical momentum diffusion) [20].

Figure 4.

Arakawa C staggered grid configuration.

The wind waves are resolved by an empirical model WACCAS [21, 22] based on the JONSWAP spectrum. The wave height, wave period, and wave direction on the center of the mass of the slick is used for the computation of the Stokes’ drift and the vertical viscosity-diffusivity. The oil slick simulation and tracking are done in Lagrangian coordinates [23, 24], by replacing the oil mass by a big number (10⁴) of particles or parcels, each one of which represents a group of oil droplets of similar size and composition. The oil spill simulation by a large number of passive (not interacting with the existing hydrodynamics) particles, transported in a Lagrangian frame of reference, is the most effective procedure, used in most contemporary models as: (a) it can describe oil spill geometries and shapes at subgrid level, (b) the model results are not distorted by numerical errors, like numerical diffusion, and (c) the fraction of the initially introduced particles (at the spill source) reaching the coast can be considered as the probability of pollution from a minor spill [25]. According to the tracer technique, it is theoretically possible to track the motion of all the “particles” (parcels) of an oil plume as they are transported by the sea. The density, the buoyancy velocity, the volatility, and other properties are statistically distributed to the particles. Their trajectories are tracked in space and time. Each one is advected according to the local (interpolated) water velocity and it is diffused by means of random walks related to the local diffusivity induced by currents and waves.

Evaporation and emulsification are processes quantified by means of commonly used functions relating the processes to the oil properties (to the fraction density, the temperature, the vapor pressure, the wind and water velocity, etc.). The Mackay approach is used [26, 27] for the computation of the evaporated volume fraction from each oil component. Emulsification is activated for oil fraction with densities and under wave height/length values beyond critical ones. The evaporation and emulsification processes result in: (a) loss of contained oil volume in the evaporated particles and (b) creation of a mixture of water in oil and formation of a floating mousse. Another important oil weathering process accounted in this model is beaching, which is better described by the duration of trapping on the coastal boundary of a beached particle as a function of the coastal morphology (from rocky to flat sandy beach). This process gives crucial information on the oil quantities retained on the coast and the consequent environmental effects.

The spatial and temporal distributions of values like oil concentration (or volume per unit area) are deduced from the number of particles and the contained volume inside each mesh of the discretization grid used for the hydrodynamic analysis. It is worth noting that the time step Dt, used for the integration of the hydrodynamic model, differs by orders of magnitude from the one used for the tracking of the oil slick.

3.1.2 Operational application

The produced operational code was updated and adapted to the environment of a coastal basin of NW Aegean Sea, namely Thermaikos bay. The commercial port of Thessaloniki hosts among others a busy oil terminal. Although the terminal operates under the strict International Maritime Organization (IMO) and European Union (EU) regulations for oil terminal operational safety, the port area is at risk for a potential oil slick accident.

The transport and fate of the slick are closely related to the hydrodynamic conditions in the bay, the bay morphology, and the oil properties. The hydrodynamic circulation in the bay is experiencing mainly N-NW winds, whereas the circulation patterns are regulated by the bay bathymetry and the coastal topography. The bay configuration and bathymetry are shown in Figure 5, along with the initial oil spill location.

Figure 5.

Thermaikos bay configuration and bathymetry, and initial oil slick conditions.

The model was tested for hypothetical scenarios of oil spill under NW and NE winds of 45 knots, to assess the risk of pollution in the coastline of the touristic areas of Chalkidiki (east) and Katerini (west), respectively. The oil is presumably spilled in the port of Thessaloniki (Figure 5), as a result of a hypothetical accidental spill during a tanker operation. An oil spill of 1000 parcels, representing a total volume of 10,000 tn, is released instantaneously on the sea surface. Total simulation run was for 72 hrs (3 days), providing output results every 3 hrs. The evolution of the slick trajectories, for the first 24 hrs and for the two wind-generated circulation scenarios, is depicted in Figures 6 and 7, respectively, documenting the difference in the coastal impact of the slick for the two cases.

Figure 6.

Evolution of the oil slick under strong NW wind.

Figure 7.

Evolution of the oil slick under strong NE wind.

In the case of strong NW wind, the touristic coast opposite to the bay is damaged within some hours after the spillage, while in the case of the NE wind, the oil slick travels along the bay toward the open sea giving time to the relevant authorities for intervention, blocking, and final cleaning. The dispersion effect is also revealed showing that the oil slick diameter is moderately increasing during this time. Figure 7 depicts more clearly the spreading of the oil spill and its statistical center of gravity movement (noted by red circles), as it is advected for many kilometers in Thermaikos bay. The oil slick dispersion after 3 hrs is shown in dark blue color, whereas the spreading of the plume at the end of a day (after 24 hrs) is in light blue. The spiral movement of the oil’s statistical center of gravity in Figure 6, near the shore, may be attributed to the inertial currents frequently observed in the sea, and points out that the largest part of the oil reaches the shore (beaching) within the initial time steps (after 9 hrs). The model results show that a contingency emergency plan should be, in this case, activated in less than 12 hrs, in order to avoid the social, economic, and environmental damage in the area.

3.2 DIAVLOS-NASOS

The universities of Thessaloniki and Athens cooperated, in the framework of the project DIAVLOS [28, 29, 30], to produce a simulation model that could trace an oil spill in case of an accident in the Alexandroupolis Gulf, which is considered of high risk, after the construction of the Burgas-Alexandroupolis oil pipeline that has been signed by the Greek, Bulgarian, and Russian governments (on March 15, 2007), and is expected to carry huge oil quantities (35–50 Mtn/yr) that will be loaded to tankers and further transported by sea to international markets.

3.2.1 Model construction

A 48-hr oil spill dispersion forecasting system was developed and implemented in the wider northern Aegean Sea. The system was based on wind, wave, and ocean circulation models, coupled with the operational systems ALERMO and SKIRON of the Department of Physics of the University of Athens, and an oil spill dispersion model by the Department of Civil Engineering of the Aristotle University of Thessaloniki. The basic components of the model are presented hereafter.

3.2.1.1 Meteorological model

The meteorological model is part of the weather forecasting system SKIRON developed by the University of Athens, which is based on a Limited Area Model of Eta/NCEP, which provides high-resolution (of 1°) weather forecasts for 120 hrs [31]. In this work, the non-hydrostatic model Eta of the operational forecasting system SKIRON has been used, with an even higher resolution (of 1°/10) since the description of non-hydrostatic phenomena is necessary for the model to be applied for such a high resolution. The available data are: air temperature and humidity at 2 m, wind velocity and direction at 10 m, short- and long-wave radiation, atmospheric pressure at the sea level, and precipitation.

3.2.1.2 Hydrodynamic model

The hydrodynamic forecasting model, which was developed for the area under study, is based on the well-known numerical model Princeton Ocean Model (POM) [21], a model that has been widely used for the simulation of open and coastal sea circulation. The Department of Physics at the University of Athens developed a POM-based sea circulation forecasting model, named ALERMO, under the European programs MFSPP and MFSTEP, which was made public in 2004 (http://www.phys.uoa.gr). Its resolution is 2 nm and it is coupled with the Mediterranean forecasting system MFS-OGCM (INGV, Italy), which also assimilates satellite data and field measurements (XBTs, CTDs, Profiling Float). ALERMO uses an innovative variational initialization technique (VIFOP) for its initial conditions and provides initial and boundary conditions to other coastal forecasting systems (e.g., Greece, Cyprus, and Israel). The developed hydrodynamic model has been applied to the area of interest (North Aegean: 38.7°–41.1° N, 22.5°–27.1° E) with horizontal resolution of 1°/60 (1 nm) and 25 σ-layers in the vertical. The open boundary conditions, as well as the initial conditions, are provided by ALERMO. Finally, the momentum, heat and water fluxes at the air-sea interface are calculated by using the weather forecasting data provided by the meteorological model. A series of sensitivity tests have been conducted, with respect to the atmospheric conditions and the inflow of water from the Black Sea [32].

3.2.1.3 Wave model

The wave model WAM [33, 34] has been adapted and applied to the area under study, with a resolution of 1°/20 by the University of Athens. This model is a third-generation numerical model, which has been widely used, as it best describes the evolution of the wave spectrum in the sea. WAM makes a distinction between deep and shallow waters, depending on the sea depth at the area, where the wave equations are being applied.

3.2.1.4 Oil spill dispersion model

The oil spill transport and dispersion model, named NASOS (North Aegean Sea and Oil Slick), is based on the 3D model that has been developed by the civil engineering department of the Aristotle University of Thessaloniki [16, 17]. The produced operational code was updated and adapted to the area of interest (38.7°–41.1^° N, 22.5°–27.1° E), that is, the North Aegean Sea. This model assumes that the oil slick is described by a big number (equal to 10⁴) of particles (or “parcels”), each of which represents many cubic meters of oil, and has individual physicochemical properties. The distinction of the oil spill in these “parcels” is supported by the fact that oil contains a variety of hydrocarbon components with different physicochemical properties, and therefore makes the monitoring of the evolution of each “parcel” easier. This model recognizes four different oil types, enough to cover the variety of physicochemical properties and their effect on the environment. In terms of the simulation of the fate and transport of the oil spill in the sea, the model recognizes the following processes: initial spreading and transport, horizontal and vertical diffusion and dispersion, evaporation, emulsification, beaching, and sedimentation. The Mackay approach [26, 27] is again used for the computation of the evaporated volume fraction from each oil component. Emulsification is activated for oil fraction with densities and wave height/length ratios beyond critical values. The evaporation and emulsification processes lead on the one hand to the deficit of contained oil in the evaporated particles, and on the other hand to the creation of a mixture of water and oil, and formation of a floating mousse. Beaching, which expresses the entrapment of oil on the beach, is expressed by the duration of particles trapped on the coastal boundary. It strongly depends on the coastal morphology (rocky to flat sandy beach).

3.2.2 Operational application

This study’s objective was to produce an operational tool to forecast the fate of a possible oil spill in Alexandroupolis gulf and contribute to manage environmental crises from such pollution, after the Burgas-Alexandroupolis pipeline would be constructed. It was clearly very important to make this oil spill dispersion forecasting system user-friendly and easily accessible by all the involved authorities (e.g., the Ministry of Mercantile Marine/Oil Spill Prevention Department, Prefecture of Evros, local port authorities, etc.).

The required input data are the initial time and location of the oil spill, the oil quantity, and quality characteristics (oil type). These data can be introduced to DIAVLOS forecasting system’s Website on-line by the user. The information is processed by the wave-atmospheric-hydrodynamic forecasting system. The oil spill dispersion model runs for 48 hrs and provides hourly images in less than 3–4 minutes. The system can also provide an updated forecast if updated input information is available (in less than 48 hrs). The model computes the hydrodynamic flow and the 3D oil spill evolution and returns to the user the following information (in percentages): surface and subsurface dispersion, evaporation, emulsification (or water content), and beaching of the oil. The output of the model, as Figure 8 shows, contains a graphical representation of the hypothetical oil spill accident area, as well as the horizontal dispersion of the oil, depicted by different colors (surface, subsurface, or beached). Additionally, the University of Athens provides 60-hr high-resolution atmospheric and oceanic forecasting at the region, to facilitate operations to contain oil spill spreading, and beaching and cleaning operations. Wind, wave, and oceanic circulation forecasts are available on DIAVLOS Website by the Ocean Physics and Modelling Group and the Atmospheric Modelling and Weather Prediction Group of the University of Athens. The password-protected and interactive Website (http://diavlos.oc.phys.uoa.gr) was created and offered to public in May 2008, combining efficiency and simplicity.

Figure 8.

Graphical and statistical results of the oil spill operational transport model (DIAVLOS).

The system was also tested against field observations of special drifting floats that monitor the fate and evolution of oil spills. These “smart drifters,” were especially designed for this project by MARAC Electronics Co. (Greece) in cooperation with the University of the Aegean, co-flowing and transmitting along with the oil slick. The model proved to give satisfactory results and in most cases the forecasting error is quite small, allowing the operational use of the system.

3.3 BSB Net-Eco—SLICKNEW

This latest study utilizes another upgraded version of PARCEL and 3D Sea and Oil Slick (SOSM) models, developed by the Aristotle University of Thessaloniki [16, 17]. The newly formalized model was tested in the Aegean Sea [35] and adjusted to the characteristics to the Black Sea, and more particularly to the Azov Sea [36].

3.3.1 Model construction

The model suite comprises of a Lagrangian (tracer) model for the transport and physicochemical evolution of an oil slick [37]. The input requirements of the model are surface wind velocities, air temperature, vertical and horizontal diffusion coefficients, surface currents, wave characteristics (height, period, direction), and, in terms of the oil transport, the initial coordinates of each parcel, its initial volume and mean density and droplet diameter, as well as evaporation rate and parameters relevant to the oil type and identity. The POSEIDON system [19], which is utilized in this modeling effort, provides information on wind speed and direction, atmospheric pressure, air temperature, wave parameters, current speed and direction, water temperature, salinity, dissolved oxygen, chlorophyll-a, and radiation. The observational basis of the system is a network of oceanographic 11 buoys (7 SEAWATCH, 3 Wavescan and 1 deepwater SEAWATCH module), operating in the Aegean Sea since June 1999. The suggested oil spill model utilizes the POSEIDON wave and wind datasets, to produce the sea velocity fields due to currents and waves. The wave, generated near surface, velocity field is computed from the classical Stokes’ drift formulae [20], based on the local values of wave height H_s and wave period T_s.

$$U_m\left(z\right)={\left(\frac{\pi ·H_s}{L_0}\right)}^2\frac{C_0}{2}\frac{cosh(2k\left(H_s+z\right)}{sinh{\left(k·H_s\right)}^2}$$E4

where $L_0=g·T_s^2/2\pi$, $C_0=g·T_s/2\pi =L_0/T_s$, and $k=2\pi /L_0$ (π = 3.14). The horizontal diffusion coefficient D_h is estimated adopting the Smagorinsky formula [21, 38] as follows:

$$D_h=C·\mathit{dx}·\mathit{dy}·{\left[{\left(\frac{\partial u}{\partial x}\right)}^2+0.5{\left(\frac{\partial v}{\partial x}+\frac{\partial v}{\partial y}\right)}^2+{\left(\frac{\partial v}{\partial y}\right)}^2\right]}^{0.5}$$E5

where C the horizontal diffusion coefficient. This equation is used to estimate the velocity of each particle, if it is selected via a random number from a sample following the uniform distribution over a range {−Ur, +Ur}

$$U_r=\sqrt{\frac{6·D_h}{\mathit{\Delta t}}}$$E6

The generation of that sample of random velocity values is based on Monte Carlo sampling, a very common and powerful procedure in simulation theory [37]. The vertical displacement of the oil “particles” is considering the vertical diffusion due to currents and waves [39]:

$$D_V=D_{\mathit{Vc}}+D_{\mathit{Vw}}$$E7

$$D_{\mathit{Vc}}=\lambda ·h·W$$E8

$$D_{\mathit{Vw}}=0.028·\frac{H_s^2}{T_s}·e^{\frac{4\mathit{\pi z}}{L_0}}$$E9

where λ = 0.001; h the water depth; and $W=\sqrt{W_x^2+W_y^2}$ for W_x, W_y the given wind velocities in the x, y axes, respectively. Consequently, the vertical velocity due to buoyancy and diffusion of the oil “particles” is given, similar to the horizontal velocity [37], by

$$W_r=\sqrt{\frac{6·D_V}{\mathit{\Delta t}}}$$E10

With respect to the input requirements of the slick model, it also requires bathymetry data of the selected area, and a file containing the characterization of the coastal meshes, according to their oil-holding capacity and the open sea boundaries. Based on all that, new horizontal positions of the oil “particles” are estimated; some are “trapped” on the beach, others may be vertically displaced due to buoyancy and diffusion; a fraction of heavy classes of oil may be emulsified over a certain wave curvature, whereas a light oil fraction in sea or on coast may be evaporated. Among the oil weathering processes that take place in the sea is evaporation. It affects the surface oil particles, in sea or on coast. A complete review of various approaches in estimating the evaporated oil is presented in [40]. Thus, adopting the empirical equation of oil evaporation representative of the oil type “Barrow Island, Australia,” the evaporation formula used in this model is described by the following equation:

$$\%\mathit{Ev}=\left(4.67+0.045·T\right)ln\left(t\right)$$E11

where T is the air temperature (°C) and t the time (in minutes). Another important oil weathering process is emulsification, which is expressed as a function of the wind velocity W expressed as $W=\sqrt{W_x^2+W_y^2}$ and the temperature T [26]:

$$\mathit{Em}=\frac{1-e^{-0.0000056\left(1+W^2\right)·T}}{1.25}$$E12

Different oil products react differently to these processes. Lighter oil fractions tend to evaporate, whereas heavier fractions tend to emulsify. Therefore, the model, taking that into account, can simulate different oil types according to their density and buoyancy velocity. The processes of photo-oxidation and biodegradation are not considered in this version of the model, as their effects are more significant at a later stage of a spill’s evolution (see Figure 3). All particles are considered to spread at a single location, while they can be released all at the same time (instantaneous discharge), or in sequence over a given period of time (continuous discharge of specified duration).

3.3.2 Operational application

The produced operational code was adapted to the environment of a coastal basin of NW Aegean Sea (39.96°–40.66° N, 22.50°–23.40° E), namely Thermaikos bay, which also includes Thessaloniki’s gulf [41]. The area that has been selected to be studied is shown in Figure 9 [42], with a horizontal discretization of 1°/30 (Dx = 2800 m, Dy = 3700 m). The commercial port of Thessaloniki provides a good study area, as it is enclosed in Thermaikos bay with an important coastline and a busy oil terminal. Although the terminal operates under strict IMO (International Maritime Organization) and EU (European Union) regulations for oil terminal operational safety, the port area is a potential oil slick accidental source [38].

Figure 9.

Thermaikos bay and ships’ route to and from the port of Thessaloniki.

The transport and fate of the slick are closely related to the hydrodynamic conditions in the bay, the bay morphology and the oil properties. The hydrodynamic circulation in the bay is sustained mainly by the strong NW to NE winds and the circulation patterns are regulated by the bay bathymetry and the coastal topography (Figure 10). Following a thorough consideration, the period of data selected to be used from the POSEIDON dataset was a winter week of January (21–28/01/2014) which exhibited significant winds, able to transport the oil slick and produce substantial results. The data were given for 29 time steps (6 hrs. each). The POSEIDON dataset for the Aegean Sea offers 3D hydrodynamic forecasts (POM) nested to POSEIDON Med forecast, for the depth levels of 5, 15, 30, 50, 80, and 120 m. The depth in the area under study, as shown in Figure 10, is less than 100 m. One hypothetical oil spill scenario has been selected to portray the model’s simulation capabilities and results, at the east side of Thermaikos bay, near the touristic areas of Chalkidiki (Figure 11). An oil spill of 5000 parcels, representing an estimated total volume of 10,000 tn, was hypothetically released instantaneously on the sea surface. Total simulation run was for 180 hrs (7.5 days). Figure 11 presents snapshots of the spatiotemporal evolution of the oil slick at a 36-hr step after the release. The touristic coast of Chalkidiki is damaged within some hours after the spillage (clearly after 2–3 days). The dispersion effect is also revealed showing that the oil slick diameter is moderately increasing during this time.

Figure 10.

Thermaikos bay bathymetry.

Figure 11.

Snapshots of the spatiotemporal evolution of an oil spill hypothetically released near the coast of Chalkidiki, at the east side of Thermaikos bay, after: (a) 36 h, (b) 72 h, (c) 108 h, (d) 144 h, and (e) 180 h.

Acknowledgments

Part of this research has been funded by the General Secretariat of Research and Technology (GSRT) of Greece, under the title “Development and application of crisis management tools aiming at the prevention from oil spill pollution in the area of Burgas – Alexandroupolis pipeline (N.E. Aegean Sea) – DIAVLOS.” Another part of this research has been produced with the assistance of the European Union, through the BSB JOP 2013-2017 (under the project with the acronym “BSB Net-Eco”).

Thanks

I would like to express my gratitude to all my colleagues who contributed to this work, my co-authors to the published articles and conference papers, but most of all to Emeritus Professor Ch. Koutitas.