Determination of the oil spill removal area by oil particle tracking in a harbor

A hydrodynamic model is the foundation of trajectory prediction of oil spill particle. The widely used EFDC model is used to provide the dynamic flow field and oil particle movement. This model is developed by John Hamrick (Hamrick 1992, 1995; Hamrick & Wu 1997). It is one of the most respected models and has been widely used by universities, the government and environmental consulting in a wide variety of fields (Shen and Haas 2004; Xu et al. 2008; Gong et al. 2009; Wang et al. 2010; Gong et al. 2012; Zhao et al. 2013). In this model, turbulence is included in the second moment turbulence closure model developed by Mellor and Yamada (1982) and modified by Galperin et al. (1988). Turbulent viscosity and diffusivity is related to the turbulent intensity, the turbulence mixing length and the Richardson number. The controlling equations of the EFDC model can be found in Hamrick (1992).

The macroscopic motion of oil spill particles is simulated by the particle tracking model of the Lagrange method. The particle diffusion process caused by turbulence is random. It is simulated by the random walk method, namely the motion of particles in a turbulent flow field is similar to the random movement of molecules.

The particle tracking model with the random walk of particles is written as follows:

where, t is time, (x,y,z) is Lagrangian coordinates of a particle, V→=(u,v,w) ${\overrightarrow{V}} = (u,v,w)$ is flow velocity, u, v and w stand for the velocity in the direction of x, y, z, respectively. D_H and D_z are the horizontal and vertical diffusion coefficients, respectively, dt is time step and p with the mean value of 0.5 is a random number from a uniformly distributed random variable generator. The above equations are solved by the fourth-order Runge-Kutta method.

Tianjin Nangang is selected as the study area. It is an important part of the Tianjin port. It is located in the east of Bohai Bay (Fig. 1). The development and construction of the Tianjin Nangang industrial zone is an important strategy of the Tianjin city. The eastern harbor reclamation project has been established and the removal of oil spills leaked by ships is a key issue in preserving the good quality of water and in the development of tourism. The study area is surrounded by the north and east breakwaters of Nangang. The wave action is very weak in sheltered waters and the impact of surface waves is not taken into account.

Figure 1

A multi-nested grid model was used to perform an accurate analysis. It consists of the Tianjin port model and the Nangang model presented in Figs 2 and 3. In the computational domain, the seabed is almost flat and the depth of water from the shore gradually increases to about 20 m at the open-sea boundary. The ports and waterways are about 15 m deep.

Figure 2

Figure 3

The Tianjin port model has 49221 horizontal grids. The sufficient grid resolution is provided to effectively present the coastline as well as sufficient details of a specific area to support the detailed modeling. Curvilinear grids are used in the domain, and the spatial resolution ranges from 50 m in the nearshore zone to 1000 m in the far-shore zone. In order to ensure the accuracy of the surface wind flow and the bottom friction stress, three layers are selected in the vertical direction.

The initial conditions of the Tianjin port model are water level and flow velocity. As the water level and flow velocity follow rapidly towards the external dynamic response, the cold start mode sets the initial water level as the mean water level and the initial velocity to zero. The water elevation at the open boundary is given by the nearest neighbor interpolation. The area adjacent to the sea is mostly the muddy coast and the seabed is smooth. After several debugging processes, the roughness is selected as 0.022 for the best fit of data. The horizontal momentum diffusion coefficient is calculated according to Smagorinsky’s formula and the mixed constant level is 0.1. The horizontal eddy viscosity coefficient is the same with the turbulent diffusivity of momentum. The model time step is 2 s based on the CFL condition, which is determined by the grid dimension and the velocity in the adjacent sea area. For a model run of 100 days, the computer runtime is about 4.5 days on a Dell T5500 workstation with 12 processors.

The nested Nangang model has 17522 horizontal grids with the size of 50 m. The same parameter settings, such as the vertical roughness, bottom roughness and the time step, are used for the Tianjin port model. The horizontal momentum diffusion coefficient is also calculated by the Smagorinsky formula. The water level is provided by the results of the Tianjin port model. The initial simulation starts with calm water conditions without a wind feld. The time step of the nested model is also 2 s. For a run of 100 days of the nested model with tracing three thousand oil particles, the computer runtime is about 5 days on a Dell T5500 workstation with 12 processors.

The water level and the magnitude and direction of velocity are compared with the feld monitoring data. Stations monitoring the water level and flow velocity are shown in Fig. 1. Figs 4-6 show that the numerical results are consistent with the measured data during low and high tide periods. Furthermore, the numerical results are also in line with the measured data during low and high tide periods at the two other velocity monitoring stations, which are not shown in this paper for the sake of brevity. It is evidenced that our numerical model can well simulate the flow field.

