Interactive comment on “ Development of a plume-ingrid model for industrial point and volume sources : application to power plant and refinery sources in the Paris region ” by Y .

The authors extend the sub-grid plume treatment in the Polyphemus air quality modeling platform to include particulate matter and treat volume sources, and they apply the plume-in-grid model to power plant and refinery sources in the greater Paris region. The topic is a good fit for GMD, the analyses appear sound, and the manuscript is generally well written. I would recommend publication after the comments below are addressed.

due to the instantaneous emission dilution in the source grid cell and subsequent transport and diffusion processes, and incorrect concentrations of secondary pollutants due to chemical transformations involving under-or overestimated concentrations of emitted species (e.g., Seigneur et al., 2006;Karamchandani et al., 2011). The errors can be reduced by coupling a plume or puff model with the Eulerian model to create a multi-scale model, typically referred to as a plumein-grid (PinG) model. PinG modeling has been used for ozone (O 3 ) since the 1980s (Karamchandani et al., 2011). It was later extended to particles . These models were evaluated for O 3 and PM 2.5 (fine particles with an aerodynamic diameter less than 2.5 micrometers). PinG modeling, compared to a standard grid-based Eulerian model simulation, typically leads to a spatial redistribution of O 3 concentrations due to its interaction with nitrogen oxides (NO x ), but a negligible effect on the domain-wide O 3 budget (e.g., Seigneur et al., 1981;Kumar and Russell, 1996;Karamchandani et al., 2002). On the other hand, PinG modeling for fossil fuel-fired power plants led to lower concentrations of secondary sulfate and nitrate particles in a simulation over the eastern United States, because of lower oxidant concentrations in the plumes .
The performance of a PinG model (or its subgrid-scale model) has been evaluated with measured in-plume concentrations.  compared concentrations simulated with a subgrid-scale puff model (SCICHEM) used in a PinG model (CMAQ-APT) with aircraft measurements of power plant plumes. Modeled plume excess concentrations (i.e., the differences between the plume and background concentrations) of NO x and O 3 showed good agreement with measured concentrations about 4 km downwind of the sources. Karamchandani et al. (2006) introduced an approach to calculate subgrid-scale concentrations using a similar PinG model (CMAQ-MADRID-APT) that accounts for uncertainty in the wind direction by calculating the concentration at a receptor site along an arc centered at the source and extending 30 • on each side of the receptor site. The source impacting the site was identified based on the sulfur dioxide (SO 2 )/ NO x plume excess ratio and the backtrajectory. The subgrid-scale concentrations showed variability in SO 2 concentrations within the grid cell of a factor of three and reproduced the plume concentration observed at the receptor site. The Eulerian model (i.e., without PinG treatment) underestimated this observed concentration by a factor of two.
Previous studies using PinG models have focused mostly on the simulation of elevated point source emissions such as fossil fuel-fired power plant and cement plant stacks. Power plants burning coal and fuel oil emit NO x and SO 2 , but emit negligible amounts of volatile organic compounds (VOC). NO x emissions, which principally consist of nitric oxide (NO) (> 90 %) consume oxidants such as O 3 and hydroxyl radical (OH) near the source. Since the combination of NO x and VOC emissions can lead to the formation of oxidants, it is of interest to apply a PinG model to study the impact of refineries where both NO x and VOC emissions can be important. However, VOC emissions result mostly from leaks and should be treated as fugitive emissions distributed spatially over a finite volume of the industrial site. Thus, the PinG model needs to be modified to account for volume sources in addition to point sources.
First, the development of the PinG model for gaseous and aerosol species due to emissions from point and volume sources is presented. Then, simulations are conducted to evaluate the effect of PinG modeling on two major types of industrial sources: a fossil fuel-fired power plant and a refinery. The simulation domain, the episode, the model configuration, and the treatment of sources are described. Model performance is evaluated by comparison to measurements and the impact of the PinG model on air pollutant concentrations due to emissions from the industrial sources is presented.
The Polyphemus Gaussian puff model treats the transport and the dispersion of puffs as well as the gas-phase chemical reactions occurring in the puffs in interaction with the ambient background. This model has been described by Korsakissok and Mallet (2010a). A brief description of the model is presented here to present the context for the addition of aerosols in the Gaussian puff model and the new treatment of volume sources.
Each puff transports all gaseous and particulate chemical species. The concentration distribution of a chemical species in the puff is assumed to be Gaussian and can be written as follows: where σ x and σ y are the Gaussian standard deviations on a horizontal plane, respectively in the wind direction and perpendicular to the wind direction, σ z is the Gaussian standard deviation in the vertical direction, and Q is the mass of a species in the puff. That mass is defined as the species emission rate from the source multiplied by the time interval between two puffs released from that source. In general, the Gaussian standard deviations for a point source, which is a single, identifiable source of air pollutant emissions, are given by empirical formulae: Geosci. Model Dev., 7, 569-585, 2014 www.geosci-model-dev.net/7/569/2014/ where σ x turb , σ y turb and σ z turb are the plume dispersion coefficients due to atmospheric turbulence, and σ y pr and σ z pr represent the added dispersion due to the plume rise. The term d 2 s 4 represents the initial standard deviation characterized by the diameter of the source d s (Korsakissok and Mallet, 2009). An algorithm was added to model the emissions of a volume source, which is a 3-D source of air pollutant emissions. Two algorithms have been used typically to treat volume sources in air pollutant dispersion models: one using a virtual source (US EPA, 1995) and one using initial plume dimensions (US EPA, 1995;Cimorelli et al., 2004). The former algorithm uses a virtual source located upwind of the actual volume source such that the plume dimensions at the location of the actual source correspond to the dimensions of that source. The location of the virtual source varies as a function of wind direction, wind speed, and atmospheric stability. The latter algorithm simply uses initial values of the standard deviations that are consistent with the dimensions of the source. The two approaches are mathematically identical and only differ in their numerical implementation. Here, the algorithm using the initial dimensions of the plume (or puff) is used. Following US EPA (1995), the initial dimensions of the puff for a volume source are given by where x, y and z are, respectively, the length, the width, and the height of the volume source.
The term d 2 s 4 of Eq.
(2) is therefore replaced by σ yi . Additional terms σ xi and σ zi are added to σ x and σ z , respectively. The plume rise is negligible for a volume source because fugitive emissions have little impact on the ambient temperature of the volume source domain. Therefore, we obtain the following formulations for the puff standard deviations of a volume source: Following the modeling approach used in the SCICHEM reactive puff model (Karamchandani et al., 2000), the concentrations of the chemical species in each puff are treated as perturbations ( c) of the background concentrations, where c p is the concentration in the puff and c b is the concentrations modeled by the Eulerian host model, i.e., the background concentration. To calculate the concentrations of the species in the puff, we use the procedure described by Korsakissok and Mallet (2010a): 1. Let a chemical reaction in a puff be A + B → P, where A and B are gas-phase or particulate species. The concentrations of the species A and B in the puff at time t are c p A (t) and c p B (t). These puff concentrations correspond to the background concentration c b added to the perturbation of the puff c.
2. The new concentrations of the species A in the puff due to the chemical reaction is computed for the time step t. The reaction rate of the species A in the puff is ), where k is the rate constant, (1) represents the background chemistry, (2) is the interaction between the puff and the background species, and (3) is the chemistry specific to the puff.
3. The new background concentration due to the chemical reaction is computed separately from the previous equation for the time step t. The reaction rate for species A in the background is 4. The new background concentration is subtracted from the new puff concentration to obtain the new perturbation of the puff. The reaction rate representing the perturbation is then ).
The new perturbation of the puff corresponds to the terms (2) and (3) in Eq. (6).
The original formulation of the Polyphemus PinG model pertained only to gaseous species. Therefore, the model was also modified to treat particulate species. To that end, the following modules were included in the Gaussian puff model to simulate the concentrations of particles in the puff: the CB05 mechanism (Yarwood et al., 2005) for the gas-phase chemistry, the SIREAM model  for the particle size distribution, the ISORROPIA model (Nenes et al., 1998) for inorganic aerosol species, and the H 2 O model (Couvidat et al., 2012) for organic aerosol species. Aqueous-phase chemistry is simulated using a model, which treats SO 2 oxidation and the gas/liquid and ionic equilibria of major watersoluble species .
The aerosol-and aqueous-phase models were implemented in the Gaussian puff model in the same way as they were implemented in the Eulerian host model for the sake of consistency. In the PinG model, the perturbations of the puff concentrations were calculated at each time step first for the gas-phase species and next for the particulate species.

