Large eddy simulation of radiation fog: impact of dynamics on the fog life cycle

. Large eddy simulations (LESs) of a radiation fog event occurring during the ParisFog experiment are studied with a view to analyse the impact of the dynamics of the boundary layer on the fog life cycle. The LES, performed with the Meso-NH model at 5 m resolution horizontally and 1 m vertically, and with a 2-moment microphysical scheme, includes the drag effect of a tree barrier and the deposition of droplets on vegetation. The model shows good agreement with measurements of near-surface dynamic and thermody-namic parameters and liquid water path. The blocking effect of the trees induces elevated fog formation, as actually observed, and horizontal heterogeneities during the formation. It also limits cooling and cloud water production. Deposition is found to exert the most signiﬁcant impact on fog prediction as it not only erodes the fog near the surface but also modiﬁes the fog life cycle and induces vertical heterogeneities. A comparison with the 2 m horizontal resolution simulation reveals small differences, meaning that grid

Abstract. Large eddy simulations (LESs) of a radiation fog event occurring during the ParisFog experiment are studied with a view to analyse the impact of the dynamics of the boundary layer on the fog life cycle. The LES, performed with the Meso-NH model at 5 m resolution horizontally and 1 m vertically, and with a 2-moment microphysical scheme, includes the drag effect of a tree barrier and the deposition of droplets on vegetation. The model shows good agreement with measurements of near-surface dynamic and thermodynamic parameters and liquid water path. The blocking effect of the trees induces elevated fog formation, as actually observed, and horizontal heterogeneities during the formation. It also limits cooling and cloud water production. Deposition is found to exert the most significant impact on fog prediction as it not only erodes the fog near the surface but also modifies the fog life cycle and induces vertical heterogeneities. A comparison with the 2 m horizontal resolution simulation reveals small differences, meaning that grid convergence is achieved. Conversely, increasing numerical diffusion through a wind advection operator of lower order leads to an increase in the liquid water path and has a very similar effect to removing the tree barrier. This study allows us to establish the major dynamical ingredients needed to accurately represent the fog life cycle at very high-resolution.

