Atmospheric Contamination of Coastal Cities by the Exhaust Emissions of Docked Marine Vessels: The Case of Tromsø

: Docked ships are a source of contamination for the city while they keep their engine working. Plume emissions from large boats can carry a number of pollutants to nearby cities causing a detrimental effect on the life quality and health of local citizens and ecosystems. A computational ﬂuid dynamics model of the harbour area of Tromsø has been built in order to model the deposition of CO 2 gas emitted by docked vessels within the city. The ground level distribution of the emitted gas has been obtained and the inﬂuence of the wind speed and direction, vessel chimney height, ambient temperature and exhaust gas temperature have been studied. The deposition range is found to be the largest when the wind speed is low. At high wind speeds, the deposition of pollutants along the wind direction is enhanced and spots of high pollutant concentration can be created. The simulation model is intended for the detailed study of the contamination in cities near the coast or an industrial pollutant source of any type of gas pollutant and can easily be extended for the study of particulate matter.


Introduction
Plume emissions from industrial activities and large boats can carry a number of pollutants to nearby cities causing a detrimental effect on the life quality and health of local citizens and ecosystems [1,2]. The main pollutants are waste products of combustion processes, mainly carbon monoxide (CO), nitrogen and sulfur oxides (NO x , SO x ) in gas phase, as well as solid particulate matter and black carbon. They are suspended in the released exhaust gases. Carbon dioxide (CO 2 ) is the main by-product in exhaust plumes. While it does not have a direct health effect, its influence on climate warming is well known.
Docked ships are a source of contamination for a city while they keep their engine working [3][4][5][6]. Emission from manoeuvring ships can also be substantially high. However, docked vessels become an important source of contamination when they must keep their engines working to maintain their vital functions while in the harbour. The deposited fraction of particulates and gas traces at street level is an indicator of the air pollution in the city. The exhaust gas and the suspended particulates are transported by the wind over the habited areas [7,8]. The interplay between the wind and the local orography will determine the transport distance and the fraction reaching the ground level [9,10].
The pollutants can be divided into two types: Gas molecules and larger particulate matter [11]. Each has different transport mechanisms. While gas pollutants are transported as part of the main fluid, solid particles are subject to drag by the fluid and they tend to suspend by gravity. Both types of particles are emitted from an exhaust at a relatively large concentration and they mix and dilute thanks to buoyancy effects and lateral diffusion. This last effect is especially relevant in increasing the diffusion range of the gas phase pollutants. On the other hand, particulate matter is mainly driven by buoyancy and the drag forces. This has an effect on their distribution close to the emission point. An accurate assessment of the contaminant distributions can be achieved using real scale atmospheric simulations including the local terrain orography and building distribution. A real scale model should consider as input parameters the atmospheric condition (temperature, wind and pressure), the topography of the studied area and the position and size of the main buildings. To study the contamination coming from docked vessels, the height of the exhaust and the temperature of the plume must also be considered.
In this work, we use computational fluid dynamics simulations of flowing air on a topographic model of the sea side in the city of Tromsø. The emission from docked vessels and the interplay between the general wind flow and the exhaust gas have been modelled in detail. We expect that this result will help to asses the effect of docked vessels in the contamination of the city of Tromsø and to introduce mitigation measures. Special attention has been given to the dispersion of a trace gas, CO 2 in this study, that is representative for the dispersion of pollutants emitted by the vessel and their distribution at street level. For that purpose, the emission from a vessel situated in the Tromsø harbour has been considered. The parameters that influence the final distribution of the emitted pollutants have been varied. Special attention has been given to the influence of the wind speed and direction, the ambient and outflowing gas temperatures and the chimney height. In Section 2, the used computational technique are described. The main results are given in Section 3 and the most important results are discussed in Section 4. Finally, the main conclusions are summarized in Section 5.