Figure 4

Figure 5

Figure 6

In order to determine the optimal removal area of an oil spill, the particle tracking model is introduce to simulate the trajectory of oil particles based on the flow field of the Nangang model. The possible oil spill area is divided into six regions based on the characteristics of a harbor shown in Fig. 3. Regions I, III and V correspond to a coal or ore berth on the west side of the harbor. Regions II, IV and VI correspond to an LNG reservoir, a petrochemical berth and a universal berth on the east side of the harbor. The drift process of an oil spill is described by a movement of oil particles. Particles are put in the center of each grid of a region to represent the movement of an oil spill within a given region. At each time step, the new position of oil particles is computed by considering the transport induced by currents, wind and turbulent dispersion.

Fig. 7 shows the trajectory pattern of oil particles in different regions. Particles from regions I, II and III have a wide sweeping area. In addition to moving in the offshore direction, some of them can even move toward the west basin. Particles from regions IV, V and VI are exchanged mainly offshore. A small part of particles from region III and a few particles from regions I, II and IV migrate to the bottom of the basin. Particles from regions I, II, III and IV migrate from the harbor in a short time while those from regions V and VI take a long time. This is due to the fact that regions V and VI are characterized by low water exchange capacity.

Figure 7

It has been found that many of the particles pass through the south breakwater before they move offshore. Thus, it makes sense that the oil removal area is located in the vicinity. In this approach, other possible areas are also selected to study the effect of oil removal in different oil spill areas.

By considering the trajectory of oil particles, the sweeping area of an oil spill is determined for a given removal area. The position of different removal areas is shown in Fig. 8. By comparing the oil removal rate of different areas, the area with the highest efficiency in oil removal is selected.

Figure 8

According to the trajectory of particles, six areas are selected for the potential location of the oil removal area. In this method, particles moving through the removal areas are assumed to be removed completely from the computational domain and the number of removed particles is counted by a Matlab program which can acquire data on the movement of particles from the particles tracking models in EFDC. The ratio of the removed particles to the total particles is referred to in this paper as the removal rate.

Each of oil removal areas 1, 2, 3 and 5 has a length of 300 m and a width of 50 m. Oil removal area no. 4 has a length of 50 m and a width of 300 m. Oil removal area no. 5 is a combination of oil removal areas no. 3 and 4.

Of the above 6 oil spill removal areas, particles from region III have the largest sweep range and they are used as an example to show the oil removal result of different oil removal areas. Fig. 9 shows the trajectory of oil particles from region III after 100 days of the oil removal operation for different oil removal areas. The oil removal rate of oil removal areas 2, 3, 4 and 6 signif cantly reduces the oil particles from regions III. The oil removal rate of removal area no. 6 has the largest oil removal rate.

Figure 9

Removal area no. 5 is located in the northwest of the harbor. If this area is selected as an oil removal area, only oil particles from regions I, II and III moving to the west are removed and those from regions IV, V and VI are removed in small numbers. Therefore, the oil removal rate of removal area no. 5 is low.

Oil removal area no. 1 near the north breakwater is the least efficient oil removal position. It can only remove a few oil particles from regions I, II and III and it cannot remove the oil spill from regions IV, V and VI.

The study area is located in the monsoon climate zone. The wind shows different characteristics in different seasons. The wind has a strong effect on the movement of an oil spill. Table 1 shows the direction, speed and probability of prevailing wind in each season. It has been found that the probability of prevailing wind in each season is about 15%. This means that the prevailing wind occurs in about 60% of the days in a year. Therefore, the seasonal prevailing wind is considered as the simulation wind of the season.

The oil removal rate is defined as the ratio of the number of removed oil particles to the number of initial total oil spill particles in the harbor. Table 2 shows the oil removal rate at different removal areas after 100 days of the oil removal operation. It has been found that oil removal areas no. 1 and 5 always show poor performance. Removal area no. 2 is easily affected by wind. Its removal rate increases in spring and summer; however, it significantly decreases in autumn and winter, especially in winter. Compared to no-wind conditions, removal areas no. 3, 4 and 6 with a seasonal wind show high removal rates. The removal ratio of removal areas 3, 4 and 6 is about 80% throughout the year. Considering the fact that removal area no. 6 overlaps with removal areas no. 3 and 4 and removal area 3 is on the cross-section of the harbor gate, removal area 4 is selected as an optimal area for the oil spill removal in the harbor.