Introduction
Despite a long-standing interest in understanding fog processes, uncertainties still exist in the physical mechanisms driving fog variability. Forecasting fog remains a challenge because of the diversity of mechanisms involved during the fog life cycle and their interactions: local flow, turbulence, radiation, microphysics, aerosols and surface effects. Several field experiments have been carried out since the 1970s and have contributed to the important progress made in understanding fog processes. Noteworthy studies include campaigns at Cardington in the UK (Roach et al., 1976;; Fog-82 in Albany, New York (Meyer et al., 1986); Lille 91 in France (Guedalia and Bergot, 1994); a campaign in the Po valley in Italy (Fuzzi et al., 1998) and ParisFog in France (Haeffelin et al., 2010).
Many important features of fog have also been characterized using one-dimensional (1-D) modelling (Bergot et al., 2007;Tardif, 2007;Stolaki et al., 2015). However, to study some aspects of the characteristics of a fog layer, it has become necessary to explicitly simulate turbulent motion in 3-D as shown by Nakanishi (2000), who was the first to use a large eddy simulation (LES) for fog. LES is a turbulence modelling technique in which most of the energy-containing eddies are explicitly resolved while eddies smaller than a certain cutoff size, usually taken equal to the grid spacing, are parameterized by a turbulence scheme. Since then, Porson et al. (2011) have explored the static stability in a fog layer, and Bergot (2013) has shown the various organized structures occurring in a fog layer, which cannot be resolved in 1-D. Thanks to these studies, the dynamical characteristics Published by Copernicus Publications on behalf of the European Geosciences Union.
of radiation fog are more clearly identified during the three stages of the fog life cycle defined by Nakanishi (2000): the onset, development and dissipation phases. During the formation phase, small banded structures, identified by Bergot (2013) as Kelvin-Helmholtz (KH) billows, occur in the middle of the fog layer in dynamical and thermodynamical fields. They are sometimes associated with a burst of turbulent kinetic energy (TKE) (Nakanishi, 2000;Bergot, 2013) but this is not always the case (Porson et al., 2011). During the development phase, the main dynamical processes relocate to the top of the fog layer and are associated with the maximum of TKE and horizontal rolls (Bergot, 2013). During the dissipation phase, coupled processes between the ground and the top of the fog layer explain the spatial variability in fog (Bergot, 2015) but the link between dynamics and microphysics has not been explored specifically in these LES studies.
The quality of the LES depends on the horizontal and vertical resolutions. Beare and MacVean (2004) demonstrate that simulations in stable conditions converge at a 2 m horizontal resolution. Very high vertical resolution is also essential for representing the divergence of the radiative fluxes in the first few metres above the surface and therefore to produce the radiative cooling necessary for the formation of fog (Duynkerke, 1999;Tardif, 2007).
So far, most fog LES studies have considered homogeneous canopies. Only Bergot et al. (2015) have accounted for the effect of surface heterogeneities, such as buildings, on radiation fog. Other studies, such as those by Zaïdi et al. (2013) or Dupont and Brunet (2008), have considered the impact of forests on turbulence structures but not for fog situations. In this study, we will explore an LES of a fog case that was observed during ParisFog and was strongly influenced by trees.
Few fog LES studies are based on sophisticated 2-moment microphysical schemes, which allow the impact of aerosols on the radiation fog life cycle to be represented. Maalick et al. (2016) studied the effects of aerosols on radiation fog with an LES but in a 2-D configuration that could present some limitations for the dynamical patterns of the fog layer. Additionally, most of the studies using 1-or 2-moment microphysical schemes fail to reproduce realistic liquid water contents (LWC) as they tend to overestimate values near the ground. For instance, Zhang et al. (2014) simulated N c = 800 cm −3 and LWC = 0.4 g m −3 , and Stolaki et al. (2015) simulated N c = 250 cm −3 and LWC = 0.34 g m −3 near the surface, both a in 1-D configuration. These values are outside the range found by Mazoyer et al. (2016) for the same site. So the question of a possible missing mechanism arises, the inclusion of which might improve the modelling of microphysical fields. Some aspect of deposition that relates to the interaction with the ground surface is important as already shown by Price and Clark (2014) with measurements and von Glasow and Bott (1999) or Zhang et al. (2014) with 1-D simulations.
The goal of this study is to better understand the physical processes dominating the fog life cycle at a complex site.
LES modelling at very high resolution (1 m vertically and 5 m horizontally) is used with surface heterogeneities (barrier of trees) and a 2-moment microphysical scheme. In order to establish the main ingredients driving the fog life cycle, and to evaluate how dynamics affects the evolution of fog, sensitivity simulations are conducted. To our knowledge, this is the first time that an LES study of radiation fog has been performed at such high-resolution with a sophisticated microphysical parameterization scheme while considering the effect of heterogeneities such as forests on the fog dynamics.
Section 2 presents the measurement set-up and the observed case and describes the numerical model. The reference simulation is analysed in Sect. 3, and Sect. 4 is devoted to sensitivity tests. Finally, some conclusions are drawn and perspectives suggested in Sect. 5.
2 Experimental design and model description

Measurements set-up
The selected fog event was observed on 15 November 2011 during the ParisFog field campaign (Haeffelin et al., 2010) at the SIRTA (Site Instrumental de Recherche par Télédétection Atmosphérique) observatory (48.713 • N, 2.208 • E). The objective of the ParisFog campaign during three winters from 2010 to 2013 was to better understand the radiative, thermodynamic, dynamic and microphysical processes occurring during the fog life cycle. The site where the instrument platform was installed was a semi-urban area with mixed land cover including forest, lake, meadows and shrubs next to a built-up area. As shown in Fig. 1a, the instrumented zone was located near a forest area. Zaïdi et al. (2013) demonstrated the impact of the tree barrier on the observed flow when the wind was blowing from the barrier of trees over the instrument location, as in our case study. This fog case has previously been studied by Stolaki et al. (2015) using a 1-D model.
Temperature and humidity sensors were located at heights between 1 and 30 m on an instrumented mast, with uncertainties of 0.2 K in temperature and 2 % in relative humidity (RH). Wind speed was measured by two ultrasonic anemometers at 10 and 30 m a.g.l. (above ground level) on the same mast. Radiative fluxes were measured at a height of 10 m with 5 and 4 W m −2 uncertainties for downward and upward fluxes, respectively. Additionally, radiosondes were launched by Météo-France twice a day from Trappes (48.7 • N, 2 • E), situated 15 km to the northwest of SIRTA.
Aerosol particle measurements were performed using a scanning mobility particle sizer (SMPS) measuring dry aerosol diameters between 10.6 and 496 nm every 5 min, and by a CCN (cloud condensation nuclei) chamber that gave the CCN number concentration at different supersaturations from 0.1 to 0.5 % (Roberts and Nenes, 2005). An RPG-HATPRO water vapour and oxygen multi-channel mi-Atmos. Chem. Phys., 17, 13017-13035, 2017 www.atmos-chem-phys.net/17/13017/2017/ crowave profiler was used to measure the liquid water path (LWP) with an error of up to 20 g m −2 (Lohnert and Crewell, 2003). Measurements of dewfall and fog-droplet deposition were not taken.

Presentation of the observed case
Radiation fog formed at 02:00 UTC on 15 November 2011 and dissipated at the ground around 10:00 UTC on the following morning. Conditions favouring fog were due to a ridge at 500 hPa centred over the North Sea and anticyclonic conditions near the surface. One of the features of this event was the initial formation of a cloud layer at 150 m a.g.l., followed 30 min later by fog occurring at the surface. As underlined by Stolaki et al. (2015), this characteristic is very common at SIRTA and 88 % of the radiation fog events during the field experiment followed a similar pattern. However, these events are not classified as stratus lowering as they were followed rapidly by formation of fog at the surface. A delay of 30 min between the formation at 150 m height and at the ground seems too short to be a stratus lowering, which is mainly driven by the evaporation of slowly falling droplets that cool the sub-cloud layer (Dupont et al., 2012). This suggests that this type of radiation fog could be linked with, and specific to, the configuration of the SIRTA site. The fog case is presented following the three phases of the fog life cycle defined by Nakanishi (2000). Before the fog onset, between 22:00 and 02:00 UTC, the surface boundary layer was stable and a near-surface cooling was observed, inducing an increase in RH (Fig. 2). Between 00:00 and 01:30 UTC, the RH near the ground remained nearly constant around 97 %. Wind at 10 m height was light (speed around 1.8 m s −1 ) as was TKE, with small variability (Fig. 3). At 02:00 UTC, the attenuated backscatter coefficient measured by the lidar increased significantly at 150 m a.g.l. (not shown), revealing the formation of liquid water at this height, while the RH at the surface remained at 97 %. The cloud base height progressively subsided over the next 30 min, at which point it reached the ground. During this time, the near-surface temperature decreased by about 1 K in a stable stratification layer. At 02:30 UTC, the appearance of fog at the ground was associated with a temperature homogenization in the first 30 m, called temperature convergence by  and corresponding to a neutral stratification. The downwelling longwave (LWD) radiation flux increased progressively to 325 W m −2 during the development of the fog layer (Fig. 4).
During the fog development and mature phases, between 02:00 and 07:00 UTC, the near-surface layer remained quasineutral and potential temperature at the different levels remained constant. The temporal variability in the 10 m wind speed and TKE was greater from this period. Around 04:00 UTC, TKE at 10 m height increased significantly, from 0.4 to 0.7 m 2 s −2 , and then presented some variability around this value. The vertical gradient of TKE between 30 and 10 m remained positive. The sodar indicated that the fog-top height reached a maximum of 300 m a.g.l. during the mature phase (Stolaki et al., 2015;Dabas et al., 2012).
At the beginning of the dissipation phase, starting at 07:00 UTC, the surface temperature increased slowly (less than 0.5 K in 2 h) and then more rapidly after 09:00 UTC. At 10:00 UTC, the downward shortwave (SW) fluxes exceeded 100 W m −2 , while near-surface temperature had increased by 1 K compared to the pre-sunrise values. TKE at 30 m decreased during 08:00 to 10:00 UTC, while at 10 m TKE remained approximately constant.
The LWP value measured by the profiler was maximum around 07:30 UTC, at the beginning of the fog dissipation phase, with 70 g m −2 (Fig. 5). The non-zero values (5 g m −2 ) before the fog onset were within the error range of the measurement.

Model description 2.3.1 Presentation of the model
The non-hydrostatic anelastic research model Meso-NH (Lafore et al., 1998) (see http://mesonh.aero.obs-mip.fr) is used here in an LES configuration. The LES is based on a 3-D turbulent scheme with a prognostic TKE (Cuxart et al., 2000) and a Deardorff mixing length (Deardorff, 1980). The atmospheric model is coupled with the ISBA surface scheme (Interaction between Soil Biosphere and Atmosphere; Noilhan and Planton, 1989) through the SURFEX model (Masson et al., 2013). This scheme simulates the exchanges of energy and water between the land surface (soil, vegetation and snow) and the atmosphere above it. It uses five prognostic equations for deep temperature, deep soil water content, surface temperature, surface soil water content and water interception storage by vegetation.
In order to consider the impact of trees at the measurement site, we used the drag approach developed by Aumond et al. (2013) for a vegetation canopy. Both this study and Zaïdi et al. (2013) have shown that the drag approach gives better results than the classical roughness law when reproducing the turbulence downstream of a forest area. The drag approach consists of introducing an additional term into the momentum and TKE equations as follows: with α = u, v or TKE, where u and v are the horizontal wind components; C d is the drag coefficient, set to 0.2; and A f (z) is the canopy area density, representing the surface area of the trees facing the flow per unit volume of canopy. A f (z) is the product of the fraction of vegetation in the grid cell by the leaf area index (LAI) and by a weighting function representing the shape of the trees, as presented in Aumond et al. (2013). The trees introduced in the simulation domain for the land surface scheme correspond to Atlantic coast broad leaved trees. The model includes a 2-moment bulk warm microphysical scheme (Khairoutdinov and Kogan, 2000;Geoffroy et al., 2008), that considers droplet concentration N c and mixing ratio r c as prognostic variables for the fog. An additional prognostic variable N ccn is used to account for already activated CCN, following the activation scheme of Cohard et al. (2000c). The aerosols are assumed to be lognormally distributed and the activation spectrum is prescribed as the following: where N ccn is the concentration of activated aerosol, F (a, b; c; x) is the hypergeometric function, C (m −3 ) is the concentration of aerosols, and k, µ and β are adjustable shape parameters associated with the characteristics of the aerosol size spectrum such as the geometric mean radius (r) and the geometric standard deviation (σ ), as well as solubility of the aerosols (ε m ) and temperature (T ) (see below for the values in our case study). S max is the maximum of supersaturation for that grid box at a time step corresponding to dS dt = 0. The evolution of the supersaturation S includes three terms accounting for the effects of a convective ascent of vertical velocity w, the growth of droplets by condensation for the newly activated droplets, and radiative cooling, as in Zhang et al. (2014): where φ 1 (T ), φ 2 (T , P ) and φ 3 (T ) are functions of temperature and pressure. Following Pruppacher et al. (1998) and Atmos. Chem. Phys., 17, 13017-13035, 2017 www.atmos-chem-phys.net/17/13017/2017/ after simplification, S max can be diagnosed by with B the beta function and ρ w the density of water. Thus, the aerosols potentially activated are exactly those with a critical supersaturation lower than S max . The number of aerosols actually activated in a time step is the difference between the number of potentially activated aerosols and the number of aerosols previously activated during the simulation. The condensation/evaporation rate is derived using the Langlois (1973) saturation adjustment scheme. Cloud droplet sedimentation is computed by assuming Stokes law for the cloud droplet sedimentation velocity and assuming that the cloud droplet size distribution n c (D) fits a generalized gamma law: where λ is the slope parameter, depending on the prognostic variables r c and N c : α and ν are the parameters of the gamma law, and ρ a is the density of dry air. They were adjusted using droplet spectra measurements from the FM-100 database of our case study and were set at α = 1 and ν = 8. These parameters are also used for the radiative transfer. In addition to droplet sedimentation, fog deposition is also introduced to represent direct droplet interception by the plant canopies. In the real world, it results from the turbulent exchange of fog water between the air and the surface below, leading to collection (Lovett et al., 1997). In numerical weather prediction models (NWPs), this process is not usually included, e.g. in the French NWP model AROME (Seity et al., 2011), the physics of which comes from Meso-NH. As fog deposition is a newly introduced process, only a simple formulation is considered here as a first step in order to perform a sensitivity study. The fog deposition flux F DEP is predicted at the first level of the atmospheric model (50 cm height) for grassy areas, and over the 15 m height for trees, in a simplistic way following Zhang et al. (2014): where χ = r c or N c , and V DEP is the deposition velocity. In a review based on measurements and parameterization, Katata (2014) showed that V DEP values ranged from 2.1 to 8.0 cm s −1 for short vegetation. Here V DEP is assumed to be constant, equal to 2 cm s −1 . A test of sensitivity to this value is presented below. Water sedimentation and deposition amounts are input to the humidity storage of the surface model. A more complete approach in a further study would include a dependence of V DEP on momentum transport as in von Glasow and Bott (1999) and also on LAI. The radiative transfer is computed with the ECMWF radiation code, using the Rapid Radiative Transfer Model (RRTM; Mlawer et al., 1997) for longwave (LW) radiation and Morcrette (1991) for shortwave radiation. Cloud optical properties for LW and SW radiation take account of the cloud droplet concentration in addition to the cloud mixing ratio. For SW radiation, the effective radius of cloud particles is calculated from the 2-moment microphysical scheme, the optical thickness is parameterized according to Savijärvi et al. (1997), the asymmetry factor is from Fouquart et al. (1991) and the single scattering albedo from Slingo (1989). For LW radiation, cloud water optical properties refer to Savijärvi et al. (1997).

Simulation set-up
For the reference simulation (denoted as REF), the horizontal resolution is 5 m over a domain of 1 km × 1 km and 126 vertical levels are used between the ground and the top of the model at 1500 m. The vertical resolution is 1 m for the first 50 m and increases slightly above this height. Momentum is advected with a 4th order centred scheme (denoted as CEN4TH), whereas scalar variables are advected with the PPM (piecewise parabolic method) scheme (Colella and Woodward, 1984). The time step is 0.1 s. The domain of simulation is presented in Fig. 1b. It has a tree barrier 15 m high and 100 m wide perpendicular to the wind direction and the rest of the domain surface is composed of grass. The lateral boundary conditions are cyclic. The radiation scheme is called every second.
The simulation begins at 23:20 UTC on 14 November 2011 before any fog has formed and it lasted for 12 h. Temperature, humidity and wind speed vertical profiles were initialized with data from the radiosonde launched from Trappes. Meteorological conditions at Trappes can differ slightly from those at the SIRTA site. Therefore wind, temperature and humidity were modified in the nocturnal boundary layer up to 400 m a.g.l. to fit the data recorded at the 30 m meteorological mast at the SIRTA site, as illustrated in Fig. A1. The soil temperature and moisture were given by the soil measurements, corresponding to a surface temperature of 276 K and a soil moisture of 70 %. Following the profiles from soundings, a geostrophic wind of 8 m s −1 was prescribed, without any other forcing. To generate turbulence, a white noise of 0.5 K was applied in the first 100 m.
It was also necessary to characterize the aerosol size spectrum for Eq. (2). The supersaturations reached in fog were lower than 0.1 % meaning that the CCNC (cloud condensation nuclei counter) measurements were not directly usable, as shown by Hammer et al. (2014) and Mazoyer et al. (2016).
However, by using the kappa-Köhler theory and the SMPS observations, the aerosol concentrations at supersaturations under 0.1 % can be retrieved if the aerosol hygroscopicity (κ) at these supersaturations is known. This method, proposed by Mazoyer et al. (2016), was applied to our case study in the hour before fog onset. Thus, above 0.1 % supersaturation, the activation spectrum was found from observations and below 0.1 % it was computed. This computed activation spectrum is fit according to Eq. (2) (Fig. A2a), which corresponds to the size distribution of aerosol particles (C = 2017 cm −3 , σ = 0.424, r = 0.1, ε m = 1) in red in Fig. A2b. This does not match the measured distribution (in black) or the lognormal distribution fit on the accumulation mode (in blue) because the Cohard et al. (2000c) formulation was not developed for fog with low supersaturation. Deducing the activation spectrum from measurements provides the exact solution.

The reference simulation
The performance of the REF simulation will first be examined, based on a comparison with observed values of thermohygrometric, dynamic and radiative parameters near the ground and LWP. Considering that the REF simulation reaches good agreement with observations, the vertical evolution and horizontal variability in the simulated fog will be characterized during the different phases of the fog life cycle. It should be emphasized that observations localized at one point will be compared to simulated fields averaged over a horizontal area located downstream of the tree barrier (blue contour area of Fig. 1b), which is representative of the measurement area. We will see that the simulation domain is divided into four parts with significant differences among them but with similar characteristics within each one.

Comparison to measurements
Figure 2a and c show the time series of near-surface observed and simulated temperature and RH. At the initialization of the simulation, near-surface temperatures are in agreement with the observations while RH is very slightly underestimated. During the cooling before fog onset, the model develops a layer that is too stable, especially in the first 5 m, between 00:00 and 01:00 UTC. The convergence of temperature is simulated with a 30 min delay with respect to the observations.
Based on RH near the surface, the fog starts to appear around 02:00 UTC. Between 04:30 and 09:00 UTC, simulated and observed temperature are in fairly good agreement, with a quasi-neutral near-surface layer. The fog starts to dissipate from the ground at 09:00 UTC, approximately 1 h ahead of the local observation. This time discrepancy induces a slight overestimation of near-surface temperature, which is less than 0.5 K at 11:00 UTC. Nevertheless, the negative temperature gradient near the surface representative of the  Fig. 1b). The grey shaded areas represent the error for the observational curves. development of the convective boundary layer is quite well reproduced after the beginning of the dissipation.
Dynamical fields at 10 and 30 m are fairly well reproduced by the model (Fig. 3 in red): the 10 m wind speed (Fig. 3a) is in good agreement with observation throughout the simulation. Until 03:00 UTC, a quasi-linear increase in TKE is produced by the model with a higher TKE at 10 m a.g.l. than at 30 m contrary to observations (Fig. 3b). Around 03:00 UTC, a more sudden increase in TKE occurs, as in the observations but 30 min before and with a lower magnitude. Then the simulated TKE remains almost constant around 0.7 m 2 s −2 from 04:00 UTC onwards, with a slightly higher variability than before. The model develops similar TKE values at 10 and 30 m, while observed values are higher at 30 m.
Considering the radiative fluxes (Fig. 4), the increase in the LWD flux associated with fog onset is simulated with a delay of 30 min, meaning that there is a delay in the simulated formation of fog at elevated levels. This delay is corroborated by the LWP evolution (Fig. 5). After that, the LWD flux of 325 W m −2 is correctly reproduced, indicating that the temperature and the optical thickness of the fog are fairly well simulated. Observations develop a difference of 8 W m −2 between longwave upwelling (LWU) and LWD during the fog life cycle, but the model fails to reproduce this difference, leading to a slight underestimation of LWU. If the measurements do not contain any errors, this probably means that the Atmos. Chem. Phys., 17, 13017-13035, 2017 www.atmos-chem-phys.net/17/13017/2017/ radiative properties of the simulated surface are not perfectly represented. A test on the emissivity of the surface (1 instead of 0.96) had no impact on the radiative fluxes, suggesting that the soil temperature was probably underestimated. The simulated LWP presents a maximum of the same magnitude than the observation, but occurring 2 h before. After sunrise (06:59 UTC), the downward and upward SW fluxes are overestimated by up to 15 W m −2 . LWD is slightly underestimated in a similar way due to the advanced dissipation time, as well as the LWP. The comparison between the REF simulation and observation for the set of parameters shows fairly good agreement, even though there are a few discrepancies. The main discrepancies concern the fog life cycle with an underestimation of the effect of elevated fog formation, and an advance of 1 h in the dissipation time. These elements are probably partly due to the semi-idealized representation of the SIRTA surface in the simulation, and also to the comparisons with point observations, given the horizontal variability that we will see below. They are felt to be acceptable and we can therefore consider that the REF simulation can be used to explore the processes driving the fog life cycle and to conduct sensitivity tests.  Fig. 1b). The grey shaded area represents the error for the observational curve.