Materials and Methods
The study of the pollutant distribution is based on a computational fluid dynamics (CFD) simulation of the exhaust gas distribution over and within the city of Tromsø. The simulation domain is constructed using a topographic model of the city harbour area. The coordinates of the topography were obtained from satellite data [12], which allows to download the longitude, latitude and altitude along selected paths. The longitude and the latitude from the raw data have to be converted to Cartesian points (x, y, z) and written in comma separated value (CSV) format. The position data is imported into a computer aided design software [13]. The positioning and size of the largest buildings in the harbour area were defined in the same way as the topography of the terrain and incorporated into the ground model. The structures of concern were storage facilities of the port, a high school, an office building, part of the hospital and a prison. Finally, three cylinders in the harbour area, representing the chimneys of a docked vessel, were added. The collection of points was then converted into a 2D model of the terrain within Solidworks. The resulting ground topographic model has a diameter of 800 m, shown in Figure 1, and describes the main buildings near the harbour area together with the main features of the terrain. The 3D geometry is completed by extruding a cylinder on top of the terrain model with a total height above sea level of 300 m. The height of the ground model ranges between 0 m at the sea and 34 m at the highest point. The typical height of the buildings ranges between 5 m and 10 m. An air column of up to 300 m was added for the fluid simulation. The resulting simulation domain has a cylindrical shape truncated below by the non-flat terrain.
Several chimney heights were considered: 5 m, 10 m and 20 m. For the temperatures of the emitted gas at the exit of the chimney values between 10 C and 200 C were used. The exit speed of the exhaust gas, 25 m/s, and the fraction of emitted CO 2 , 40%, were set at typical values for large vessels. The environment temperature varied between 0 C and 20 C.
The 3D model was meshed using ANSYS Workbench [14] to a 1.3 10 6 tetrahedral cells mesh. The boundaries were meshed using a triangular mesh and they were divided into four regions: Top boundary, side boundaries, ground and exhaust of the chimney. The exhaust region will be used to define the plume inlet conditions. The side boundaries will accommodate the inlet/outlet conditions of the wind and the atmospheric boundary layer conditions. The ground will be distinctive in that it will require no-slip and wall conditions for the simulated quantities. Finally, the top boundary will be used to ensure a smooth vertical distribution of temperature and pressure. The mesh was refined around the main features of the terrain. The size of the mesh elements and the degree of refinement close to the main features are illustrated in Figure 2.

Computational Fluid Dynamics
In CFD the air column is described with the Navier-Stokes equation with the gravitation g as an external force The magnitudes solved for are the fluid velocity u and the atmospheric pressure P. The fluid has a kinematic viscosity ν and its density is ρ. The gases are compressible and ρ is modified by the combined effect of the temperature T and pressure P. They are linked by the equation state of ideal gases P = nRT.
R is the molar gas constant 8.315 J/(K mol) and n is the molar density. The total fluid density ρ and the fluid velocity u are linked by the mass transport equation The atmosphere of Earth is composed of nitrogen (78%) and oxygen (21%) mainly. Several other gases appear at a lower concentration (noble gases argon, neon, helium and krypton, carbon dioxide CO 2 , methane NH 4 and water vapour). In this work, we are interested in the transport and diffusion of trace gases. For this purpose, we separated the gas phase into two components. The main component describes dry air with a mass density of 1.2 kg/m 3 at sea level and a trace component expelled from the vessel exhaust. For the trace component CO 2 was chosen throughout the work. Within the ideal gas approximation, the behaviour of the trace gas is representative of other gas phase pollutants also.
The wind velocities were defined as the wind speed and direction in the upper part of the simulation domain. The speed values considered for this work ranged between 1 m/s and 15 m/s. Four directions, north was (0, 1, 0), were considered starting from u = ( 1, 1, 0) and varying the orientation of the direction vector by 90 . The high Reynolds number Re = uL/ν 1000 corresponded to the turbulent regime. The κ Reynolds Averaged Simulation (RAS) turbulence model [15] was used to describe the turbulence.

