`Journal of Applied MathematicsVolume 2012, Article ID 681585, 15 pageshttp://dx.doi.org/10.1155/2012/681585`
Research Article

## Numerical Simulation of Oil Spill in Ocean

1Department of Civil and Environmental Engineering, Hanyang University, 222 Wangsimni-ro, Seongdong-gu, Seoul 133-791, Republic of Korea
2Department of Civil Engineering, Kyungnam University, 449 Wolyoung-dong, Masan 631-701, Republic of Korea

Received 4 January 2012; Revised 29 March 2012; Accepted 2 April 2012

Copyright © 2012 Yong-Sik Cho et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The spreading of oil in an open ocean may cause serious damage to a marine environmental system. Thus, an accurate prediction of oil spill is very important to minimize coastal damage due to unexpected oil spill accident. The movement of oil may be represented with a numerical model that solves an advection-diffusion-reaction equation with a proper numerical scheme. In this study, the spilled oil dispersion model has been established in consideration of tide and tidal currents simultaneously. The velocity components in the advection-diffusion-reaction equation are obtained from the shallow-water equations. The accuracy of the model is verified by applying it to a simple but significant problem. The results produced by the model agree with corresponding analytical solutions and field-observed data. The model is then applied to predict the spreading of an oil spill in a real coastal environment.

#### 1. Introduction

In an ocean environment, a spilled oil due to an unexpected accident may cause a serious damage. Thus, an accurate prediction of behaviors of split oil is very significant to keep coastal environmental system. The area of oil spreading can be predicted numerically by solving the proper equations governing the flow field and the associated mass transport phenomenon. The most reasonable choice is probably an advection-diffusion equation consisted of both advection and diffusion. The equation is difficult to solve numerically because it is represented by both hyperbolic- and parabolic-type partial differential equations. A standard split operator approach is a plausible and practical choice for solving the advection-diffusion-type partial differential equation [1]. The hyperbolic (advection) and parabolic (diffusion) components of the equation are solved separately using methods that properly and simply simulate the physical behaviors. Numerical solutions obtained from the advection and diffusion equations are then combined.

An oil spill accident is regarded as a kind of a disaster because it causes not only fatal destruction of the marine environment but also enormous cost of the disaster prevention and the damage compensation. Thus, it is important to accurately predict the spread range of the spilled oil as an early stage countermeasure against a disaster. In the last three decades, many investigators have studied the transport processes of oil spills based on the trajectory method [24]. Those methods have been applied in river-lake system [57] and seas [811]. Some commercial oil spill models, such as, COZOIL [12], NOAA [13], OILMAP [14], WOSM [15], have been used to determine the oil movement and distribution in the water body. However, transport of oil spills has been conducted considering tidal currents simultaneously in only few researches. Furthermore, oil spills processes in a field of strong tide and tidal currents have not been investigated.

In this study, the spilled oil dispersion model has been established in consideration of tide and tidal currents simultaneously. The Hebei Spirit oil spill that occurred on December 7, 2007 is the largest oil spill accident occurred in the Yellow Sea. In the Yellow Sea, tidal currents are very strong and should be considered to investigate related coastal processes. Therefore, the accuracy of the model for predicting tide and tidal currents is very important to investigate oil spreading in this area. Subsequently, by computing the diffusion distribution on the Hebei Spirit oil spill considering tidal currents simultaneously, verification has been made through comparison of the diffusion distribution and the field-observed data between at 8 p.m. on December 7, 2007 when the spilled oil flowed into the whole Manripo and Sindu-ri seashore and at 11 a.m. on December 11, 2007 when the satellite photographs exist. Also, this study aims to prepare for a possible accident in the future and provide the basic materials for establishing the disaster prevention of the oil spill pollution.

In the following section, the governing equations are described first. The numerical model is then presented and detailed description is followed. Numerical simulations of tide and tidal currents are conducted to test the applicability of the spilled oil spreading model in the Yellow Sea, and the verification has been made comparing numerical results with the field observed data. Also, the model is verified conducting the numerical calculations on the simple diffusion distribution problem. The model is then applied to predict the Hebei Spirit oil spill in a real coastal environment. Finally, concluding remarks are made.

#### 2. Governing Equations

The behaviors of tide and tidal currents may be described by the following nonlinear shallow-water equations [16, 17]: In (2.1)–(2.3), represents the free surface displacement and is the total water depth defined as with being a local still water depth (Figure 1), and are volume flux components along the - and -axis directions defined as and with and being the depth-integrated velocity components in the - and -axis directions, and and are bottom frictional effects of the - and -axis directions, respectively.