Vertical evolution
First the vertical evolution of the fog is analysed. Figure 6 represents the time variations in vertical profiles of r c and N c , the radiative cooling rate and the vertical velocity in the updraughts. The vertical and temporal variations in simulated N c can be studied as the LWP is realistic, but values of N c must be carefully considered as a first comparison to near-surface measurements clearly shows an overestimation of simulated values. Figure 7a, c and d represent the time variation for total TKE (resolved plus subgrid) and dynamical and thermal production of TKE for the REF simulation, all averaged over the horizontal area downstream of the tree barrier. A first feature is that subgrid kinetic energy is 1 order of magnitude lower than resolved kinetic energy (not shown). This means that the 5 m horizontal resolution allows an LES approach as most of the eddies are resolved.
The evolution of r c serves as a basis for decomposing the fog life cycle into the three phases: formation, between 02:00 and 03:00 UTC, until the fog becomes optically thick; development, between 03:20 and 08:20 UTC, until r c at upper levels of the fog layer begins to decrease, and dissipation from 08:20 UTC (Fig. 6a).
Before the fog onset and during the formation phase, the TKE is small and spread over a 30 m layer that deepens slowly because of the tree barrier (Fig. 7a). TKE mainly occurs by dynamical production, which presents maxima at two levels: near the surface and at 15 m height due to the trees (Fig. 7c). Thermal production is negative because of the thermal stratification (Fig. 7d). Radiative cooling near the ground (Fig. 6c) and mixing by the tree drag effect are the ingredients that allow fog to appear at the same time over a 30 m deep layer (Fig. 6a). Then the mixing by the tree barrier causes the fog layer to develop vertically at greater heights (Fig. 6a). Hence, the effect of elevated formation is reproduced, even though the height of fog onset is underestimated (3) for supersaturation evolution with the two source terms depending on vertical velocity and radiative cooling, activation of fog droplets during the fog formation is mainly produced by radiative cooling at the top of the fog layer ( Fig. 6b and c). At the beginning of the development phase (around 03:00 UTC), when the fog depth reaches approximately 80 m, it becomes optically thick to LW radiation. At that time, TKE increases significantly by dynamical production (Fig. 7a and c), in agreement with the findings of Nakanishi (2000), which indicates a dynamical change. The optical thickness of the fog layer causes strong radiative cooling at the top of the layer (greater than 5.5 K h −1 in absolute value; Fig. 6c), and r c values increase in the upper part of the fog layer. Hence, the fog top becomes the location of the dominant processes. Radiative cooling induces small downdraughts and buoyancy reversal. In addition to the vertical ve-locity of the updraughts, now higher than 0.2 m s −1 throughout the fog layer, a second maximum of droplet concentration of 1100 cm −3 occurs in the upper part of the fog layer around 03:20 UTC. The sudden optical thickening corresponds to the increase in surface LWD to 320 W m −2 (Fig. 4) and to maximum cooling at the ground (Fig. 2a). At the same time, temperatures converge in the vertical levels near the ground, showing the effect of fog on the stratification as analysed by .
During the development phase, the top of the fog layer is characterized by vertical wind shear inducing positive dynamical production of TKE, while small values of positive thermal production appear at the top due to buoyancy reversal. In the lowest 40 m of the fog layer, the drag effect of the trees induces values of kinetic energy higher than 0.6 m 2 s −2 . The maximum of r c continues to increase in the upper part of the fog layer until 05:00 UTC, reaching 0.37 g kg −1 at 120 m (Fig. 6a). At the same time, LWD surface fluxes remain constant while the fog layer continues to deepen and the LWP continues to increase until 05:00 UTC (Fig. 5).
Around 05:00 UTC, a change occurs in the development of the fog layer: it continues to thicken but at a slower rate, Atmos. Chem. Phys., 17, 13017-13035, 2017 www.atmos-chem-phys.net/17/13017/2017/ while the LWP begins to decrease in the simulation. This change of growth at the top of the fog layer is associated with a warming in the fog layer (not shown) and a decrease in the maximum radiative cooling near the top which spreads over a greater depth (Fig. 6c). This also corresponds to an increased number of resolved updraughts and downdraughts near the top (Fig. 6d). The variability in the fog depth also becomes stronger, in connection with fog-top waves as we will see below. This change of growth seems to be linked to the fact that the fog layer reaches the top of the nocturnal boundary layer, meeting stronger temperature, humidity and wind gradients. This increases the top entrainment process, limiting the deepening of the fog layer. With the decrease in top radiative cooling, cloud droplet concentration becomes more homogeneous in the fog layer, except near the ground where it decreases by deposition. The cloud mixing ratio also begins to decrease near the ground (Fig. 6b). The beginning of the dissipation phase in the simulation (around 08:20 UTC) is preceded by the beginning of solar radiation, and a divergence between surface LWU, which starts to increase, and surface LWD, which starts to decrease (Fig. 4). The dissipation of the fog begins at the surface, and the fog lifts into a stratus layer. The radiative heating of the surface induces the convective structure of the fog as vertical velocity in the updraughts increases ( Fig. 6c and d) and thermal production of TKE becomes significantly positive (Fig. 7d). Additionally, after sunset, downdraughts at the top of the fog layer increase the amount of solar radiation reaching the ground and this feeds the heating at the base of the fog layer. Hence, near the ground, both thermal and dynamical effects contribute to the production of TKE, and to a deepening of the TKE layer to 60 m. The height of the fog top continues to increase as it is driven by radiative and evaporative cooling, which induces vertical motions and top entrainment. Although the mixing ratio decreases at all levels, droplet concentration increases sharply when the fog layer lifts from the surface (Fig. 6b). As the cloud evolves into a stratus layer, droplet activation is no longer induced by radiative cooling at the top of the fog layer but by updraught vertical velocity at all cloud depths, and especially www.atmos-chem-phys.net/17/13017/2017/ Atmos. Chem. Phys., 17, 13017-13035, 2017 near the stratus base. The stronger vertical velocity activates more droplets for the same water content. Droplets become smaller and more numerous, preventing the droplet sedimentation process and limiting the decrease in LWP. Moreover, the deposition process is no longer active as there are no cloud droplets at the surface. We will now consider the horizontal heterogeneity of the fog layer.

