Spill , Transport and Fate Model ( STFM ) : Development and Validation

The Spill, Transport and Fate Model (STFM) is a new computational tool for modeling oil spills in Brazilian waters, developed at the Institute of Astronomy, Geophysics and Atmospheric Sciences at the University of São Paulo. The STFM was designed to assist in the analysis of environmental impact and water quality, for academic purposes and for environmental licensing studies. This work presents the formulation of the model, performance evaluations and validation of the results. These are necessary phases for the model to be made available to the public.


INTRODUCTION
The Spill, Transport and Fate Model (STFM) is a new trajectory and weathering model for handling oil spills. It is under development at the Institute of Astronomy, Geophysics and Atmospheric Sciences, University of Sao Paulo (IAG/USP) for application in marine studies and environmental impact assessment (Zacharias et al., 2018). In Brazil, the STFM is registered at the Instituto Nacional da Propriedade Industrial (INPI − Brazil) under contract BR512021002447-8. It is a three-dimensional model based on Lagrangian elements (LE) that uses atmospheric and coastline data provided by the Weather Research and Forecasting (WRF) Model, and also uses current, temperature, salinity, and depth data provided by the Hybrid Coordinate Ocean Model (HYCOM).
The STFM was successfully tested in the mysterious oil spill that occurred off the northeast coast of Brazil in 2019. The model was used in association with scenario trees, allowing the first estimate of the original volume of spilled oil and computationally confirming that the oil trajectory was subsurface (Zacharias et al., 2021a;2021b).
There are several oil spill models in the literature which are used to assess accidents at various stages of both oil production and transport. These models may range from simple parametric calculations to advanced three-dimensional numerical models. They could be coupled with meteorological, hydrodynamic, and wave models in order to forecast the transport and weathering of the oil in high resolution. However, the oil spill model will only be as accurate as the quality of the atmospheric and ocean models in which it is embedded (Keramea et al., 2021).
The development and testing of these models constitute an important area of research even eleven years after the Deepwater Horizon events (French-McCay et al., 2021a;2021b). Even in the case of the Brazilian-2019 "Mysterious Oil Spill", in which the oil drift was subsurface, the models still needed further development to provide more definitive answers (Lessa et al., 2021;Zacharias et al., 2021a;2021b).
In Brazil, some progress has also been made in this area with the development of ECOS (Easy Coupling Oil System) model, registered at the Instituto Nacional da Propriedade Industrial (INPI − Brazil) under contract: BR 51 2013 00013 (Stringari et al., 2013;2014;Marques et al., 2017;Lopes et al., 2019) and Coupled Model for Oil Spill Prediction (CMOP), developed by the National Institute for Space Research (INPE) in São José dos Campos, Brazil (Tessarolo and Innocentini, 2016;Tessarolo et al., 2018;Barreto et al., 2021). Along with the STFM (Zacharias et al., 2018;2021b), these are the three most-developed models currently in Brazil.
The development of Brazilian models that can be integrated with operational hydrodynamic platforms such as REMO (Oceanographic Modeling and Observation Network, https://www.rederemo.org/, last visit on: 08/28/2021) and BSO (Brazilian Sea Observatory, https://brazilianseaobservatory.org, last visit on: 08/28/2021) is an important and necessary step for the development of a marine disaster prevention system (Franz et al., 2021).
The behavior of an oil spill in the marine environment is initially governed by ocean and atmospheric conditions, which determine the transport and spread of the slick and by the physicochemical properties of the oil that determine the weathering (Fay, 1971). The most common model sets involve spreading, advection, diffusion, evaporation, emulsification and dispersion, and they do not consider processes, such as oil dissolution, photo-oxidation, biodegradation and vertical mixing ( Figure 1). Currently, uncertainty estimates and timely response to oil spills are lacking in the new generation of oil spill models. Further improvements should emphasize better parametrizations of oil dissolution, biodegradation, entrainment, and prediction of oil particle-size distribution (Keramea et al., 2021;Mohammadiun et al., 2021). The objective of this work is to present the latest updates of the STFM and its performance evaluation, implemented after the simulations of the oil spill that occurred in 2019 near the northeastern Brazilian coast. This will enable the STFM to be made available in the future, allowing users to have prior access to the STFM performance results compared to the reference models. This is a fundamental step for confirmation of the method quality developed by the STFM.

MATERIAL AND METHODS
The STFM was designed to connect the concepts of particle and box-model, which means that each simulated Lagrangian Element (LE) is also a small box-model with its own mass balance, calculated independently for each LE. This design was particularly useful in simulating the 2019 oil spill off the Brazilian coast, to represent hundreds of small slicks that spread along the shore (Zacharias et al., 2021b).
The position of each Lagrangian Element (LE) is given by the vector sum of the internal and external velocities of the oil slick (Equation 1) (Lynch et al., 2015).
Where P (x, y, z), VEf, Vspd, and Vdif define the instantaneous position of the particle in space, resulting speed of the winds and currents, gravitational spreading speed, and the turbulent diffusion velocity parameterized by a random motion equation, respectively.
Advection is the tridimensional transport of the LE caused by winds, currents, and vertical transport due to breaking waves and buoyancy. The STFM uses the friction velocity from WRF to calculate the atmospheric drag coefficient and solve the advection as Equation 2 (Lynch et al., 2015).
Where ksea and kbuoy are the transport coefficient constants (=1.0) for the ocean current and buoyancy, respectively; kwind is the coefficient of transport by wind, Vsea, Vwind and Vbuoy are the speed vectors for currents, winds and buoyancy, respectively.

The buoyancy speed is based on Stokes' law (Equation 2a
). This velocity is particularly important for undersurface spills or leaks , Equation 2b.
Where g is the gravitational acceleration, ρo is the oil density, ρw is the water density, vw is the water viscosity, di is the droplets diameter randomly, z is the LE depth (0 for surface and negative for undersurface depths), e 10z is an inertial wind speed damping function that quickly decays to zero as the LE sinks, u* and U10 are the friction velocity and wind speed module at a height of 10 m.
Spreading is the first effect solved by the STFM when the oil is released. It consists of a horizontal evolution of the slick, resulting from the action of gravity and buoyancy forces, facing the viscosity resistance (Fay, 1971;Dodge et al., 1983;Lynch et al., 2015).
The initial area (Equation 3a) evolves very quickly (a good number of minutes) compared to the oil slick's lifetime (Equation 3b). After a certain moment, turbulent diffusion becomes predominant and spreading reduces to zero, following the increase in viscosity. (3b) Where k1 and k2 are Fay's constants (0.5 and 0.75 respectively), A is the area of the slick, V is the volume of spilled oil, and M is the mass of the oil remaining on the slick. Turbulent diffusion is the mass transport inside the oil slick through random and chaotic motions. This effect is usually parameterized in different ways within the Lagrangian models to represent the subgrid scale effects of atmospheric and hydrodynamic models.
The STFM solves the Brownian motion of the LE introduced by Langevin's theory, using the Ornstein-Uhlenbeck process (Equation 4) resulting in less abrupt changes in the velocity and position of the Lagrangian particles (Gillespie, 1996;Lynch et al., 2015).
Where R and c are positive constants called the relaxation time and the diffusion constant, respectively; and G(t) is the ''Gaussian white noise,'' which may be defined as the dt→0 limit of the temporally uncorrelated normal random variable with mean 0 and variance 1/dt.
The dissolution process consists of the loss of soluble fractions from the oil slick to the water column. The oil components' dissolution rates are much lower than the other weathering effects. Analyzing the mass loss of the slick, dissolution is practically negligible (0-2% loss) 5 Spill, Transport and Fate Model (STFM) … Rev. Ambient. Água vol. 17 n. 1, e2789 -Taubaté 2022 when compared to its evaporation (30-90% loss); this is the reason why many models presented in the literature do not include dissolution in their formulations (Keramea et al., 2021). However, some specific components, such as benzene, may be of special interest because of their higher dissolution rate and their levels of toxicity to marine life.
Hydrocarbons that ultimately dissolve in water are mostly aromatic and dissolved in the sub-surface interaction between oil droplets and water. The STFM discretizes dissolution by elements of interest (Equation 5) to assess the residual impact on water quality resulting from dissolved fractions, which cannot be removed even after oil clean-up (Huang and Monastero 1982;Cohen et al., 1980;Lynch et al., 2015).
Where Sd is the dissolution rate (g/h/m 2 ) of the "i" component, Kd is the mass transfer speed due the dissolution (0.0002 m min −1 ), fs is the fraction of the surface covered by crude oil, T i is the mass content of component "i" in the spilled oil (dimensionless), E i is the increase of dissolution of component "i" in comparison to the average dissolution of the spilled oil (dimensionless), S0 is the oil initial solubility (30 g/m 3 ), α is the decay constant (min -1 ), t is the time after the spill in minutes.
The values of E and S0 are not easily found in the literature, so some reference values have been provided in Table 1 and Table 2.   (Lynch et al., 2015).
The dispersion process of the oil slick consists of the gradual entrance of oil fractions (droplets) in the vertical column of water, and, due to weathering of lighter fractions or the collision/coalescence process with suspended sediment, or the use of dispersants, these droplets no longer return to the surface.
In probabilistic approaches or preliminary simulations, STFM uses a simplified model that integrates vertical dispersion, sedimentation and sinking in the water column Zacharias et al., 2021a). In the usual simulations, STFM uses the "classic" formulation given by Equation 6 (Delvigne and Sweeney, 1988 Where dQ is the mass rate entangled in the water column by the diameter class D (kg/m 2 /s), E is the energy of the waves dissipated per unit area for a single breaking wave (J/m 2 ), Fwc is the fraction of oceanic coverage with waves, Scov is the fraction of oceanic coverage with oil and C0 is the empirical constant of the oil.
The STFM uses a very simple droplet size distribution with 20 classes based on Delvigne and Sweeney (1988) for heavy crude oils, and on Cekirge et al. (1997) for orimulsion or fuel residual oils in plumes.
Evaporation is the main process in oil weathering, affecting mainly the lighter components, and also responsible for the natural removal of oil spilled on the ocean surface. It is calculated by the STFM using the empirical set of equations obtained by Fingas (Equation 7), based on the atmospheric exposure time of the oil slick and the oil temperature (Fingas, 2016).
Where Ka is the evaporation coefficient A of the oil, Kb is the evaporation coefficient B of the oil, T is the oil temperature (°C) and t is the time in minutes.
The formation of emulsions (water immersed in the oil) results in a significant increase in the oil slick volume, changing its dimensional characteristics, such as area and thickness. Weathering is very important in all oil spill models; however, it is also one of the most difficult effects to be properly predicted, with 30 ~ 50% overestimation (Keramea et al., 2021). The basic model is given by Equation 8 Fingas, 2016).
STFM has been compared with the results obtained by two reference models: ADIOS2 (Lehr et al., 2002) and GNOME (Zelenke et al., 2012), and with other algorithms previously tested Wang et al. (2005), Stringari et al. (2013;2014) and Zadeh and Hejazi (2012), using the methodology previously described and analyzed on the Brazilian coast Zacharias and Fornaro, 2020) and using a model performance evaluation (Chang and Hanna, 2004). The performance evaluation of the model consists of the observation of the distance between the model's obtained results and the expected values, so that the less the distance and systematic errors, the better the performance.
The environmental data (coastline, meteorological and hydrodynamic fields) used in this study have been provided by GNOME Online Oceanographic Data Server (GOODS), which serves as a user-friendly interface to the NOAA global database, allowing the extraction of environmental information and its exportation in GNOME input file formats (https://gnome.orr.noaa.gov/goodslast visit on Aug/29/2021).
The Global Forecast System (GFS) is a complete global-scale weather forecasting system developed by the National Center for Environmental Prediction (NCEP), providing data at resolutions ranging from 0.25º to 2.5º. The GOODS platform has access to 0.5º and 1º spatial 7 Spill, Transport and Fate Model (STFM) … Rev. Ambient. Água vol. 17 n. 1, e2789 -Taubaté 2022 resolution data; therefore, due to the need for compatibility with GNOME, data with 0.5º spatial resolution, provided every 6 hours, have been used.
STFM has used the WRF (Weather Research & Forecasting Model) atmospheric data. The WRF simulation covers the area 20° S, 50° W to 10° N, 20° W, with a 1-h time interval and 0.15 o horizontal resolution (Skamarock et al., 2019). The WRF was driven by initial and boundary conditions from the National Center for Environmental Prediction (NCEP) database from the Global Forecast System (GFS) historical archive (dataset name: ds084.1) at a 6-h interval and 0.25 o horizontal resolution (NCEP et al., 2015).
GOODS extracts surface current data directly from Global HYCOM, which operationally simulates ocean conditions with 1/12º horizontal resolution (ie, one grid element every 7 km, approximately, for the domain of this study), being able to solve the oceanic circulation in turbulent scale (eddy), presenting the vortices and meanders of the main currents in the region (Chassignet et al., 2007).
The Hybrid Coordinate Ocean Model (HYCOM) is a computational model of general ocean circulation that numerically solves primitive hydrodynamic equations. The main evolutionary characteristic of HYCOM (in comparison with its predecessors) is the use of hybrid vertical coordinates, establishing vertical levels by isobars (levels of the same pressure), isopycnals (same density) or by sigma coordinates (height with respect to the background level) thus obtaining better results in coastal (shallow waters) and oceanic simulations (Chassignet et al., 2007).
The initial set of simulations has evaluated the movement of the oil drift models, including advection, spreading and diffusion. Three spill locations have been used in different areas of Santos Basin (Zacharias and Fornaro, 2020).
The set of simulations has considered 3 spill points, 5 different models and 6 initial dates, totaling 90 numerical simulations (Table 3). The model performance evaluation has been carried out using the methodology given by Chang and Hanna (2004). The weathering and physicochemical analysis of the models have been based on the comparison with reference model ADIOS2 (Lehr et al., 2002), using crude oils from the Arabian series (Extra Light, Light, Medium and Heavy), obtained from ADIOS2 (Table 4).  (Lehr et al., 2002) and Oil Spill Science and Technology (Fingas, 2016).

RESULTS AND DISCUSSION
The STFM model performance evaluation has presented satisfactory results in oil slick full-movement simulations (Table 5), comparing those with all reference models, using the previously established criteria (Zacharias et al., 2018).
The variations in the Model Performance Evaluation results are within the established limits, which means that all tested models have presented similar and concordant predictions of the oil slick's movement. Weathering results have been consistent among all tested models (Table 6). Evaporation and emulsification are always the most important factors in weathering. In general, the models use equations that are very similar to each other to estimate these two effects, so the results of the model performance evaluation are very similar. The pseudo-components method is the natural evolution for all these models; however, it faces the lack of adequate information about the composition of each type of oil.
Density and viscosity variations resulting from changes in the oil composition are affected by evaporation and emulsification. Again, without good oil speciation, parameterizations are used to represent these effects, and the tendency is for uncertainties to accumulate and estimates 9 Spill, Transport and Fate Model (STFM) … Rev. Ambient. Água vol. 17 n. 1, e2789 -Taubaté 2022 of these variables to be less accurate (Table 6). Oil dispersion in the water column is parameterized by the wind, resulting in a very easy method to be implemented. However, the formation of waves and droplets tends to be more complex. This part of the weathering will undergo substantial improvement with the assimilation of wave models as input data (Table 6).
In general, the STFM has presented excellent results of model performance evaluation when compared with other models in the same category, meeting the necessary criteria for model validation.

CONCLUSION
The STFM code has been tested against the other developed codes (Wang, Stringari and Zadeh) and reference models (GNOME and ADIOS2), with validation criteria proposed by Chang and Hanna. The model works according to the proposed premises and meets the objective of simulating events of oil spill on the ocean surface, successfully determining the drift trajectory and weathering of the spilled oil.
The physicochemical equations used in the STFM can be easily coupled with hydrodynamic models (such as MOHID or SISBAHIA), including all proposed weathering (evaporation, composition of residual oil fractions, emulsification, alteration of density and viscosity, dissolution and toxicity).
The STFM has demonstrated the ability to handle different types of oil (from extra-light to heavy oils), covering the whole spectrum of commonly used oils and spillage possibilities. The presented model has great potential for evolution and use in the existing coastal, estuarine and water quality models.
The results have shown that the STFM is fully calibrated and validated for simulations on the Brazilian coast. Evidently, future simulations with STFM will have the same quality as the data provided by HYCOM and WRF, so for time-response applications, it is really important that HYCOM and WRF are properly calibrated as well.