Atmospheric Boundary Layer
The simulation of an atmospheric air column required setting up carefully the boundary conditions of the circulating air. The wind flow along the computational domain was stabilized using the atmospheric boundary layer model developed by Richards and Norris [16] with the generalized formulation used by Hargreaves et al. [17] and Yang et al. [18]. The velocity is fixed and horizontal in the atmospheric boundaries with a vertical profile κ = 0.41 is the von Karman constant, z g (x, y) is the topographic ground height, the value of which at each point depends on the topographic height, and z 0 is a surface roughness coefficient that was set up to 5 m to match the buildings' height. The friction velocity U is a reference value calculated as U = κu re f / ln((z re f + z 0 )/z 0 ). The velocity u re f was chosen equal to the wind speed in the upper limit of the simulation domain. For the temperature of the air at sea level (z 0 m) T 0 = 0 C, 10 C and 20 C were chosen. The temperature in the atmospheric lateral boundaries was fixed to a linear profile T = T 0 0.01z. The turbulence model parameters were set up as C µ = 0.09 being the turbulent viscosity coefficient.

Numerical Solutions
The simulations were run using the OpenFOAM [19] CFD software based in the field operation and manipulation (FOAM) C++ class library for continuum mechanics. OpenFOAM uses the finite volume numerical method to integrate the Navier-Stokes equation.
The simulated gas mixture is composed by atmospheric air and a trace component that in this simulation was identified with CO 2 . The trace gas was expelled by the docked vessel through its chimney. In the plume outlet, 40% of the outgoing gas mass corresponded to CO 2 . The CO 2 fraction at each cell was updated every time step using the mass conservation equation (Equation (3)). Buoyancy effects were included and the thermodynamical model of the fluid takes included the enthalpy h f and the specific heat capacity c p . The transport properties were described through the dynamic viscosity µ and the Prandtl number Pr. Their values for each component of the mixture are given in Table 1. The boundary conditions set for the simulations are listed in Table 2.  (1)-(3) are solved using a combination of the PISO (Pressure-Implicit with Splitting of Operators) algorithm [20] and the SIMPLE (Semi-Implicit Method for Pressure Linked Equations) algorithm [21]. This choice allows to increase the size of the time step to speed up the simulation. The transient simulation was run for 1000 s simulation time in parallel in 32 cores. The simulation time was chosen to ensure that the wind flow and the CO 2 distribution reached a stationary state and it was short enough to disregard disturbances on the breeze circulation coming from tides and insolation.
The study addressed the effect of wind speeds and directions, atmospheric temperature, chimney height and exhaust gas temperature. Wind speed and directions cover the most usual atmospheric conditions in the city of Tromsø. On the other hand, chimney height and exhaust gas temperature describe the most important variables of the exhaust emission. The values chosen for the parameter study are given in Table 3. Table 2. Boundary conditions. All the boundary conditions used are detailed. The solver solves the quantity P ρgz instead of the total atmospheric pressure. Its boundary value is set up so that the gradient matches the velocity boundary condition.

Results
The simulation of the atmospheric boundary layer poses a challenge for CFD. The main problem lies in the inherently turbulent nature of the atmospheric flow and cloud formation, further complicated by the orography of the terrain and the buildings [22,23]. A realistic wind flow requires taking into account many details like changes in transport properties with temperature and pressure, as well as the boundary conditions in contact with the terrain or in the top of the simulation domain. The thermodynamical properties of the gas mixture composing the air has also to account for the heat transfer properties with diverse heat sources, the solar insolation or the content of water vapour.
When the simulation is performed for stable atmospheric conditions and it is constrained to the lower part of the atmospheric boundary layer, the surface layer, a stable wind flow can be obtained using a k turbulence model [24] and by setting the boundary conditions proposed by Richard and Hoxey [16] for the velocity, the pressure as well as the turbulent kinetic energy k and its dissipation rate .
Between the lower and the upper boundary of the simulation domain the thermodynamical quantities vary by only 3.5% for the pressure, 1% for the temperature and 2.5% for the molar density n of an ideal gas. The narrow thermodynamical conditions allow using a constant value of the transport and thermodynamic properties.
In contrast, the simplifying assumption of incompressibility cannot be applied to atmospheric gas. The equation of state for the perfect gas was used to connect the pressure to the thermodynamical properties of the fluid (temperature, enthalpy and specific heat capacity) as well as the mass density. The relevant parameters are given in Table 1. To ensure a stable run, an appropriate choice of boundary conditions and initial conditions is crucial. To constrain U and the turbulence model parameters (k, ) in the lateral boundaries the ABL condition was used. The ABL values correspond to a stable dry atmosphere with a defined wind flow at heights well above the orographic obstacles. The temperature was constrained to be equal to T 0 at ground level and had a vertical gradient of 0.01 K/m. At the top boundary, all quantities had zero gradient boundary conditions, except for U which had a fixed value corresponding to the free wind flow. On the bottom boundary, no slip condition was used for U and wall functions for the turbulence parameters. The rest of the magnitudes had zero gradient conditions.