Horizontal variability
To better characterize turbulent structures and the impact of trees on the fog layer, the horizontal variability in the fog layer is examined. Figure 8 presents horizontal and vertical cross sections of wind speed, cloud mixing ratio, potential temperature and TKE at 02:10 UTC during the formation phase. The tree barrier tends to block the flow upstream. It enhances the turbulence by wind shear downstream, accelerating the flow near the ground and creating longitudinal structures in the direction of the wind. Ascents occur upstream and small subsidences downstream (up to 2 cm s −1 , not shown). The subsidences bring warmer and dryer air from above to the ground. Therefore, structures of stronger wind near the ground downstream of the trees coincide with structures of warmer, clear air as they delay fog formation. The fog forms at the surface upstream of the trees, and 500 m downstream, while it appears first at elevated levels over the intermediate area between the trees and downstream (Fig. 8d). The fog takes about 1 h to cover the entire domain at ground level. Thus, heterogeneity of the surface vegetation explains heterogeneities in fog onset over the SIRTA site, as well as the fog property of developing first at elevated levels.
After the formation phase, the base of the fog layer is at the ground over the whole domain. These results are in agreement with the effects of buildings on fog studied by Bergot et al. (2015) who found a 1.5 h period of heterogeneity of fog formation over the airport area.
During the development phase, as shown in the vertical cross sections of Fig. 9 at 06:20 UTC, horizontal rolls appear at the top of the fog layer and are associated with dynamical production of TKE by shear. They are aligned almost perpendicularly to the mean wind direction (not shown). These structures correspond to KH instability, previously observed Atmos. Chem. Phys., 17, 13017-13035, 2017 www.atmos-chem-phys.net/17/13017/2017/ by Uematsu et al. (2005) and modelled by Nakanishi (2000) and Bergot (2013). They have depths corresponding to about one-third of the fog layer height, as in Bergot (2013), and a horizontal wavelength of the order of 500 m. These horizontal rolls explain the oscillations at the top of the fog layer visible in Figs. 6 and 7. They become well marked from 05:00 UTC when the increase in depth of the fog layer begins to slow down, as the fog layer reaches the top of the nocturnal boundary layer, meeting stronger wind gradients. The horizontal rolls induce strong horizontal variability in the cloud mixing ratio near the top of the fog, with larger values in the ridges of the fog-top rolls, and smaller ones in the troughs (Fig. 9a). Local updraughts occur upstream of the crest of the wave, and downdraughts downstream (both up to 1.2 m s −1 ; Fig. 9d). The maximum of droplet concentration occurs near the top of the fog layer (Fig. 9b) in the radiative cooling layer (Fig. 9c), and preferentially upstream of the crest of the wave rather than downstream, in the ascent area, where the droplets are preferentially activated and transported. These extrema of droplet concentration do not appear in Fig. 6 as they are hidden by the spatio-temporal average. Deposition velocity equal to 8 cm s −1 DX2 Horizontal resolution = 2 m WE3 3rd order WENO advection for momentum WE5 5th order WENO advection for momentum WENO: Weighted Non-Oscillatory, Shu (1998).
Inside the fog layer, the radiative cooling is negligible, while vertical velocity presents strong spatial heterogeneities. Maxima of supersaturation appear to be strongly correlated with vertical velocity (Fig. 9e), with values up to 0.25 %, which are probably overestimated, although this cannot be confirmed as measurements of supersaturation peaks are not available beyond the surface. However, droplet concentration variations are smooth, and do not show a strong correlation with the maximum supersaturation, because of the pre-existing droplets. Near the ground, maximum simulated values of supersaturation lie around 0.1 % while Hammer et al. (2014) and Mazoyer et al. (2016) reported observed supersaturation peaks lower than 0.1 %.
During the dissipation phase, heterogeneities remain at the top of the fog layer, but the signature of KH waves disappears (not shown). The dissipation of fog at ground level takes about 20 min, and, as noted in Bergot et al. (2015), does not reveal a clear effect of surface heterogeneity.
Having characterized vertical and horizontal heterogeneities of the fog during its life cycle, sensitivity tests are now presented to identify the sources of variability and their impact on the fog life cycle.