Description of the simulations
The model presented in the previous section was applied to the Paris region (Greater Paris) for a three-week simulation. Two industrial sources are treated explicitly with the PinG representation: a coal-fired power plant and a refinery. The simulation configurations are described below.

Simulation setup
The air quality simulations were carried out to evaluate the impact of the PinG model on the concentrations of O 3 and PM 2.5 . Three modeling domains are used with one-way nesting. The largest domain covers western Europe and part of eastern Europe with a horizontal resolution of 0.5 • × 0.5 • . The first nested domain covers France with a resolution of 0.1 • ×0.1 • and the smallest domain covers Greater Paris with a resolution of 0.02 • × 0.02 • . The vertical resolution consists of 9 levels up to 12 km with a finer resolution near the surface. The altitudes of the vertical upper boundary of the grid cells are 40 m, 120 m, 300 m, 800 m, 1500 m, 2400 m, 3500 m, 6000 m, and 12 000 m. The smallest domain is presented in Fig. 1. The Weather Research and Forecast (WRF) model version 3.3 with the Advanced Research WRF (ARW) dynamics solver was used to simulate the meteorological fields over Greater Paris (Skamarock et al., 2008). The regular latitude-longitude map projection is used for the three simulation domains with one-way nesting. Horizontal grid spacing of the coarse domain is 0.5 • and those of the two nested domains are 0.1666 • and 0.0555 • , respectively. The largest 0.5 • domain covers Europe and the smallest domain covers Greater Paris. The vertical resolution consists of 28 levels up to 100 hPa (about 16 km). The descriptions of the two coarser domains are given by Royer et al. (2011).
Two simulations over Greater Paris were carried out from 4 July to 29 July 2009 with six days of spin-up to initialize the simulation.
1. "Reference" simulation: the two industrial sources are treated in a standard way by the Eulerian model, i.e., their emissions are released into the appropriate grid cells.
2. "Plume-in-Grid (PinG)" simulation: the two industrial sources are treated with the subgrid-scale puff model.
For the PinG simulation, the time interval between puffs must be selected to reproduce the continuous plumes with sufficient accuracy as well as to restrict the computation time.  Korsakissok and Mallet (2010b) investigated the sensitivity of PinG modeling to time intervals between puffs. They compared a simulation conducted using a time interval of 1 s to simulations using greater time intervals. The impact of using different time intervals was negligible for time intervals below 50 s; the error was about 10 %, with a time interval of 200 s. Here, the time interval between puffs was set to 200 s, which corresponds to a third of the time step for the Eulerian model. Korsakissok and Mallet (2010b) tested two criteria for the puff injection into the Eulerian host model: a time criterion and a criterion based on the puff horizontal size (2λ y σ y ) reaching the grid cell horizontal size. λ y is a constant for the effective size of puffs; it is set to 2 in this study (i.e., 95.4 % of the puff material is included) following Seigneur et al. (1981) and Korsakissok and Mallet (2010b).
The sensitivities of modeling results to the puff injection criteria have been reported by Korsakissok and Mallet (2010b). They are summarized here. The use of the time criterion (time elapsed since the emissions from the source), with a 1 h injection time, improved the performance of the PinG model compared to that of the Eulerian model for a domain with a grid cell size less than about 60 km. The size criterion (puff size commensurate with the grid size) gave better results than the time criterion for a domain with a fine grid cell size (< 25 km); however, the results for a domain with coarser grid cell size were poor because of a greater injection time (more than 6 h), which may lead to large errors in puff trajectories. Therefore, the time criterion gave better results than the size criterion for a domain with a coarse grid cell size (> 50 km).
Geosci. Model Dev., 7, 569-585, 2014 www.geosci-model-dev.net/7/569/2014/ The time criterion using 1 h of puff travel time allows one to minimize the number of puffs, which minimizes the computational time. However, the subgrid-scale treatment is then used only near the sources and one may not benefit fully from it. These two injection criteria are used jointly here and the puffs are injected into the corresponding grid cells as soon as one of the two criteria is met. The sensitivity of the PinG model simulation to the puff injection criteria is investigated below, using only the puff size criterion, i.e., using the subgrid-scale treatment over longer distances from the source. The Gaussian puff model formulation uses similarity theory for the parameterization of the Gaussian standard deviations and the integrated injection for the injection method. In the integrated injection method, puff mass is distributed over the grid cells covered by the puff horizontal and vertical extents (here, 2σ y and 2σ z ). The detailed parameterizations are given by Korsakissok and Mallet (2010b).
For the comparison of the concentrations calculated in the Reference and PinG simulations, the perturbations of concentrations in non-injected puffs in the PinG simulation are taken into account by injecting the puffs in the corresponding grid cells at the moment of the comparison (and solely for the purpose of the comparison). The integrated injection method is used for this temporary merging of the puffs into the host model.