Figure 1: A definition sketch of the free surface displacement and local water depth.

The study is focused to investigate the oil spreading in a limited area in the Southern part of the Korean Peninsula and thus the Cartesian coordinate system is then employed. The boundary conditions along the offshore are provided by a large-scale numerical model [18]. The bottom frictional effects may play a significant role in the very shallow zone. Thus, the Manning’s empirical formula is employed in the study given as: In (2.4), is the Manning’s roughness coefficient. In this study, the Manning coefficient of has been used.

It is well known that the change of sea level is very complicated and affected by many physical phenomena such as waves, tides, and others. Among these, the tide may be the most dominant force for the oil spreading in the Yellow Sea. In this study, therefore, the tide is only investigated in a very limited area without considering any interaction with other physical phenomena.

In this study, a source or sink can be described to a following equation [19]: in which represents the concentration of the oil in the sea. The source or sink is assumed to be nonlinear and is linearized using the Newton-Raphson iterative technique [20] as: in which the superscript () represents a previous step and the sedimentation rate, can be derived as the following relation [21]:

By applying an operator splitting approach to (2.1), the following advection equation can first be obtained by: and the diffusion-reaction equation can then be derived by in which , , , represent the diffusion coefficients.

#### 3. Numerical Model

The continuity and momentum equations given by (2.1)–(2.3) are solved by the finite difference approximations. The detailed description of the finite difference method can be found in literatures [2, 3] and is not repeated here. Firstly, the advection equation can be discretized by using the SOWMAC (second-order wave equation method for advective calculation) scheme based on the characteristic method and is given as [22] in which is the Courant number defined as , and superscript and subscript are time level and computational point, respectively. and are weighting factors of the finite difference scheme.

An unnecessary wave propagates upstream if the initial concentration distributions at both and are not consistent with the downstream advection progress. To overcome this difficulty, the values of the concentration at time step can be employed as shown in Figure 2 instead of time step . Drawing the characteristic curves which go downstream through points and , time step can be introduced corresponding to time when each characteristic curve intersects the (-)axis or -axis. According to the concept of characteristics for a pure advection, the concentrations are conserved between two time steps and are given as follows:

Figure 2: Estimation of concentration at a previous time step.

Using (3.2), (3.1) can be rewritten as The optimal values of and in (3.3) do not seem to be constant, but instead are dependent on the Courant number . To obtain the functional relationship of and to , the Taylor series expansion on (3.3) may be employed [22].

The final expression of the advection equation can be written as follows: in which the coefficients are given by The values of , , and are defined as: in which AINT is one of the intrinsic functions in FORTRAN that carries out the function of truncating decimals.

Secondly, the diffusion-reaction equation can be written as: The diffusion-reaction equation is discretized by a three-level locally implicit scheme with reference to the computational grid, as shown in Figure 3 [23].

Figure 3: Computational mesh for the diffusion-reaction term.

Equation (3.7) contains two dimensions and is different from the form derived in Hobson et al. [23]. However, their numerical scheme can be easily extended to the -direction in a similar manner and the final form of discretized equations using the finite difference approximation is given as: in which is used and the diffusion coefficients are assumed as a constant value, . In (3.8), and are the artificial values which are determined to give the correct amount of diffusion and avoid any time lag problems. In general, a whole family of diffusion schemes can be obtained by changing the values of and . The optimum scheme of this type is the one which exhibits no time lag and is the most numerically stable.

In this study, (3.8) is solved by a two-dimensional three-level locally implicit scheme with second-order accuracy [23]. Equation (3.8) can be rewritten for as follows: in which , , , and are given as:

Since a first-order sink and zero-order reaction cannot generate numerical stability alone, this study is concerned only with the presence of a first-order term. Thus, the zero-order sedimentation rate, , is assumed to be zero, and the first-order sedimentation rate of a suspended material, , is described as: and a first-order sink is defined as: in which is a coefficient for calculation of the concentration at the bottom and, for convenience, is assumed to be 1, and is the sedimentation rate, assumed to be 0.08 mm/sec.

#### 4. Verification of the Numerical Model

The numerical model developed in this study is employed to simulate the change of tide and tidal currents near the accident area of the Hebei Spirit oil spill occurred on December 7, 2007. The Hebei Spirit oil spill is the largest oil spill accident occurred in the Yellow Sea located between Korea and China. In general, tide effects play a dominant role in coastal processes in the Yellow Sea, which is well known for the strongest tidal currents in the world. Therefore, the accuracy of the model for predicting tide and tidal currents is very important to investigate oil spreading in this area. To verify the model, free surface elevations and tidal ellipses at several points are computed and compared with the corresponding observed data.