Sensitivity study
In order to better characterize the physical processes dominating the fog life cycle, sensitivity tests were conducted in a second step. The resulting simulations and their differences relative to the REF simulation are summarized in Table 1.

Impact of trees
To evaluate the impact of trees on the fog life cycle, a simulation called NTR is run, in which the tree barrier was replaced by grass. Deposition on the grass was considered over the whole domain. Figure 3a shows that, without trees, the 10 m wind speed is overestimated over the measurement area. As in REF, but 30 min earlier, the model develops a sudden increase in TKE around 02:30 UTC at the beginning of the development phase. This change is linked to the increase in the optical thickness and not to the turbulence induced by the trees (Figs. 3b and 7b). After this period, TKE is underesti-mated and remains stronger at 10 m height than at 30 m, contrary to observation. This means that the drag effect of trees is responsible for the observed stronger TKE at 30 m height. The fact that the REF simulation develops very similar TKE at 10 m and 30 m a.g.l. probably means that the representation of surface heterogeneities is still underestimated. This can be explained by the broad range of surface covers present in reality, in addition to the trees (lake, small buildings, etc.) but not included in the simulation.
The main differences in dynamics between NTR and REF appear first on total TKE, with a thinner layer of TKE values greater than 0.5 m 2 s −2 and smaller maxima (Fig. 7b). Before the fog formation, the too thin layer of turbulence near the ground in NTR limits the supply of warmer air from above. This induces an overestimation of the vertical temperature gradient before the fog, and emphasizes the cooling in the low levels of 2 K less than in REF (Fig. 2b). Figure 10a presents the temporal evolution of cloud mixing ratio vertical profiles during the NTR simulation, to be compared to Fig. 6a for REF. Figure 11a and b show instantaneous vertical cross sections of potential temperature at the fog formation with REF and NTR. The stronger cooling in NTR homogenizes the fog formation at the ground and prevents elevated fog formation. The consequence is that the onset of fog in NTR occurs almost 2 h earlier than actually observed and than in the REF simulation (Fig. 2d).
During the formation and development phases, the fog layer is thinner in NTR than in REF. This is due to the formation at the ground and the absence of mixing without trees, thus limiting the vertical development. The maximum of cloud mixing ratio in NTR is increased compared to REF, due to the absence of warming by entrainment. It leads to largely overestimated cooling near the ground in comparison to observations (Fig. 2b). Inside the fog layer, despite the increase in r c , the positive temporal evolution of N c , called the production of N c is not higher than in REF (Fig. 10b), as smaller vertical velocities and higher cloud mixing ratio production compensate for the stronger cooling in the activation process. During the development and the mature phases, the LWP is also largely overestimated with NTR (Fig. 12a).
During the development phase, 500 m wavelengths of KH waves are more smooth and regular without trees and this is noted during the whole phase. This is shown on kinetic energy spectra applied to vertical velocity over the whole fog depth, computed according to Ricard et al. (2013) and presented in Fig. 13. The spectra of REF and NTR present two main differences: firstly the TKE variance is smaller with NTR at wavelengths shorter than 200 m. This means that the flow presents fewer fine-scale structures without the tree drag effect. Secondly, the peak of variance at 500 m wavelength, corresponding to the KH waves, is more pronounced in NTR.
To summarize, the absence of the tree barrier produces an unrealistic simulation, as it causes the fog onset to occur too early (almost 2 h in advance). It also induces cooling that is too strong in the low levels, and a large overestimation of the LWP. The absence of trees also modifies the signature of the KH waves at the top of the fog layer, with a more regular pattern and fewer small-scale heterogeneities. The impact of the deposition process will now be examined more precisely.