Simulation Run
In a first step, the simulation is run for 100 s without exhaust gas to let the wind flow stabilize; the resulting velocity distribution is shown in Figure 3a). In the upper part of the simulation domain the velocity field was parallel to the general wind flow. Closer to the ground, the wind flow was disturbed by the orography and the city buildings. In a second step, the simulation was run with the outflow from the exhaust included for at least 500 s so that the new air flow reached a stationary regime. The atmospheric wind flow was then disturbed by the ascending air + CO 2 plume, see Figure 3b). The dynamic pressure, P ρgz, shows the low pressure zone created by the plume on the general air flow. CO 2 gas was transported upwards by the fast outgoing plume and subsequently drifted following the general wind flow. At the same time, sidewards from the plume an outdiffusion of the pollutant occurred. Part of it reached the ground level where an excess of CO 2 could be built. The pollutant excess is defined as a 10% excess over the background level of 400 ppm. This happened mainly right below the plume except for the width of The concentration of CO 2 in the city only increases when the wind blows inland, north-west. The wind blowing in the opposite direction, towards the sea or south-east, also deposit pollutant gas on the sea. In both directions, the gentle slope of the terrain near the sea increases the deposition rate both inland and on the sea. In line with this, the winds blowing parallel to the coast, north-east and south-west, yield a much lower deposition of the pollutant. The main reason for this is that the exhaust plume travels upward higher.

Wind Speed
The speed of the incoming flow is a determining factor in the range of dispersion of the pollutant. In general the largest concentration at ground level is always detected right below the plume. However, the sidewards distribution and the peak concentration value of the deposited gas molecules can vary remarkably. The widest distribution of excess pollutant is obtained for the lowest wind speed. While the highest concentration is reached downwind, the concentration sideways remains high in the whole simulation domain. In addition, there is also a noticiable deposition of pollutant upwind. For a medium value, 5 m/s, and the highest value, 15 m/s, the gas deposition remains detectable in narrower areas. It has a width of 120 m in the first case and 200 m in the second at the maximum deposition concentration distance from the chimneys.
The case with wind speed of 3 m/s is different from the rest in the fact that the deposition of the pollutant is lower by an order of magnitude. This anomaly is better understood by noting that the plume, see the side view in Figure 5, flows higher and therefore less pollutant reaches the ground. Larger wind speeds increase the drag of the gas molecules away from the chimney diminishing the time for lateral dispersion. This effect is the main reason why the concentration of pollutant is radically smaller when going from wind speeds of 1 m/s to 3 m/s. At larger wind speeds, the maximum height reached by the plume is shortened by the stronger pull exerted by the main wind flow. The interplay between these factors, the first being dominant for low wind speeds and the latter for wind speeds above 5 m/s, explains the observed trend.  ( 1, 1, 0). The chimney height is 10 m, the ambient temperature is 10 C and the exhaust gas temperature 110 C in all cases. The terrain height is represented in a grayscale. Lighter colors correspond to lower heights and the sea has the lightest color. The chimneys are three white spots near the coastline. The upper panels show the CO 2 concentration isosurfaces (0.1%, 1000 ppm) representing the plume distribution within the simulation domain.