Figure 4 shows the computational domain and sea level observation points. Sea bottom topography is also shown in the figure. Numerical simulations of tide and tidal currents were conducted during 16 days (December 7 ~ 22, 2007) and the harmonic analysis for computed results is presented. Table 1 displays a comparison of observed data and computed results using harmonic analysis. Computed results represent both amplitudes and phase differences very well. Figure 5 shows a comparison of numerically obtained free surface profiles and corresponding observed measurements at points of T-1, T-2, and T-3. The agreement between the numerical solutions and field measurements is quite reasonable.

Table 1: Observed data and computed results using harmonic analysis.
Figure 4: Computational domain and sea bottom topography.
Figure 5: Comparison of observed and computed free surfaces of tides.

Tidal currents are very strong in the Yellow Sea and those may play an important role in oil spreading. In Figure 6, tidal ellipses at points of C-1 and C-2 are computed and compared with field-observed data to verify the accuracy of the numerical model. The numerical model represents well-strong tidal currents observed in this area.

Figure 6: Comparison of observed and predicted tidal ellipses.

#### 5. Application of the Model

The spilled oil dispersion model as established in previous sections is applied to predict the diffusion distribution on the Hebei Spirit oil spill. Subsequently, by conducting the numerical calculations on the diffusion distribution considering tide and tidal currents simultaneously, the prediction has been made through the comparison of the diffusion distribution observation between at 8 p.m. on December 7, 2007 when the spilled oil flowed into the whole Manripo and Sindu-ri seashore and at 11 a.m. on December 11, 2007 when the satellite photographs exist (Figure 7). In Figure 7, a target area of the numerical simulation is also described.

Figure 7: Satellite image of the Hebei Spirit oil spill (December 11, 2007).

Figures 8(a) and 8(b) show distribution of simulated tidal currents field, 8(a) flood phase and 8(b) ebb phase, where the Hebei Spirit was stranded on duration of a spring tide. The results show that the tide current reaches the maximum speed in the northern area from the accident and mainly flows not perpendicular to the coastline but parallel to that. As a result, strong invasion of the oil spill to the coastline is not represented in the numerical simulation.

Figure 8: Distribution of tidal currents near the target area.

Figures 9(a) and 9(b) show numerically obtained diffusion distribution at 8 p.m. on December 7, 2007 and 11 a.m. on December 11, 2007. Figure 9(a) shows diffusion distribution after 13 hours since the Hebei Spirit had been stranded. The south-western currents on a flood phase and north-eastern currents on ebb phase are well reproduced and those are well reflected in the shape of oil distribution. Therefore, the numerical result shows similar distribution of the oil spreading with the satellite photograph (Figure 7). Since there is lack of quantitative observed data near the target area, numerical results are compared with a satellite image only in a qualitative viewpoint. On the other hand, by comparing the satellite photograph and Figure 9(b), numerical results do not represent strong invasion of the oil spill to the coastline. It may be occurred because effects of wind wave were neglected in the numerical simulation. Nevertheless, it is obvious that the model represents well the whole processes of the oil spreading, and the model can be applied to predict the oil spreading in a real coastal environment.

Figure 9: Computed diffusion distribution of the Hebei Spirit oil spreading.

#### 6. Concluding Remarks

In this study, a numerical model is employed to solve a two-dimensional advection-diffusion-reaction equation. The model is based on a standard split operator (fractional step) approach. Thus, the hyperbolic (advection) and parabolic (diffusion) portions of the equation are solved separately by using techniques describing properly the physical behaviors of each. In the model, the advection step is solved by using the SOWMAC scheme, while the diffusion-reaction step is done by a three-level locally implicit scheme.

The numerical model is first applied to an idealized problem to verify its accuracy. The model produces results agreeable with corresponding analytical solutions and field observed data. The model is then employed to investigate the behavior of the oil spreading in a coastal environment. The model yields reasonable results. The employed model could be used to forecast the behaviors of oil spreading in various practical situations.

#### Acknowledgment

This research was supported by Basic Science Research Program through the National Research Foundation of Korea(NRF) funded by the Ministry of Education, Science and Technology (2011-0015386).

#### References