Impact of deposition
Three simulations were carried out to better characterize the role of the deposition process, all keeping the tree barrier. The first one, called NDT, removed only deposition over trees compared to REF. In this case, trees acted as grass for deposition. This was done by activating deposition only at the first level of the model. The second one, called NDG, removed deposition altogether. The third one, noted DE8, considered a deposition velocity V DEP of 8 cm s −1 over grass and trees, which is the upper bound given by Katata (2014) instead of 2 cm s −1 as in REF.
NDT very slightly increases the LWP during the fog life cycle (Fig. 12a). Conversely, removing deposition everywhere with NDG has a considerable impact. The onset of fog occurs at the surface and not at 30 m height and almost 2 h earlier than in observations and in the REF simulation (Fig. 10c). During the development phase, there is no longer a vertical gradient of r c and N c (Fig. 10c and d).
The fog layer is deeper throughout the life cycle, and therefore the LWP is largely overestimated with a maximum be- tween 05:00 and 06:00 UTC of about twice the observed value (Fig. 12). Due to the larger amount of cloud water near the ground, the dissipation at the ground is delayed by more than 1 h. In contrast, DE8 induces a significant reduction of the LWP, and the onset of fog near the ground coincides relatively well with observation. The formation of fog at elevated levels is more pronounced, and r c over the whole fog depth is reduced during the development phase compared to REF ( Fig. 10d and e). This means that the deposition process is highly sensitive to the deposition velocity. Zhang et al. (2014) have already shown that including a deposition term in simulations seems to have some effect on the droplet concentration in the layer near the ground and consequently on visibility. However, the effect they found was less pronounced than the one seen here. A possible explanation is that both u * , the friction velocity, and the mean volumetric diameter of droplets used in their parameterization were underestimated. In our case, the deposition process, even with a simple parameterization, appears to be essential to correctly simulate the fog life cycle and to approach the observed LWP. Neglecting this process modifies the fog life cycle in terms of onset and dissipation times. The elevated fog formation, which is a climatological characteristic of the SIRTA site, is the result of two effects: the tree drag effect, which mixes the lowest levels, and the deposition process, which erodes the near-surface water content. We will now examine the impact of the horizontal resolution on the simulated fog life cycle.