Temperature of the Plume
Another important factor determining the final distribution of CO 2 in the ground is the temperature of the outflowing gas at the chimney exit. Higher outflow temperature enhances the upwards buoyancy of the exhaust gas due to its lower density. In contrast, at lower temperatures the exhaust gas will get dragged downwards. When the temperature of the exhaust gas is equal to the ambient temperature T c = 10 C, see Figure 6, the amount of deposited pollutant is strongly enhanced. When the exhaust gas temperature is T c = 210 C the concentration decreases by 2 orders of magnitude and the side range of the deposited gas decreases to less than half.
!""#$ !"#$%&"$' ()*%&"$' ( " *+#,#,$- Figure 6. Exhaust gas Temperature. Concentration of CO 2 at ground level for wind speed u = 5 m/s and north-west ( 1, 1, 0) direction. The chimney height is 10 m and the ambient temperature 10 C in all cases. The terrain height is represented in a grayscale. Lighter colors correspond to lower heights and the sea has the lightest color. The chimneys are three white spots near the coastline. The upper panels, show the concentration of CO 2 within the plume. The 0.1%, 1000 ppm, surface for CO 2 was plotted in all cases.

Chimney Height
The chimney height (see Figure 7) has a negligible effect on the deposited gas within the simulated area. Moreover, the ambient temperature, shown in Figure 8, has a small effect on the deposited pollutant. In all cases, the region where a significant deposition of CO 2 is observed remains largely unchanged and the concentration of the deposited pollutant is also similar.

Discussion
The wind blowing directly inland, north-west ( 1, 1, 0), is found to be the most relevant to understand the transport of gas phase pollutants into the city. In all cases, the deposition of CO 2 is the largest right below the center of the plume. When the wind speed is very low the concentration downwind stays relatively constant within the simulated domain, as can be seen in Figure 9. The drag of the main wind flow is low and the pollutant disperses mostly by diffusion. It must be noted that the concentration upwind has a sharp maximum near the chimney and decays fast at larger distances as shown in Figure 5. For slightly larger wind speeds, the deposited gas increases up to a distance of around 200 m. At longer distances it remains roughly constant until the boundary of the simulation domain. When the wind speed is 5 m/s or higher, the pollutant still has a peak concentration, but at longer distances from the chimney, the concentration of pollutant decreases, especially for the largest wind flow speed. The peak CO 2 concentration is the lowest at moderate speeds, 3 m/s, and increases steadily as the wind speed increases.

Conclusions
An accurate CFD model of the harbour area of Tromsø was built in order to model the deposition of CO 2 gas emitted by docked vessels within the city. It includes buoyancy and fluid compressibility effects as well as the thermodynamics of perfect gases. Stabilization of the general wind flow was achieved by forcing the atmospheric boundary layer conditions in the lateral boundaries of the domain. The diffusion of CO 2 out of the plume into the surrounding air is also modelled and it is found to have noticeable effects, especially at low wind speeds. The influence of the wind speed and direction, vessel chimney height, ambient temperature and exhaust gas temperature were studied. The wind direction, as expected, determined the direction where the pollutant gas will be deposited. More interesting is the trend found when varying the wind speed blowing inland. The deposition range was the largest when the wind speed was very low. In this case the diffusion effects dominate over the transport effects and the pollutant gas has a wide distribution range extending several hundred meters downwind and in a lateral direction. At moderate speeds, the plume is transported far from the city area before the deposition happens. Interestingly, at even larger speeds, the height of the plume is shortened increasing the deposition below the plume. The lateral distribution is 100-200 m wide in most cases. The simulations show that the interplay between the main wind flow, the local orography and the city buildings can create pockets of high pollutant concentration.
This simulation model can be extended to study other types of pollutants or particulate matter. They can be of high interest for a detailed study of the contamination in cities near coasts or an industrial pollutant source. A next step to increase the accuracy of this model would be validation of the results using field data.