1. B. Lin and R. A. Falconer, “Tidal flow and transport modeling using ultimate quickest scheme,” Journal of Hydraulic Engineering, vol. 123, no. 4, pp. 303–314, 1997.
2. D. Mackay, S. Paterson, and K. Trudel, “A mathematical model of oil spill behavior,” Environmental Protection Service, Fisheries and Environment Canada.
3. J. C. Huang, “A review of the state-of-the-art of oil spill fate/behavior models,” in Proceedings of the Oil Spill Conference, pp. 313–322, Washington, 1983.
4. M. L. Spaulding, “Oil spill trajectory and fate modeling: state-of-the art review,” in Proceedings of the 2nd International Oil Spill Research and Development Forum,, London, UK, 1995.
5. H. T. Shen and P. Yapa, “Oil slick transport in rivers,” Journal of Hydraulic Engineering, vol. 114, no. 5, pp. 529–543, 1988.
6. H. T. Shen, P. D. Yapa, and M. E. Petroski, “Simulation of oil slick transport in Great Lakes connecting channels,” Report Nos. 86-1 to 4, vols. I-IV, Department of Civil and Environmental Engineering, Clarkson University, Potsdam, NY, USA, 1986.
7. P. D. Yapa, H. T. Shen, and K. S. Angammana, “Modeling oil spills in a river-lake system,” Journal of Marine Systems, vol. 4, no. 6, pp. 453–471, 1994.
8. S. A. Lonin, “Lagrangian model for oil spill diffusion at Sea,” Spill Science and Technology Bulletin, vol. 5, no. 5-6, pp. 331–336, 1999.
9. X. Chao, N. J. Shankar, and H. F. Cheong, “Two- and three-dimensional oil spill model for coastal waters,” Ocean Engineering, vol. 28, no. 12, pp. 1557–1573, 2001.
10. S. D. Wang, Y. M. Shen, and Y. H. Zheng, “Two-dimensional numerical simulation for transport and fate of oil spills in seas,” Ocean Engineering, vol. 32, no. 13, pp. 1556–1571, 2005.
11. S. D. Wang, Y. M. Shen, Y. K. Guo, and J. Tang, “Three-dimensional numerical simulation for transport of oil spills in seas,” Ocean Engineering, vol. 35, no. 5-6, pp. 503–510, 2008.
12. M. Reed, E. Gundlach, and T. Kana, “A coastal zone oil spill model: development and sensitivity studies,” Oil and Chemical Pollution, vol. 5, no. 6, pp. 411–449, 1989.
13. J. A. Galt, G. Y. Watabayashi, D. L. Payton, and J. C. Petersen, “Trajectory analysis for the Exxon Valdez: Hindcast study,” in International Oil Spill Conference (IOSC '91), pp. 629–634, Washington, DC, USA, 1991.
14. E. Howlett, K. Jayko, and M. L. Spaulding, “Interfacing real-time information with OILMAP,” in Proceedings of the 16th Arctic and Marine Oil Spill Program Technical Seminar, pp. 517–527, 1993.
15. V. S. Kolluru, M. L. Spaulding, and E. L. Anderson, “A three-dimensional oil dispersion model using a particle based approach,” in Proceedings of the 17th Arctic and Marine Oil Spill Program Technical Seminar, pp. 867–894, 1994.
16. Y. S. Cho and P. L. F. Liu, “Crest-length effects in nearshore tsunami run-up around islands,” Journal of Geophysical Research C, vol. 104, no. 4, pp. 7907–7913, 1999.
17. Y. S. Cho and H. J. Kim, “Sea level change due to typhoons in Korean peninsula,” Journal of Coastal Research, vol. 22, no. 3, pp. 678–682, 2006.
18. B. H. Choi, H. M. Eum, and S. B. Woo, “Modeling of coupled tide-wave-surge process in the Yellow Sea,” Ocean Engineering, vol. 30, no. 6, pp. 739–759, 2003.
19. P. L. Brezonik, Chemical Kinetics and Process Dynamics in Aquatic Systems, Lewis Publishers, 1994.
20. S. V. Patankar, “Numerical Heat Transfer and Fluid Flow,” Hemisphere Publishing Coporation, 1980.
21. D. L. Stefanovic and H. G. Stefan, “Accurate two-dimensional simulation of advective-diffusive-reactive transport,” Journal of Hydraulic Engineering, vol. 127, no. 9, pp. 728–737, 2001.
22. T. Komatsu, K. Ohgushi, and K. Asai, “Refined numerical scheme for advective transport in diffusion simulation,” Journal of Hydraulic Engineering, vol. 123, no. 1, pp. 41–50, 1997.
23. J. M. Hobson, N. Wood, and P. J. Mason, “A new finite-difference diffusion scheme,” Journal of Computational Physics, vol. 125, no. 1, pp. 16–25, 1996.