Sensitivity to effective resolution
In order to assess the impact of spatial resolution on the fog life cycle, a 2 m horizontal resolution simulation (called DX2) was carried out using the same momentum advection scheme as in REF (CEN4TH). According to Skamarock (2004), kinetic energy (KE) spectra deduced from simulations allow the effective resolution to be set up as the scale at which the model starts to depart from the theoretical slope, which is −3 for vertical velocity spectra applied to stable turbulence. Mean KE spectra applied to the vertical wind component reveal an effective resolution of the order of 4-5 x for simulations with CEN4TH (DX2 and REF), in agreement with Ricard et al. (2013), namely 8 and 20 m, respectively (Fig. 13).
In two other tests performed on the wind transport scheme, keeping the 5 m horizontal resolution, the CEN4TH scheme was replaced by the WENO (Weighted Non-Oscillatory; Shu, 1998) scheme at 3rd order (called WE3) or 5th order (called WE5). These spatial schemes, associated with an explicit Runge-Kutta temporal scheme, allow time steps 10 times larger than CEN4TH associated with a leap-frog temporal scheme, but they were run here with the same small time step (0.1 s) for comparison. Due to the upstream spatial discretization, WENO schemes are implicitly diffusive and are therefore characterized by a coarser effective resolution, especially WE3 because of its lower order. Figure 13 shows that the effective resolutions are 35 m (i.e. 7 x) and 70 m (i.e. 14 x) for WE5 and WE3, respectively.
WE3 significantly reduces the top entrainment and the supply of warmer, dryer air from above. This emphasizes the cooling near the surface (Fig. 11c) as the diffusive contribution of the advection operator dissipates small updraughts and suppresses part of the resolved KE variance, in particular that present at the top of the fog layer. This induces an overestimation of the thermal gradient near the surface before the fog, and leads to cooling that is too strong by 1 K during the fog (not shown). The consequences of the increased cooling are that the LWP is largely overestimated throughout the fog life cycle, and the dissipation is delayed (Fig. 12b). Considering the LWP, WE3 tends to be closer to the NTR simulation, meaning that a diffusive transport scheme significantly diminishes the tree drag effect.
In contrast, the differences between WE5 and REF are very small: only the LWP is higher with WE5 during the dissipation phase due to a slightly deeper fog layer. This underlines the less diffusive behaviour of WE5 and its higher accuracy compared to WE3.
Thus the jump in the effective resolution with the diffusive WE3 scheme affects the fog life cycle significantly, while the smaller deviation with WE5 has almost no impact. Increasing numerical implicit diffusion seems to have almost the same effect as removing the drag effect of trees. This also underlines the importance of the numerical schemes for correct handling of the cloud edge problem (Baba and Takahashi, 2013).

