Ch 4. When pollution meets forest air

OVERVIEW

    • In last chapter we played around with the OH oxidation of isoprene with the presence of NOx and compared the results with a chamber study. This condition is in fact quite relevant in the real world (well, sorta). Plants and trees emit a lot of isoprene, which will under go various processes in the air. On the other hand, we humans emit a lot of pollutants into the atmosphere (i.e. NOx, CO, etc). We see from the last chapter that isoprene oxidation pathways could be quite different when you have different NOx levels. How would pollution affect the photo-oxidation of biogenic VOCs? This is of great interest.
    • In addition to chemical reactions, a lot of other processes occur in the atmosphere too. For example, emissions and depositions. What's even more interesting is that, in the ambient, many parameters will have diurnal variations! In this chapter we will try to mimic some of these processes in our very simple model.

1. From chamber to ambient

  • Yupp, we already set up the model in the last chapter, let's grab it!
  • We already know how to model chamber (well, sort of). Imagine you got a gigantic chamber that covers a fairly large forest, which emits isoprene. And then you injected a bunch of pollutants into the chamber...

2. Mixing Layer & emissions & depositions

  • This is an important concept but we will not talk too much about it. From the box model point-of-view: we put all chemical species in a box and, we assume everything is well mixed inside this box. In real world this is never the case; but under certain circumstances, this is not the worst assumption.
  • But how big is the box? Well, in our previous examples, this doesn't matter (because the model operates on the volume basis). But if we want to mimic the emission and deposition processes, this becomes important. Emissions and depositions are often expressed as fluxes (in molec/cm2/s), i.e. the number of molecules crossing the interface (ground) and enter / leaving the box. If we assume everything is well mixed, then this amount of molecules crossing the interface are evenly distributed vertically, and hence can be translated into volume rates (in molec/cm3/s), which is what the model solver needs.
  • Therefore, for emissions, we have:

EmissionRate (molec/cm3/s) = EmissionFlux (molec/cm2/s) / MixingLayerHeight (cm)

From time to time, you may find emission fluxes measurements for different species, under different conditions.

  • Deposition is treated slightly different, as often times deposition velocity (in cm/s) is reported. Therefore:

DepositionRate (molec/cm3/s) = (DepositionVelocity (cm/s) / MixingLayerHeight (cm)) * Concentration (molec/cm3)

  • The height of a well mixed mixing layer varies: usually lower at night and higher during daytime. Let's assume it's 500 m at night, and reaches 1500 m at local noon. Again this is vaguely representative for a good sunny day in the mid/low-latitude.

3. Diurnal variation

  • Many variables have diurnal variations, such as SZA (and photolysis frequencies), temperature, and mixing layer height we just talked about. For a very rough approximation, let's assume these are all sin-shaped.
  • I wrote this snippet to mimic the diurnal variation shape (don't laugh):

 Function DiurnalProfile(HourOfDay) Variable HourOfDay Variable Time_SunRise = 6, Time_SunSet = 18 Variable OutputValue If (HourOfDay<6 || HourOfDay>18) OutputValue = 0 Else OutputValue = Sin( Pi*(HourOfDay-Time_SunRise)/(Time_SunSet-Time_SunRise) ) Endif Return OutputValue End

This returns a diurnal profile shape on the right: sun rises at 6:00 and sets at 18:00.

4. SZA: diurnal variation

  • Inside the FOR loop, find the section where j-values are defined. Then replace the original SZA line with the following two lines:

 SZA = 90 + (20 - 90) * DiurnalProfile(Mod(i*TimeInterval, 86400) / 86400 * 24) SZA = SZA * Pi / 180 // this converts degree into radius

The first line sets SZA to 90 degrees at night, starting to decrease at 6:00, reaching daily miminum at 12:00 at 20 degrees, then starting to increase and reaching 90 at 18:00. This section should look something like this:

5. Mixing layer height: diurnal variation

  • Inside the FOR loop, after the j-values section (or before, whatever), insert these two lines:

 MBL_m = 500 + (1500-500) * DiurnalProfile(Mod(i*TimeInterval, 86400) / 86400 * 24)

this sets MBL_m at 500 m at night, which starts to increase at 6:00, reaching 1500 m at 12:00, then starts to decrease and reaches 500 m at 18:00. Remember you'll have to declare this before the FOR loop.

6. Dry depositions and emissions

  • This is a bit more complicated but I'll walk you thru. Firstly, before the FOR loop, define a bunch of 1-dimensional waves of the same length, which is the number of species. Depositions and emissions may very A LOT for different species.
// ==================== // Dry depositon & emissions // ==================== Make/O/N=(DimSize(Species,0)) DryDepRate_molec_cm3_s = 0 Make/O/N=(DimSize(Species,0)) DryDepVelocity_cm_s = 0.1 DryDepVelocity_cm_s[339] = 0.5 // HNO3, ref... DryDepVelocity_cm_s[366] = 0.5 // N2O5, ref... Make/O/N=(DimSize(Species,0)) EmissionRate_molec_cm3_s = 0 Make/O/N=(DimSize(Species,0)) EmissionFlux_molec_cm2_s = 0
  • For dry deposition velocities, we set 0.1 cm/s for all species, which is probably not very wrong. There are models to calculate dry deposition velocities, depending on compound and conditions. We will not cover this topic here. But, for HNO3 and N2O5, we know these two are very "sticky" so we set higher values: 0.5 cm/s is probably not the worst guess.
  • At this point, if you wanna set dry deposition velocities to other values, or for other species, you now know how to do it ;)
  • Then, inside the FOR loop, AFTER the MLH_m definition, let's calculate the dry deposition rates and emission rates:
// ============== // Set dry deposition // ============== DryDepRate_molec_cm3_s = DryDepVelocity_cm_s / (MLH_m * 100) * Concentration_Matrix[i][p]
 // ====================== // Set emission for isoprene only // ====================== EmissionFlux_molec_cm2_s[297] = 1.5e+12 * DiurnalProfile(Mod(i*TimeInterval, 86400) / 86400 * 24) EmissionRate_molec_cm3_s = EmissionFlux_molec_cm2_s / (MLH_m * 100)
  • Deposition rate (molec/cm3/s) is basically concentration (molec/cm3) divided by the deposition lifetime (s), which is calculated by dividing deposition velocity (cm/s) by the mixing layer height (converted into cm).
  • Emission rate (molec/cm3/s) is simply flux (molec/cm2/s) divided by the mixing layer height (converted into cm). Here we set isoprene emission only, which has a sin-shape and peaks at 1.5e+12 molec/cm2/s. Again if you wanna use other values, be my guest.
  • Now the most important part: we have set the values and calculated the deposition and emission rates, but in order to incorporate these into the ODE, we need to pass them to ODE_ScratchBoard.ipf. Click ODE_ScratchBoard button on the GUI. Scroll all the way to the end, copy-paste these:
// ===================== // Set deposition and emissions // ===================== Wave DryDepRate_molec_cm3_s, EmissionRate_molec_cm3_s YDOT = YDOT - DryDepRate_molec_cm3_s YDOT = YDOT + EmissionRate_molec_cm3_s

7. Set initial concentration

  • We're almost there!!! Now let's set the initial conditions.
  • Let's set RH at 60%, set O3, NO, and CO to 30 ppb, 10 ppb, and 100 ppb, and run the model!!! You should get something like this.
  • As shown, OH and HO2 both peak at around local noon. HO2 is about two orders of magnitude higher.
  • Isoprene remains zero at night but increases during day time, because of the way we set up the emissions. Its two major oxidation products, MVK and MACR, both increase during daytime too.
  • Initial NO (10 ppb) rapidly quenched by O3, followed by the increase of NO2. In the early morning, photolysis of NO2 produced some NO. O3 increases strongly with the presence of light, but on the first day only. On the 2nd day, O3 formation is minor since we're running out of NOx. Question: where does NOx go?

8. A few notes: spherical chicken in the vacuum

  • This example demonstrates how to set time-dependent variables, also briefly introduced how to set up emissions and depositions. We also get to play around with isoprene oxidation, as well as how this biogenic VOC interacts with NOx which is often in higher levels in polluted air.
  • This model (as is) is only conceptually useful, as the real ambient environment is way more complicated. Many physical and chemical processes are not covered in this chapter.
  • For example, in the evening, the mixing layer shrinks, forming a inversion layer that is essentially de-coupled with the surface, this is often called "residual layer". Chemistry in the residual layer can be quite different from the near-surface. In the early morning, the increase of mixing layer height would lead to the enhanced exchange with this residual layer, this processes is often referred to as "entrainment". To describe this process, one has to either model the boundary layer dynamics as well as the chemistry in the residual layer, or, to prescribe the chemical characteristics in the residual layer (if you have measurements). Technique-wise, it's fairly easy to add this (check the Seinfeld and Pandis textbook).
  • Advection also plays an important role in the real world. Again to do this, a thorough understanding of the "background air" is needed. Equations are out there and it's quite straightforward.
  • Last but not least, people usually would "spin-up" the model for a certain period of time. In this example, a few days should be sufficient.