Emissions
The European Monitoring and Evaluation Programme (EMEP) anthropogenic emission inventory was used for the domains covering Europe and France. Over Greater Paris, the anthropogenic emissions were generated with the 2005 inventory of Airparif, the air quality agency for the Paris region (http://www.airparif.asso.fr/). The surface emissions and aircraft emissions were spatially distributed over the Eulerian grid. There are 196 industrial sources in the domain. Those point source emissions were located in the corresponding grid cell(s) following calculation of the plume rise. Two industrial sites, the Vitry power plant and the Grandpuits refinery were selected for PinG treatment. Therefore, they are treated in the same manner as the other point sources in the Reference simulation and are treated with the subgrid-scale puff model in the PinG simulation. The locations of these two sites are indicated in Fig. 1.
The plume rise is computed in the simulations based on the characteristics of the sources given by Airparif (exhaust velocity, exhaust temperature, stack height, and stack diameter). The Briggs formulae were used for the calculation of the plume rise, as described by Korsakissok and Mallet (2009).
The temporal variation for the emissions was obtained by applying temporal factors available from Airparif. These factors are computed by source categories, which are defined by the SNAP (Selected Nomenclature for sources of Air Pollution) code. Therefore, the factors are generic for a category, and they are not specific to a given source.
Source-specific temporal profiles were used for the sources treated with the subgrid-scale model. For the power plant, the emissions were obtained from data provided by EDF, the operator of the power plant. Figure 2 presents the temporal profile obtained from EDF for the power plant. According to Total, the operator of the refinery, the refinery runs continuously (except for some maintenance operation periods) and the emissions of the refinery were considered to be constant for the entire period of simulation. These sourcespecific temporal profiles are also applied when these sources are treated in the "Reference" simulation.
Emissions of air pollutants from the oil refinery were assumed to be released into the atmosphere from the stack except VOC emissions, which were assumed to occur as fugitive emissions over the entire volume of the refinery, e.g., from tanks, pumps, compressors and valves (Bénassy et al., 2008). Therefore, the VOC emissions from the refinery are treated with a volume source and other emissions (SO 2 , NO x and PM) are treated with a point source. For the Grandpuits refinery, the dimensions of the VOC volume source given by Total are 1300 m (width) × 940 m (length) × 20 m (height).
Annual emission rates are provided in the Airparif inventory for each species and each point source (tons year −1 ). For all other sources, temporal factors (monthly, weekly and hourly) were applied to obtain a hourly factor, which represents the ratio of the emission rate for a given hour to the annual emission rate. The annual emission rates for the Vitry power plant and the Grandpuits refinery are presented in Table 1. The NO x emission rate of the Vitry power plant is about three times greater than that of the Grandpuits refinery, while the latter has a greater VOC emission rate than the former. Biogenic emissions were computed with the Model of Emissions of Gases and Aerosols from Nature (MEGAN) (Guenther et al., 2006).

Results
Figure 3 presents surface concentrations of O 3 , NO x , SO 2 and PM 2.5 averaged over the simulation period (4-29 July 2009) obtained with the Reference simulation for the Greater Paris domain. The results of this simulation are first compared to available measurements. Then, the Reference and PinG simulations are compared and differences due to the PinG treatment are discussed.

BDQA network
The model evaluation is performed using the hourly concentrations observed at the stations of the BDQA (Base de Données de la Qualité de l'Air) monitoring network. The locations of the BDQA stations, which are in the Greater Paris domain, are shown in Fig. 1. Table 2 shows the model performance statistics over the entire period and all stations for the concentrations of PM 10 (particles with an aerodynamic diameter less than 10 micrometers), PM 2.5 , O 3 , and NO 2 (nitrogen dioxide). The results are presented only for the PinG simulation because the differences in the performance statistics between the Reference simulation and the PinG simulation are negligible. This close similarity between the performance statistics of the two simulations results from the fact that the plumes from the two industrial sources impact only a few monitoring stations for short periods and have small impacts on region-wide pollutant concentrations. Similar conclusions were obtained in previous PinG modeling studies (e.g., Karamchandani et al., 2002;Korsakissok and Mallet, 2010a). Differences appear when focusing on the impacts of those specific sources, as discussed below.

Megapoli campaign
During the Megapoli campaign of July 2009 (http://www. pole-ether.fr/megapoli/index.jsp), concentrations of various species were measured that are not included in the BDQA measurements. Therefore, model performance is also evaluated here using the Megapoli campaign measurements. Two sites are available for the evaluation of hourly concentrations: SIRTA (Site Instrumental de Recherche par Télédétection Atmosphérique, Haeffelin et al., 2005) and LHVP (Laboratoire d'Hygiène de la Ville de Paris). Model performance statistics are presented in Table 3.
Model performance statistics for O 3 and NO 2 at SIRTA are similar to those at the BDQA stations and meet the performance criteria for O 3 . Gas-phase concentrations of isoprene are evaluated at SIRTA. The model reproduces well the temporal variability in the concentrations of isoprene (correlation coefficient: 0.58); however, it significantly underestimates the observed values by about a factor of three. This discrepancy can be explained by uncertainties in isoprene Geosci. Model Dev., 7, 569-585, 2014 www.geosci-model-dev.net/7/569/2014/ emissions, which are estimated to be up to a factor of three for specific times and locations when different input variables are used in the emission calculations (Guenther et al., 2006). For elementary carbon (EC) in PM 2.5 (EC 2.5 ), the correlation is better at SIRTA (correlation coefficient: 0.49) than at LHVP (correlation coefficient: 0.27). This result is consistent with a previous model evaluation for EC (Couvidat et al., 2013). However, biases and errors are greater at SIRTA than at LHVP. For organic carbon (OC) in PM 2.5 (OC 2.5 ) at LHVP, model performance statistics are similar to those for EC 2.5 at the station except for greater biases for OC 2.5 . The performance goals are met for EC 2.5 at LHVP but are slightly exceeded at SIRTA and for OC 2.5 at LHVP because of an underestimation greater than 30 %. For inorganic aerosols, sulfate is underestimated and nitrate is overestimated at SIRTA. The underestimation in the concentrations of sulfate may contribute to the overestimation in the concentrations of nitrate because of interrelated equilibria of ammonium sulfate and ammonium nitrate. This result is typical of current air quality model simulations (Solazzo et al., 2012). The concentration of ammonium is better estimated than those of sulfate and nitrate. The performance goal is met by ammonium and almost met by sulfate (a bias of −31 % versus a goal of −30 %), but is not met by nitrate. Although concentration differences between the two simulations are very small at the measurement sites and over most of the entire domain, local differences can be important near the industrial sources, which are treated in the subgrid-scale model of the PinG simulation, as discussed below.

Impacts of PinG modeling on pollutant concentrations
The spatial impacts of PinG modeling on the surface concentrations of major species are discussed below. Figure 4a and b present the differences in the surface concentrations of NO x and O 3 between the PinG and Reference simulations, respectively. The differences are averaged over the entire period of simulation (4-29 July 2009). Using the PinG model tends to decrease the concentrations of NO x near the location of the point sources. In the PinG model, NO x is not immediately dispersed in the grid cells where the sources are situated and it can remain in a plume aloft and be transported farther downwind. However, under conditions where the plumes are rapidly mixed to the ground (here, mostly with westerly winds), the PinG simulation leads to greater surface concentrations near the sources because the emitted material is less dispersed. The concentrations of NO x in the grid cells farther downwind of the sources are higher in the PinG simulation than in the Reference simulation because the plume material is released in one or a few grid cells in the PinG simulation, whereas the corresponding emitted material has been dispersed among more grid cells in the Reference simulation.
The surface concentrations of NO x are smaller with PinG near the Vitry grid cell (−1 %); however, they are generally higher with PinG near the Grandpuits grid cell (11 %). These differences can be interpreted by the different characteristics of the sources (see Table 4). The stack height is higher at Vitry (160 m) than at Grandpuits (81.5 m). Thus, the NO x plumes are emitted at a higher altitude at Vitry. Furthermore, the plume rise is higher at Vitry due to a higher exhaust velocity, although it is partially compensated for by a higher exhaust temperature at Grandpuits. In the simulations, the plumes at the Vitry power plant are emitted at altitudes ranging from 300 to 800 m, whereas the plumes at the Grandpuits refinery are emitted at altitudes ranging from 120 to 300 m. Consequently, the NO x plumes touch the ground closer to the source at Grandpuits than at Vitry. Figure 5a and b presents the vertical cross section of the differences in the NO x concentrations at Grandpuits and Vitry, respectively. As discussed above, the NO x plumes are emitted at higher altitudes at Vitry than at Grandpuits, which leads to a plume touchdown farther downwind at Vitry. Therefore, the impact of PinG modeling on the surface concentrations of NO x is greater at Grandpuits than at Vitry. Note that the plumes are transported east of the sources because the prevailing winds are from the west. The greater concentrations in the Reference simulation at the longitude of the source correspond to upwind diffusion.
The differences in the NO x concentrations between the two simulations explain the differences in the O 3 concentrations. The O 3 concentrations near the industrial sources are typically higher with the PinG simulation because the titra- the Greater Paris region (Kim et al., 2009). However, the concentrations of O 3 northeast of Grandpuits are higher with the PinG simulation where the concentrations of NO x are also higher. It is due to higher surface VOC emissions in the PinG simulation. As mentioned above, the fugitive VOC emissions at the Grandpuits refinery are treated with a volume source, i.e., VOC are emitted near the surface. The higher concentrations of VOC in the PinG simulation lead to higher concentrations of oxidants, which result in a higher ozone formation (see Fig. 4c). Nevertheless, the impact of PinG on the O 3 concentrations is small (< 1 %). Figure 4d and e presents the differences in the surface concentrations of SO 2 and sulfate between the PinG and Reference simulations, respectively. The spatial impact of PinG on the concentrations of SO 2 shows a tendency similar to that of the NO x concentrations. The SO 2 concentrations are smaller with PinG near the Vitry grid cell (−5 %), while they are higher with PinG near the Grandpuits grid cell (50 %). The spatial impact of PinG modeling on the difference in the sulfate concentrations is similar to that of SO 2 . Sulfate is formed by oxidation of SO 2 in the gas and aqueous phases. In the gas phase, the oxidation of SO 2 by OH produces sulfuric acid, which condenses to form particulate sulfate. In the aqueous phase, SO 2 is oxidized by O 3 , hydrogen peroxide (H 2 O 2 ), and oxygen (O 2 ) to form sulfuric acid, which leads to particulate sulfate when the droplets evaporate. The sulfate concentrations in the PinG simulation increase (decrease) with the increase (decrease) in SO 2 at Grandpuits (Vitry) (see Fig. 4d). Consequently, the impact of SO 2 on the formation of sulfuric acid is more important than that of the oxidants. Figure 4f and g presents the differences in the surface concentrations of nitrate and ammonium between the PinG and Reference simulations, respectively. The impact of PinG on the sulfate concentrations explains the differences in the ammonium concentrations between the PinG and Reference simulations (−2 % to 6 %), which are due to the neutralization of sulfate by ammonium. In the case where the sulfate concentrations increase (decrease) with PinG, the ammonium concentrations increase (decrease) according to the concentration of ammonia available. In the case where the sulfate concentrations decrease with PinG, the associated decrease in the ammonium concentrations leads to an increase in gas-phase ammonia concentrations, which results in an increase in ammonium nitrate formation and, therefore, in particulate nitrate concentrations. This phenomenon appears clearly east of the Vitry and Grandpuits sources where sulfate concentrations are lower (higher) in the PinG simulation compared to the Reference simulation, but nitrate concentrations are greater (lower).

Impacts of the volume source treatment in the PinG model
In the previous section, the fugitive VOC emissions at the Grandpuits refinery are treated with a volume source. For further analysis of impacts of the volume source treatment on the VOC emissions, an additional simulation was conducted in which the VOC emissions at the refinery were injected into grid cells as in the Reference simulation, but other pollutants (NO x and SO 2 ) were emitted through the stack using the point source subgrid-scale treatment. The results of this simulation (hereafter PinG-NonVolumeSource) are compared to those of the PinG simulation. Figure 6a presents the differences in the concentrations of a VOC, namely toluene. As expected, toluene concentrations are greater. Differences in oxidant concentrations are shown in Fig. 6b and c. The volume source PinG treatment leads to greater oxidant concentrations downwind of the source due to greater VOC concentrations from the refinery fugitive emissions. However, the differences are small (< 1 %) because VOC refinery emissions are a small contribution to total VOC ambient concentrations. Differences in secondary organic aerosol (SOA) concentrations between the PinG and Geosci. Model Dev., 7, 569-585, 2014 www.geosci-model-dev.net/7/569/2014/ PinG-NonVolumeSource simulations are shown in Figure 6d. SOA is formed in the atmosphere from semi-volatile organic compounds (SVOC), which are formed by oxidation of VOC emitted into the atmosphere from both anthropogenic and biogenic sources. The impact of the volume source in PinG modeling on the concentrations of VOC can lead to the differences obtained between the two simulations near the Grandpuits refinery (3 %). Although the differences in the concentrations of VOC are significant only near the Grandpuits refinery (Fig. 6a), the differences in the concentrations of SOA extend farther downwind of the source. The differences near the Grandpuits refinery are partly due to the differences in the concentrations of anthropogenic SOA formed by oxidation of VOC emitted from the sources (toluene, xylene and glyoxal oligomer). However, the impact of the PinG modeling on the anthropogenic SOA is negligible downwind of the source (not shown). In addition, the contribution of the anthropogenic SOA to the total differences in the concentrations of SOA is rather low (about 10 % of the total differences) even near the source. In fact, the differences in the SOA concentrations are due to differences in the concentrations of biogenic SOA, due mostly to the oxidation of monoterpenes. Since monoterpenes are not emitted from the Grandpuits refinery, the differences in the concentrations of biogenic SOA are due to oxidation of monoterpenes emitted in the surrounding region by oxidants (OH, O 3 and NO 3 ), which differ in concentration between the two simulations (see Fig. 6b and c). However, SOA concentration differences are small (about 0.05 µg m −3 ).

Sensitivity to the criterion for puff injection
A supplementary simulation was conducted using only the size criterion to investigate the influence of the time criterion on PinG modeling. In this simulation (hereafter PinG-SizeOnly), the subgrid-scale treatment was extended over larger distances from the source. The results of the PinG simulation conducted above with the two puff injection criteria applied jointly are compared to those of this PinG-SizeOnly simulation. Figure 7 presents the difference in the surface concentrations of NO x between the PinG and the PinG-SizeOnly simulations. The differences can be depicted by defining three spatial zones. The NO x concentrations are smaller in the PinG simulation than in the PinG-SizeOnly simulation near the sources up to a distance of about 15 km (zone 1). Beyond 15 km downwind of the sources, the NO x concentrations become greater in the PinG simulation (zone 2), while they become smaller again in the PinG simulation farther downwind (zone 3). The boundary between zones 1 and 2 is clear for the different plumes. However, the boundary between zones 2 and 3 is not as clear because it depends on the meteorological conditions influencing plume transport. The nearest boundary between zones 2 and 3 is situated at about 50 km from the sources, but it can be as far as the boundary of the simulation domain.
In zone 1, the NO x concentrations are smaller in the PinG simulation because these situations correspond to low wind speeds, where the emitted pollutants are not transported and dispersed rapidly from the sources. The emitted species staying near the sources are injected earlier into grid cells in the PinG simulation than in the PinG-SizeOnly simulation. The injected species can then be more rapidly dispersed to neighboring grid cells in the Eulerian model. Figure 8a presents the differences in the concentrations of NO x between the PinG and the PinG-SizeOnly simulations at 12:00 UTC on 4 July. Wind speeds for grid cells at an altitude of 210 m where the plumes at the Grandpuits refinery are emitted are less than 2 m s −1 . The differences are only important near the source. The concentrations of NO x are lower in the PinG simulation in grid cells where the plumes are located and they are higher in grid cells that are outside of the plume path because the dispersion following the earlier injection into the grid cells increases the NO x concentration. Nevertheless, the mean difference for the whole period is small in this zone (less than 1 %).
In zone 2, the NO x concentrations are higher in the PinG simulation. Under higher wind speeds, this zone corresponds to the distance from the source where most of the plumes are injected into the Eulerian model after one hour of the subgrid-scale treatment. NO x injected into the Eulerian model in the PinG simulation reaches the ground faster than NO x that is transported by the plumes longer in the PinG-SizeOnly simulation. The plumes in the PinG-SizeOnly simulation are treated by the subgrid-scale model until the size criterion is met. Figure 8b presents the difference in the NO x concentrations between the PinG and the PinG-SizeOnly simulations at 21:00 UTC on 4 July. At that moment, the wind speeds at an altitude of 210 m are greater than 5 m s −1 in grid cells near the Grandpuits refinery. In this figure, there is no significant difference near the source, i.e., in zone 1. Higher concentrations are simulated in the PinG simulation at distances between 15 and 30 km from the source.
In zone 3, the NO x concentrations are smaller in the PinG simulation. NO x injected into the Eulerian model in zone 2 in the PinG simulation is more easily dispersed to neighboring grid cells. On the other hand, NO x transported by plumes to this zone in the PinG-SizeOnly simulation remains more concentrated. In Fig. 8b, lower concentrations are simulated in the PinG simulation at distances between 30 and 45 km. In Fig. 7, mean concentrations of NO x between the PinG and the PinG-SizeOnly simulations differ most in this zone (up to 3 %). The differences in O 3 concentrations are negligible over the whole domain (less than 1 %). However, the differences in sulfate concentrations are important in zones 2 and 3; the sulfate concentrations are greater in the PinG simulation in zone 2 (up to 4 %) and are smaller in zone 3 (up to 5 %). These differences are due to the differences in the concentrations of SO 2 , as discussed above. The differences in SO 2 show similar tendencies to those in NO x . However, the maximum difference in SO 2 (15 %) is higher than that of NO x (3 %) because of lower background concentrations for SO 2 than for NO x . Figure 9 presents the relationship between the differences in the concentrations of NO x and the differences in the concentrations of SO 2 , between the PinG and PinG-SizeOnly simulations. The correlation coefficient is 0.86, reflecting the similar behavior of NO x and SO 2 .  Using least-square fitting, the differences in the concentrations of SO 2 are on average 1.67 times greater than the differences in the concentrations of NO x , because of different emission rates and oxidation kinetics. An indicator, I , which quantifies the influence of the time criterion on the pollutant concentrations, was calculated using the following equation: where C is the concentration in grid cells. Very low C values lead to high I values and can mislead the analysis. Therefore, a threshold value for (C PinG − C Reference ) is applied (0.1 ppb or µg m −3 ). As expected, the indicator I shows similar patterns as those displayed in Fig. 7 with negative values near the source (zone 1 typical of low wind speed conditions) and positive values farther downwind (zone 2 typical of high wind speed conditions), as depicted in Fig. 10 for NO x . The magnitude of I is a measure of the effect of the injection criterion compared to the effect of the use of a PinG treatment. It appears that the latter effect is dominant, as the former is typically 20 % or less of the latter as shown in Fig. 10; similar results were obtained for other species (not shown). Nevertheless, the effect of the injection criterion is not negligible and the use of the size criterion is recommended. However, it should be considered that trajectory uncertainties become large with long injection times, and an upper limit of 3 h for the injection time has been recommended for simulations using very coarse grids (Korsakissok and Mallet, 2010b).

Sensitivity to the horizontal grid size
The impact of PinG modeling on pollutant concentrations can be influenced by the choice of the horizontal grid size of the Eulerian host model. For example, the size criterion for the puff injection into the host model is more rapidly met for simulations with a fine resolution because the puff size reaches the grid size faster. Therefore, a sensitivity study of the impact of the grid size on PinG modeling results was conducted. The two simulations (Reference and PinG-SizeOnly) were repeated using the same modeling setup over Greater Paris, except for a coarser horizontal grid size of 0.10 • instead of 0.02 • . In this study, the resolution of 0.10 • corresponds to 7.3 km (east to west domain mean value) × 11.1 km (south to north) and that of 0.02 • corresponds to 1.5 km (east to west domain mean value) × 2.2 km (south to north).
Mean wind speeds at the two vertical levels where the plumes are emitted (210 m and 550 m) are 7 m s −1 and 7.7 m s −1 , respectively. Therefore, using the time criterion of one hour would result in most of the puffs being injected into neighboring second or third grid cells from the sources. To avoid such early puff injections, the time criterion for the puff injection was not taken into account and only the puff size criterion was used. Indeed, using solely the size criterion Geosci. Model Dev., 7, 569-585, 2014 www.geosci-model-dev.net/7/569/2014/ in the PinG modeling has been advised for horizontal grid sizes less than 25 km (Korsakissok and Mallet, 2010b). For consistency in the injection criteria, the results of the simulations for the injection criteria presented in Section 4.4 are used when the sensitivity simulations for the horizontal grid resolution are compared to the simulations with the fine grid resolution. Figure 11a presents the differences between the Reference and PinG simulations for the NO x concentrations. The NO x concentrations are smaller in the PinG simulation at Vitry and are greater in the PinG simulation at Grandpuits. These results are consistent with those of the simulation with the fine grid size. However, the maximum difference between the Reference and PinG simulations is about 1 % near the Grandpuits grid cell with the coarse grid size. It is much smaller than that obtained with the fine grid size (11 %). This result seems counterintuitive, but is due to the fact that the comparison is made with grid-averaged concentrations. The grid cell volumes are much greater in the simulations with the coarse grid size, which leads to lower pollutant concentrations with a coarse grid size in the Eulerian host model, and the impact of the PinG modeling is reduced with the coarse grid size. This result is consistent with previous studies. Using coarse grid sizes in Eulerian models leads to lower maximum NO x and O 3 concentrations because of more diluted NO x emissions (Cohan et al., 2006;Henderson et al., 2010). Also, the impact of PinG modeling on the concentrations of a passive tracer increases with finer grid resolutions (Korsakissok and Mallet, 2010b).
Two test simulations were conducted to compare characteristics of puff traveling with the fine and coarse grid size during 24 h. Puffs in the simulation with the coarse grid size can travel longer than in the simulation with the fine grid size. This longer traveling time of puffs can reinforce impacts of PinG modeling in the simulation with the coarse grid size. The traveling times in the simulations are compared in Table 5. The mean puff traveling time with the coarse grid size is twice as great as that with the fine grid size. About 80 % of puffs in the simulation with the coarse grid size could reach the boundaries of the simulation domain without being injected into the Eulerian grid. However, 70 % of puffs in the simulation with the fine grid size were injected into the host model due to the smaller grid size. Figure 11b presents the differences between the Reference and PinG simulations in O 3 concentrations. The lower O 3 concentrations in the PinG simulation at the Grandpuits refinery are consistent with those obtained in the simulation with the fine grid size and result from the titration of O 3 by NO. The O 3 concentrations farther downwind east of the Grandpuits refinery are higher in the PinG simulation in grid cells where the NO x concentrations are also higher. This result was explained above by the presence of fugitive VOC emissions, which are emitted near the surface as a volume source. However, the impact of the VOC emissions is less  Figure 11c and d presents the differences in the surface concentrations of SO 2 and sulfate between the PinG and Reference simulations, respectively. For SO 2 , the impact of the PinG model shows a tendency similar to that of the NO x concentrations. The spatial impact of PinG modeling on the difference in sulfate concentrations is similar to that of SO 2 . However, the sulfate concentrations are lower in the PinG simulation at the Grandpuits grid cell, where the SO 2 concentrations are higher in the PinG simulation. It can be explained by lower oxidant concentrations (see Fig. 11b), which lead to lower SO 2 oxidation to sulfate. Karamchandani et al. (2006) reported lower sulfate formation when using PinG modeling. In the simulations presented here, however, the effect of the oxidation kinetics is limited and surface concentration differences tend to follow SO 2 concentration differences in most parts of the domain. It is likely that the NO x emissions are not sufficiently high to alter oxidant fields over long distances because of other significant sources of NO x in the Greater Paris region (e.g., on-road mobile sources). Maximum differences in the sulfate concentrations between the PinG and the Reference simulations with the coarse grid size is lower (3 %) than that with the fine grid size (15 %).

Conclusions
The Polyphemus plume-in-grid (PinG) model was modified to include a full treatment for PM. Furthermore, a PinG treatment for volume sources was developed and a PinG simulation was conducted for the first time for industrial sources of SO 2 , NO x , and VOC.
The impact of PinG modeling over Greater Paris for two industrial sources, a coal-fired power plant and a refinery, was studied. The results show the importance of source characteristics (stack height and diameter, exhaust temperature and velocity) for the surface concentrations of primary pollutants emitted aloft (e.g., NO x , SO 2 , and primary PM) simulated with PinG. The impact of PinG on NO x leads to an impact on the formation of O 3 due to these sources; however, this impact is weak on average because the simulation domain is VOC-limited. The formation of sulfate in the subgrid-scale model is mostly influenced by the different dispersion of SO 2 in the PinG and Reference simulations and little affected by oxidant concentrations. This result suggests that the relative influence of the precursor (here SO 2 ) and oxidant concentrations in PinG modeling is sensitive to the NO x and VOC levels in the plume and the relative importance of the source treated at the subgrid scale compared to other sources (e.g., traffic) within the domains. The impact of PinG on nitrate and ammonium concentrations results from the interrelated equilibria of ammonium sulfate and ammonium nitrate. The PinG treatment of VOC fugitive emissions, using a volume source at the oil refinery, plays an important role and leads to higher surface concentrations of VOC, which subsequently lead to slightly greater O 3 concentrations. The impact of PinG on the VOC concentrations leads to negligible differences in the concentrations of anthropogenic SOA. Concentrations of biogenic SOA are also influenced by PinG modeling because of differences in oxidant concentrations; however, these differences remain small. The results obtained here reflect the impact of a single fugitive volume source of VOC in a polluted area (Greater Paris). The application of a PinG volume source approach to a heavily industrialized area with multiple VOC fugitive sources should be investigated, as impacts of the PinG approach may then be more significant.
The impact of the time criterion for puff injection into the host model is negligible for O 3 (less than 1 %). Nevertheless, the use of the size criterion is recommended for PM 2.5 although an upper limit of 3 h could be imposed for puff travel times to minimize the effect of uncertainties in puff trajectories. However, it must be noted that the injection criterion has less effect on the air quality simulation than that of the plume-in-grid treatment since the effect of the former on pollutant surface concentrations is typically 20 % or less of the effect of the latter.
When using a coarse horizontal grid size, the impacts of PinG modeling are lower than those with the fine grid size because all concentrations are more diluted in the host Eulerian model with the coarse grid resolution.
The PinG modeling results presented here demonstrate that fugitive emissions need to be taken into account in addition to stack emissions for industrial sites treated at the subgrid scale. The effect of PinG modeling on secondary pollutants is complex and depends strongly on the relative importance of the sources treated at the subgrid scale compared to other sources within the domain.