Conclusion
Large eddy simulations of a radiation fog event observed during the ParisFog campaign were performed, with the aim of studying the impact of dynamics on the fog life cycle. In order to study the local structures of the fog depth, simulations were performed at 5 m resolution on the horizontal scale and 1 m on the vertical scale near the ground, and included a tree barrier present near the measurement site, taken into account in the model by means of a drag approach. The model included a 2-moment microphysical scheme, and a deposition term was added to the droplet sedimentation, representing the interception of droplets by the plant canopies and acting only at the first vertical level above grass, and above the height of the trees.
The performance of the reference simulation was satisfactory as it gave fairly good agreement with the classical near-surface measurements and the LWP. This good performance allowed the processes driving the fog life cycle to be explored.
The formation of the fog at elevated levels and the fact that it subsided to the ground in a very short time, a frequently observed characteristic of radiation fog events at the SIRTA site, has been explained. It is a consequence of the tree drag effect when the wind meets this obstacle and the deposition effect, which reduces the formation of droplets near the surface. In contrast, the fog formed at the surface first upstream and 500 m downstream of the trees, leading to a duration of about 1 h for fog formation at the surface over the whole domain.
At the beginning of the development phase, the fog became optically thick to LW radiation, inducing a significant increase in KE by dynamical production, which was also associated with temperature convergence at low levels. The radiative cooling near the top of the fog layer was the main source of droplet activation so the droplet concentration was maximum in the upper levels of the cloud.
During the development phase, the fog layer depth grew more slowly when the fog reached the top of the nocturnal boundary layer, encountering stronger thermodynamical gradients and wind shear. Horizontal rolls at the top of the fog layer, associated with KH instabilities, became prominent. The cloud droplet concentration became quasi-homogeneous in the fog layer when averaged over time but extremes of droplet concentration occurred locally near the top of the fog in the radiative cooling layer, with maxima preferentially upstream of the crests of the waves rather than downstream in the ascent area. This indicates that vertical velocity makes up the main contribution to droplet activation at the top of the fog layer, followed by the contribution of radiative cooling. Inside the cloud layer, maxima of supersaturation were directly linked to the local updraughts, while variations in droplet concentration were smoother.
During the dissipation phase, as the fog evolved into a stratus layer, the cloud mixing ratio decreased at all levels. However, a sharp increase in the droplet concentration occurred over the whole depth of the cloud because droplets were now only activated by the convective ascents.
Various sensitivity tests allowed the main processes affecting the evolution of fog to be identified. The tree drag effect and the deposition process were considered as essential to correctly reproduce the main characteristics of the fog. The absence of the tree barrier produced an unrealistic fog simulation, with too early an onset, excessively strong cooling and a large overestimation of the LWP.
Neglecting the deposition process over the whole vegetation canopy exerted the most significant impact on the fog prediction. It overestimated LWP, prevented elevated fog formation, modified the fog life cycle and suppressed vertical and temporal heterogeneities of the microphysical fields. Conversely, increasing the droplet deposition velocity from 2 to 8 cm s −1 reduced the LWP.
Increasing the horizontal resolution to 2 m did not change the fog prediction significantly, which means that grid convergence seems to be achieved at these resolutions. Conversely, increasing the numerical diffusion with a momentum transport scheme of lower order, involving a coarser effective resolution, drastically limited the top entrainment, and tended strongly towards the solution where the tree drag effect was ignored. This underlined the importance of the properties of numerical schemes in LES, particularly at cloud edges.
This study demonstrates the feasibility and the interest of LES including surface heterogeneities to improve our understanding of fog processes. At these fine resolutions, surface heterogeneities have a strong impact, explaining part of the variability in the fog layer and making these simulations very challenging. Therefore, horizontal and vertical variabilities in the fog layer also need to be more thoroughly explored in future field experiments. The horizontal variability, especially at the onset of the fog, also stresses that one point observation may not be very representative of what happens over a coarser grid box of a numerical weather prediction model.
One of the main points of this study is that fog water deposition should not be neglected in 3-D fog forecast models, as still often occurs. It influences not only near-surface fields but also the whole fog life cycle. In this study, the deposition term was introduced quite crudely and this would need some refinement in further studies. It would need to take account of the wind speed and the turbulence, and it could also consider the hygroscopic nature of canopies. By analogy with dry deposition, it would also be better to take droplet diameter into account. Other studies have also shown that fog water deposition is strongly enhanced at the forest edge, becoming up to 1.5-4 times larger than that in closed forest canopies (Katata, 2014), so it could be interesting to simulate the edge effect of fog water deposition. It is also crucial to perform measurements of fog water deposition and dewfall during field experiments (Price and Clark, 2014).
This study has shown the great importance of some dynamical effects operating at 1st order for correct predictions of the fog life cycle. Microphysics near the ground will be further explored in a future study, and the impact of aerosols on the fog life cycle will be considered.
Data availability. The Meso-NH code in the current version is freely available from the Meso-NH site: http://mesonh.aero. obs-mip.fr/mesonh53. Data relative to the case study simulations are available on request from Christine Lac (christine.lac@meteo.fr).