A Study on Plume Dispersion Characteristics of Two Discrete Plume Stacks for Negative Temperature Gradient Conditions

The dispersion of air pollutants emitted from industries has been studied ever since the dawn of industrialisation. The present work focuses on investigating the effect of negative atmospheric temperature gradient and the plume stack orientation of two individual equal-height stacks on the vertical rise and dispersion of the plume. The study carried out upon three-stack layout configurations namely inline, 45° and non-inline, separated by an inter-stack distance of 12 times the exit chimney diameter (12 D) and 22 times the exit chimney diameter (22 D) in each case over the two temperature gradients of −0.2 K/100 m and −0.5 K/100 m. The turbulence is modelled using realisable k-ε model, a model used in the FLUENT flow solver. In the case of the inline configuration, the upwind plume shields its downwind counterpart, which in turn allows for higher plume rise at a given temperature gradient. The plume oscillates more in the case of inline than 45° and non-inline cases. Also, for a temperature gradient of −0.5 K/100 m, the plumes oscillate violently in the vertical direction, mainly because, with the initial rise of the plume, cold air from higher altitudes moves down and forms a layer of lower temperature closer to the ground. The present study is important to highlight the plume dispersion characteristics under negative temperature gradient conditions.


Introduction
Plume dispersion has been a significant field of study to understand and monitor the effects of air pollution [1,2]. Some of the major pollutants present in flue gases are carbon monoxide, sulphur dioxide, nitrogen monoxide and nitrogen dioxide. They are emitted after combustion from various industries into the atmosphere. Initially, industries were established in less populated areas, but with time, more people populated the surrounding areas of the industry. Hence, it is important to properly disperse the pollutants into the atmosphere and reduce the environmental impact of these emissions. Thus, multiple flue stacks are used to treat each pollutant separately. Consequently, after treating these pollutant gases, they are released through separate flue stacks. The pollutant dispersion is affected by various factors like flue stack height, plume exit velocity from the stack, atmospheric conditions, the distance between the flue stack exits, wind velocity and direction and stack layout. Usually, the temperature gradient is −0.6 K/100 m altitudinal increase in a large-scale normal atmosphere. But, in the region close to the ground due to convection (i.e. less than 500 m), the temperature gradient variation is highly dependent on the several local meteorological conditions. Several authors working on CFD computations take temperature gradient (α = dT/dz) as the stability criteria while dealing with pollutant dispersion close to the ground [2][3][4][5]. The atmospheric conditions are therefore categorised based on temperature gradient with respect to the altitudeviz. stable (α > 0), neutral (α = 0) and unstable (α < 0).
Most of the pollutants are emitted into the atmospheric boundary layer (ABL), of the atmosphere, and the temperature gradients and wind velocity direction adversely affect the dispersion of the contaminant. The atmospheric boundary layer (ABL) is the lowest layer of the atmosphere, which typically extends from 0.5 to 1 km in altitude in diurnal hours. The boundary layer during the night is usually shallow owing to the absence of convection. The plume dispersion into the ABL is affected by both the number of stacks and the effect of layouts. Mokhtarzadeh-Dehghan et al. (2006) were the first to study the effects of the layout on plume dispersion using multi-stack outlets of cooling towers. The stacks were placed in tandem and non-inline layouts [6]. Also, works carried out by Contini et al. (2004) have played a major role in studying the dispersion of plumes from multiple individual flue stacks placed in several orientations and various inter-stack distances with respect to the wind [7].
Numerous analytical, experimental and numerical studies have been carried out in the past to study the plume dispersion and the effect of varying the factors mentioned above upon the plume dispersion. Initially, the plume dispersion was analytically modelled using Gaussian distribution models. Briggs (1965) used dimensional analysis to develop equations to compute plume rise under various wind speeds [8]. Using the virtual stack configuration, Anfossi et al. (1978) developed a model to study the maximum plume height. Their results obtained agreed with that of the Briggs model for stacks of equal height [9].
The various aspects of dispersion of plumes from a single flue stack were also studied by several researchers [3,[10][11][12][13][14]. The work carried out by Wee and Park investigated the effect of the crosswind velocities upon the plume rise and dispersion [3]. Simulations were carried out to study the dispersion of pollutants expelled from a single flue stack under negative gradient conditions of −0.5 K/100 m, −1 K/100 m and −1.2 K/100 m. Blocken et al. (2007) compared the ASHRAE and CFD RANS models for the numerical evaluation of plume dispersion [10]. They reported that using CFD with the standard k-ε model in combination with the standard wall functions provides an acceptable prediction of concentration of plume downstream. Considerable oscillations of the plume were observed in the vertical direction in all three cases. Onbasioglu (2001) observed that a jet with higher velocities affected the reverse flow and the entrainment of mass, whereas higher exit stack temperature led to a rapid decrease in concentration inside the stack [11]. Contini et al. (2009) studied the dispersion of plume using a small-scale wind tunnel model [12]. They observed that at a Reynolds number of 2196, the exit velocity profile in the cylindrical stack is distorted relative to the full-scale model, and the overestimation of the plume rise leads to the underestimation of the ground level concentration. Kozarev et al. (2014) investigated the plume rise in calm and low-velocity conditions by using a full-scale CFD model. The influence of the emission and stack parameters, the ground level air, the ground level air temperature, the temperature gradient of the atmosphere and the surface wind velocity ranging from 0 to 1 m/s was studied [13].
Multi-stack plume dispersion study predicted the formation of a pair of counter-rotating vortices and also the mixing of plumes emanating from a cooling tower using a tandem and side-by-side arrangement of plume exits using experiments in a wind tunnel. The same was validated using CFD [4]. They concluded that tandem source arrangement showed early mixing of the plumes and thereby displayed a higher plume rise. Another work by Mokhtarzadeh-Dehghan et al. (2006) showed that the penetration of the cross wind between the plumes played a major role in the mixing of the plumes [6]. Velamati et al. (2015) studied the plume rise and dispersion under positive and neutral temperature gradient conditions and in various configurations of the stack compared with the wind direction. Also, they studied the effect that the orientation of flue stack exits had with respect to the incoming wind direction upon the pollutant dispersal behaviour [5]. The two stack configurations used are the inline and non-inline configurations, and it was found that in the inline configuration, the plume rise was higher in the inline case due to the shielding effect of one plume over the other. The vertical dispersion of the plume is higher in the inline case due to the faster inter mixing of the same, which allows it to be spread over a larger area. However, the inline configuration possesses a higher vertical momentum due to the good mixing of jets; therefore, the lateral dispersion of the inline case is lesser compared with its non-inline counterpart.
Experimental studies have been carried out using wind tunnels and towing tanks to study pollutant dispersion from multiple plume stacks. Contini and Robins (2004) conducted experiments in a towing tank to study the plume rise of plumes emanating from two separate flue stack exits and mixing in neutral cross flow and on dispersion [7]. Their study generated results that indicated the presence of largely asymmetric shapes generated by merging of two plumes having opposite vorticity. They proved that the mixing of the plumes becomes protracted as the angle between the incoming wind direction and that of the line joining the centre of the flue stack exit increased beyond 45°. The mixing of the plumes is least in the case of a non-inline stack layout, where a three-lobed structure is seen in the crosssection of the plume. It causes the centre of gravity of the combined plume to shift downwards.
Another work by Contini et al. (2006) experimentally studied the mixing of identical plumes in a turbulent boundary layer using a wind tunnel [15]. The study found that the internal turbulence of the plume dominated during the mixing process. Also, Huq and Stewart (1996) studied the effect of laminar and turbulent cross flows upon the buoyant plumes [16]. They investigated the effects of turbulent cross flows on the development of plumes using detailed analysis of flow visualisation, hot wire anemometry etc. Huang et al. (2014) studied the impact of shape and height of the upstream roof on airflow and pollutant dispersion in an urban street canyon [17]. They conducted a 2D CFD study using FLUENT code upon various shapes of buildings and for three different aspect ratios of the buildings. They also validated the same using wind tunnel results. Macdonald et al. (2002) experimentally studied the behaviour of merging buoyant plumes using scaled water flume experiments [18]. They conducted experiments for various inter-stack distances and different exhaust velocity ratios for stack pairs. They observed the enhancement of plume rise when the stacks are in line with the flow separation, and little or no improvement in case of the stacks placed perpendicular to the flow direction. They concluded that this might be due to the 'momentum shielding' of the upstream plume. Yasin (2012) investigated the flow and dispersion of gaseous emissions from the vehicle exhaust in a street canyon under changes of the aspect ratio and the wind direction using numerical methods and also validated the same using a small-scale wind tunnel study [19].
To the best of the authors' knowledge, most of the plume dispersion studies have been done in either neutral or stable environmental conditions. Only Wee and Park (2009) have performed numerical simulations for single stack in negative temperature gradient conditions [3]. And plume dispersion studies over multiple stacks have only been performed by Contini et al. (2004) [7]. This study has only been carried out for neutral temperature gradient conditions. Under normal circumstances, the atmosphere has a negative temperature gradient of around −0.6 K/100 m. The negative temperature gradient conditions have a substantial impact on the spatial distribution of the ground-level concentration of pollutants, which attracts the interest of environmentalists. Many large industries are also

Problem Description
The present work studies the effect of the negative atmospheric temperature gradient on the plume rise and pollutant dispersion for three different configurations of the multiple individual flue gas exits existing from their respective flue stacks-inline, 45° and non-inline as shown in Fig. 1. Also, both the stacks are of equal height and separated by an inter-stack distance of 12 D and 22 D. The temperature gradient prevalent in the atmosphere with respect to the altitude has been denoted by α. The values of α considered in the present study are −0.2 and −0.5 K/100 m. The exit temperature of the plume is 397.15 K, which is 100 K above the ambient air temperature of 297.15 K.

Governing Equations
The problem is investigated numerically using a commercial CFD code FLUENT [20]. The governing equations are The local mass fraction of the pollutant species CO 2 , Y i is solved using the convection-diffusion equation: In the above equations, ρ, c p , k, µ denote the density, the specific heat, the thermal conductivity and the molecular dynamic viscosity of the mixture. The turbulent viscosity µ t is closed using the realisable k-ε model. J i is the diffusion flux of the species i. The gravity ⃗ g acts in the vertically downward direction (negative y-axis).

Computational Grid and Boundary Conditions
The geometric dimensions of the whole domain are expressed in terms of the diameter of the stack (D = 1 m). The length, breadth and the height of the computational domain is 1000D × 600D × 600D as shown in Fig. 1a. The height of flue stack is 30 D. The grid is meshed using the commercial grid generation software, GAMBIT [21]. As the computational domain is very large in size, an unstructured mesh is used with mesh refinement near the stack exit to capture the effects of buoyancy and momentum that are dominant near the stack exit. The computational grid in the XZ plane is shown in Fig. 2a, and this is extruded along the y-direction. The grid along the y-direction was fine up to Y = 40D (i.e. 10D above the stack height). It is then stretched to get a coarser grid along the remaining vertical height of the computational domain. The grid near the stack exit was ensured to be fine enough to capture the effect of buoyancy flux. As the path of the plume varies, with time, grid adaptation based on the variation of the plume path has been carried out for up to two levels of refinement as shown in Fig. 2b. The mesh adaptation is based upon the gradient of the mass fraction of the pollutant. The initial grid used has a total number of hexahedral grid cells of approximately 1.069 million, while the next two levels of refinement have a cell count of 1.22 million and 2.05 million cells. Figure 3 shows the mass fraction of the pollutant along a vertical line on the symmetric plane at 400D from the plume for three different grids. The mass fraction of the pollutant for the base grid was a little overpredicting, but the first level adaptation and the second level adaptation  Figure 1a shows the computational domain along with the boundary conditions. Wind inlet (face ABCD) depicts the velocity inlet wherein the wind with an average velocity of 2 m/s enters the control volume. The air entering is pollutantfree, i.e. the mass fraction of pollutants in the incoming air is zero. The face has a uniformly varying linear temperature gradient. Mathematically, these boundary conditions at the wind inlet can be described as The top, side1, side2 of the given control volume are pressure inlets. The temperature at this surface is specified based upon the uniformly varying temperature gradient and is a constant. Mathematically, the same can be expressed as follows. Ground is modelled as no-slip, stationary wall maintained at a temperature of 297.15 K. The pressure-based solver employing the segregated algorithm was used to solve the governing equations. The QUICK scheme [22] is used for convection term discretisation, and the SIMPLE algorithm [23] is used for pressure-velocity coupling. A convergence criterion of 10 -6 for the absolute error was used for all the variables (T, u, v, w).

Solution Methodology
The velocity and the temperature of the control volume were initialised with an average velocity of 2 m/s (wind velocity with a fully developed velocity profile using uniform velocity field) and 297.15 K (ambient temperature). The time step used for all the simulations is 1 s. The plume was initially introduced into the control volume. With the initial temperature gradient of 297.15 K set, for the entire domain, the plume showed a nonperiodic behaviour. After about two oscillations which

Validation of the Computational Model
The temperature variations of the plume with height (Y) at X = 10D and Z = 1.33D have been used as a benchmark study. The comparisons are made at a line drawn at 1.33D which passes through the peak in temperature. Fig. 4, where T s is the stack exit temperature, and T a is the ambient temperature. It can be seen that the results obtained are very close to the experimental results.

Plume Dispersion from Multiple Stacks
The calculation of the plume trajectory has not been defined consistently. Contini and Robins (2004) traced the path of the plume using the average height of the plume cross-section [7]. Mokhtarzedah-Dehghan et al. (2006) employed two methods to trace the path of the plume-using the midpoint of the line drawn between the top and bottom extreme of the plume cross-section in the Y Z plane and using the midpoint of the line drawn between the top and bottom of the plume in the plane of symmetry [6]. However, in the current study, the path of the plume is taken as the locus of the points with maximum pollutant mass fraction, which is a similar approach as used by Velamati et al. (2015).
The plume rise for every 100-time steps is plotted after periodicity is obtained post time t. Figure 5 a-f depict the plume rise for inline, 45° and non-inline stack for an interstack distance of 12D, at −0.2, −0.5 K/100 m. The plume rise at a given section is the y coordinate of the point of maximum plume concentration at a given time step. The plume rise for the inline stack configuration, α = −0.2 is shown in Fig. 5a. The figure shows an analogous plume trajectory till a downstream distance of 300 m for all the time steps. After 300 m, the plume begins to oscillate. The oscillations become amplified as the plume propagates downstream. This can be seen by comparing the oscillations of the plume at a downstream distance of 600 m and 800 m. The plume oscillates with an amplitude of about 45 m at a downstream distance of 600 m, while the amplitude of oscillations increases to 110 m at a downstream distance of 800 m from the flue stack. The similar kind of behaviour of the plume oscillations is observed in Fig. 5b-d. However, the amplitude of the oscillations varies from Fig. 5a, depending on the stack layout.

Effect of the Stack Layout Configuration on Plume Characteristics
In this section, plume characteristics-shape, plume rise and plume dispersion-are studied for threestack configurations namely inline, 45° and noninline for negative temperature gradients α = −0. In general, it can be observed that for a given temperature gradient for an inline stack layout of flue exits, the plume rise was the highest compared with the 45° layout and the noninline cases. In the inline configuration of the stack exits, the plume emanating from the exit located upwind shields the plume expelled from the down-wind flue exit. This enhances the efficiency of mixing of the plumes. The shielding effect was also observed by Bornoff and Mokhtarzadeh-Dehghan (2001) and Velamati et al. (2015) [4,5]. The first author reported the same effect in case of an inline or tandem layout of a multi-flue stack in the small-scale wind tunnel study over cooling towers. Velamati et al. stated the same in the inline layout of the flue exits under neutral and positive temperature gradient conditions. It can be observed that for an inline and non-inline layout, the plume trajectory remains the same for a downwind distance of about 300 m, after which the plume is seen to oscillate. This is because until 300 m, the plume rise is governed by the initial momentum of the plume imparted to the plume, as it is expelled from the stack. After this downstream distance, the plume rise becomes dependent upon the buoyancy force imparted by the atmospheric temperature gradient. The plume rise at a given downstream distance is higher for the inline configuration compared with its 45° and non-inline counterpart.
As seen in Fig. 6, the temperature gradient begins to fluctuate beyond a downstream distance of 600 m. A plausible reason for the same is that in the inline and 45° stacks layout, there is complete and partial shielding of the plume, respectively. The shielding effect allows the plume located downwind to rise undisturbed to a higher altitude in the inline and 45° configurations than its non-inline counterpart.
The plume rising to the higher layers of the atmosphere heats the air present there due to the heat transfer between the plume and the air. This, in turn, causes the cold air to be pushed down to neutralise the temperature gradient, set up by the rising plume and the heated air. The movement of the cooler air is as observed in the first and second columns of Fig. 6 at downstream distances of 600 m. In the non-inline case, the two plumes inter-mix with the help of a pair of counter-rotating vortices, and there is no shielding in this kind of a stack layout. The two plumes emitted come into direct contact with the incoming air, and this curbs the plume rise. The maximum plume rise in the case of the noninline configuration is about 250 m, for α = −0.2 K/100 m. The lower negative value of the atmospheric temperature gradient (−0.2 K/100 m) prevents the plume from oscillating considerably. Therefore, due to the lower rise of the plumes, there is a lesser amount of heat exchange between the plumes and the atmosphere. This allows the temperature gradient to maintain itself for a larger downstream distance of up to 900 m (as observed in the third column of Fig. 6). Figure 7 depicts the pollutant mass fraction contours and temperature gradient of −0.5 K/100 m. The inline layout of the flue stacks combined with the lower value of the temperature gradient allows this configuration to have the highest plume rise of 530 m (shown in Fig. 5b) compared with the 45° and the non-inline layouts, which show a plume rise of 380 m and 480 m respectively. Also, the oscillations are higher in this case. The initial rise of the plume draws cold air from high altitudes closer to the ground. The region of cold air combined with its lowtemperature equivalent in the higher elevations causes the plume to rise and fall alternately. The temperature gradient in the mid-plane is entirely different from the atmospheric temperature gradient (−0.5 K/100 m). A lower value of temperature gradient α = −0.5 K/100 m allows for a greater value of plume rise height, for a given configuration and a given inter-stack distance (12D, in this case). This can be observed by comparing the first, second and third columns of Figs. 6 and 7 and by comparing the first and second columns of Fig. 5. The initial rise of the plume is due to the momentum imparted as it exits from the stack. As the plume moves downstream, the rise and fall are governed by the temperature gradient prevalent in the atmosphere. Figure 8 shows the concentration contours of the pollutant mass fraction, at distances of 200 m and 500 m down-stream, for the inline, 45° and non-inline layout. It can be plainly seen from the graphs that for the inline configuration, the plume rise is maximum at any given time instant. Another important feature that can be compared in the three layouts is in the shape of the cross-section of the plume. The inline layout of the stack which depicts the shielding effect shows symmetric yet narrow cross-section of the plume.
The lateral displacement of the plume is low in the inline case (as observed in the first column of Fig. 8). The 45° layout shows an asymmetric cross-section. This is because, in the 45° layout, the plume released from the upstream stack forms a rotating vortex as it moves downstream. When this rotating vortex meets the second plume being released from the downstream stack, it tends to form an asymmetric pair of counter-rotating vortices. The presence of these asymmetric counter-rotating vortices was also observed by Contini et al. (2004) in the case of the neutral atmospheric temperature gradient condition. The plume released from the upstream stack has a greater vortex strength at the plane of downwind stack. Because of this, the pollutant from the downwind stack is pulled towards the vortex generated from the upwind stack thereby increasing the asymmetric nature of the plume cross-section, as it propagates downstream. It can be clearly seen by comparing the figures in the second column of Fig. 8, for a given time instant.
The last column of Fig. 8 shows a wider and symmetric cross section of the plume, as shown in the case of the non-inline layout. The two plume exits are separated by a distance of 12D, which allows the plumes to be released simultaneously, forming a pair of counter-rotating vortices due to the intermixing of the plumes released from the two stacks. These are of equal strength at any distance downstream from the stack exit.
Another notable feature can be seen in the temperature contours in Fig. 8. As the plume initially rises, the cooler air from the higher altitude tends to move in and get stratified at lower altitudes. This layer of colder air, along with its counterpart in the upper layers of the atmosphere, causes the plume to move downwards and also sustains the oscillations that ensue.

Effect of Inter-stack Distance on Plume Characteristics
Along with the atmospheric temperature gradient and the exit stack configuration, the effect of the distance between the two flue stacks plays a crucial role in the dispersion of plumes. Figures 9 and 10 compare the plume rise for a stack layout of 45° separated by inter-stack distances of 12D and 22D. By comparing the first and second columns of Fig. 9, it can be observed that for an inter-stack distance of 12D, the 45° orientation exhibited a greater amount of oscillations than its 22D counterpart. In the first column of Fig. 9, the oscillations of the plume cause the temperature gradient to change in the higher altitudes, in the case of 12D. However, due to the lesser amount of oscillations in the case of 22D, the temperature gradient present in the atmosphere maintains itself for a larger downstream distance. Figure 11a and d show the plume rise of the 45° stack configuration in the two temperature gradients, for an inter-stack distance of 12D and 22D. Comparing  Fig. 11a, c, it is observed that for α = −0.2 K/100 m, Fig. 9 Pollutant contours and temperature gradients along the symmetric plane of the 45° configurations at inter-stack distances of 12D, 22D at α = −0.2 K/100 m for various time instants. (Contour lines represent pollutant mass fractions 5 × 10 -8 , 2 × 10 -7 , 3 × 10 -7 , 4 × 10 -7 and 5 × 10 -7 ) ◂ the plumes emitted from a 12D configuration tend to oscillate more than the 22D counterpart. This is because, in the case of 22D, a larger distance is travelled by the upwind plume, before it comes in contact with the second plume downwind. The larger distance allows the upwind plume to be diluted to a much lesser and thereby mix much more efficiently with the downwind plume. However, in the case of a stack separation of 12D, the plume traverses a lesser distance before it comes in contact with the second plume downstream. This lesser distance does not allow the plume to dilute considerably. Therefore, mixing becomes protracted, in the case of 12D. Also, the temperature gradient of −0.2 K/100 m plays a major role in confining the plumes to the lower altitudes of the atmosphere. It does not allow the plume to oscillate considerably and thereby prevents dilution caused by the oscillations.
A lower value of the atmospheric temperature gradient −0.5 K/100 m causes the oscillations to increase for any stack layout (Fig. 10). As the plume initially rises, cooler air from the higher layers of the atmosphere moves to the lower layers of the atmosphere and gets stratified there. This region of lower temperature along with its lower temperature counterpart in the higher layers of the atmosphere causes the plume to move downward and also sustain the oscillations that ensue. Even in the case of the inline and non-inline layouts, the maximum plume rise height is higher for an inter-stack distance of 12D compared with 22D at a given temperature gradient. Figure 11b and d show the plume rise for inter-stack distances of 12D and 22D at −0.5 K/100 m. It can be seen that the maximum plume rise is higher in the case of 22D compared with 12D. The maximum plume rise in the case of 12D is 400 m while that of 22D is about 550 m. Figure 12 depicts the concentration contours of the pollutant at various distances downstream from the stack along an XZ plane, which also depicts the temperature gradient of −0.5 K/100 m. The higher rise of plume causes the temperature gradient of the atmosphere to change from −0.5K/100m. This is because the plume exchanges energy with the cooler air present in the high altitudes. This exchange of energy heats the air present at the higher altitudes of the atmosphere and thereby changes the temperature gradient value. Another notable feature is the shape of the cross-section of the plume. The cross-section in the case of an inter-stack distance of 22D (second column of Fig. 12) is wider than that of 12D (first column of Fig. 12). This is due to the greater stack separation in the 22D case. Also, the plume rise at a given cross section for the 22D case is higher than the 12D counterpart.

Conclusions
The three single flue exit orientations studied were inline and non-inline and 45°, for the two temperature gradients −0.2 and −0.5 K/100 m. Based upon these, the effect of the stack configuration and the negative atmospheric temperature gradient upon the plume rise and dispersion have been enlisted below: 1. At a given atmospheric temperature gradient, the inline configuration shows the highest plume rise compared with the 45° layout and the non-inline. This is because of the shielding of the downwind jet, by its upwind counterpart which allows slower cooling of the plumes emitted by the downwind plume and also allows for higher plume rise. 2. The plume mixing occurs with the formation of a pair of counter-rotating vortices. Symmetric vortices are observed for the non-inline layout of the flue stacks, whereas asymmetric vortices are observed in the 45° layout of the flue stacks.
The maximum plume rise is higher for an inter-stack distance of 12D in the case of all three layouts, i.e. inline, non-inline and 45°. This is because of the lesser distance of 12D, travelled by the plume emitted from the upwind stack before it mixes with the downwind plume, compared with its 22D counterpart.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons . org/licen ses/by/4.0/.