chem/module_mosaic_gly.F

   1       !----------------------------------------------------------------
   2       ! SOA formation from glyoxal with complex formulation including
   3       ! * reversible formation (Kampf et al., ES&T, submitted)
   4       ! * dark/ammonium-catalyzed formation (Noziere, J. Phys. Chem., 2009)
   5       ! * OH chemistry (Ervens and Volkamer, ACP, 2010)
   6       ! * surface uptake (Ervens and Volkamer, ACP, 2010)
   7       ! Christoph Knote, ACD, NCAR, 20130326
   8       !----------------------------------------------------------------
   9
  10       MODULE module_mosaic_gly
  11
  12       IMPLICIT NONE
  13
  14       INTEGER, PARAMETER :: nspecs   = 13,                             &
  15                             igly_g   =  1,                             &
  16                             igly_r1  =  2,                             &
  17                             igly_r2  =  3,                             &
  18                             igly_nh4 =  4,                             &
  19                             igly_sfc =  5,                             &
  20                             igly_oh  =  6,                             &
  21                             ic_as    =  7,                             &
  22                             ic_an    =  8,                             &
  23                             ia_nh4   =  9,                             &
  24                             ioh_g    = 10,                             &
  25                             iph      = 11,                             &
  26                             iwater   = 12,                             &
  27                             iarea    = 13
  28
  29       ! >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  30       LOGICAL, PARAMETER :: lfast_tau1 = .FALSE.
  31       ! <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  32
  33       CONTAINS
  34
  35       ! adapted from Numerical Recipes, Second Edition, p. 706-707
  36       ! changes: added "h" (timestep) argument to DERIVS calls
  37       SUBROUTINE rk4(y, dydx, n, x, h, yout, derivs)
  38
  39       INTEGER n
  40       REAL(kind=8) :: h, x, dydx(n), y(n), yout(n)
  41       EXTERNAL derivs
  42       INTEGER i
  43       REAL(kind=8) :: h6, hh, xh, dym(nspecs), dyt(nspecs), yt(nspecs)
  44
  45       hh=h*0.5
  46       h6=h/6.
  47       xh=x+hh
  48
  49       DO i=1, n
  50         yt(i) = y(i) + hh * dydx(i)
  51       ENDDO
  52
  53       CALL derivs(xh, yt, h, dyt)
  54
  55       DO i=1, n
  56         yt(i) = y(i) + hh * dyt(i)
  57       ENDDO
  58
  59       CALL derivs(xh, yt, h, dym)
  60
  61       DO i=1, n
  62         yt(i) = y(i) + h * dym(i)
  63         dym(i) = dyt(i) + dym(i)
  64       ENDDO
  65
  66       CALL derivs(x+h, yt, h, dyt)
  67
  68       DO i=1, n
  69         yout(i) = y(i) + h6 * ( dydx(i) + dyt(i) + 2. * dym(i))
  70       ENDDO
  71
  72       RETURN
  73
  74       END SUBROUTINE rk4
  75
  76       ! Simple SOA formation from glyoxal as presented in
  77       ! Washenfelder et al, JGR, 2011
  78       SUBROUTINE glysoa_simple(dtchem)
  79
  80       USE module_data_mosaic_therm, ONLY: t_k, area_wet_a, gas, aer,   &
  81                                           jtotal, igly, iglysoa_sfc_a, nbin_a
  82
  83       IMPLICIT NONE
  84
  85       REAL(kind=8), INTENT(IN)  :: dtchem
  86       REAL(kind=8)              :: omega, gamma_gly, A, delta_gly, frac_A
  87       INTEGER                   :: ibin
  88
  89       ! mean molecular velocity of glyoxal (cm/s)
  90       omega = 1.455e4 * sqrt(t_k / 58.0_8)
  91
  92       ! aerosol uptake coefficient for glyoxal (-)
  93       ! Washenfelder et al., 2011:
  94       ! 0 - 8   x 10^-4
  95       !     2   x 10^-4 (+/- 1 x 10^-4)
  96       ! Volkamer et al., 2007
  97       !     3.7 x 10^-3
  98       ! B. Ervens, pers. comm., 2010
  99       gamma_gly = 3.3E-3
 100
 101       ! get total aerosol surface area (cm^2 / cm^3)
 102       A = 0.0
 103       DO ibin = 1, nbin_a
 104         A = A + area_wet_a(ibin)
 105       ENDDO
 106
 107       ! no aerosol surface area - no uptake
 108       IF (A > 0.0) THEN
 109         ! first order uptake, Fuchs and Sutugin, 1971
 110         ! dCg = 1/4 * gamma * A * |v_mol| * Cg * dt
 111         delta_gly = 0.25 * gamma_gly * A * omega * gas(igly) * dtchem
 112
 113         ! avoid negative concentrations
 114         delta_gly = MIN(gas(igly), delta_gly)
 115
 116         ! update partitioning
 117         gas(igly) = gas(igly) - delta_gly
 118
 119         ! distribute onto bins according to fraction of surface area
 120         DO ibin = 1, nbin_a
 121           frac_A = area_wet_a(ibin) / A
 122           ! we take the "photochemical" glysoa aerosol as surrogate
 123           aer(iglysoa_sfc_a, jtotal, ibin) = aer(iglysoa_sfc_a, jtotal, ibin) &
 124                                        + frac_A * delta_gly
 125         ENDDO
 126       ENDIF
 127
 128       END SUBROUTINE glysoa_simple
 129
 130       SUBROUTINE glysoa_complex_derivs(x, y, dt, dydx)
 131
 132       USE module_data_mosaic_therm, ONLY: conv1a,  & ! converts q/mol(air) to nq/m^3 (q = mol or g)
 133                                           p_atm,   & ! pressure (atm)
 134                                           t_k        ! temperature (K)
 135
 136       REAL(kind=8), INTENT(IN)  :: x, y(nspecs), dt
 137       REAL(kind=8), INTENT(OUT) :: dydx(nspecs)
 138
 139       REAL(kind=8), PARAMETER   :: eps      = 1.e-16 , & ! minimum allowed concentration in reservoirs
 140                                    Kh_water = 4.19e5 , & ! effective Henry's law constant of glyoxal in pure water (M atm-1) (Ip et al., GRL, 2011)
 141                                    Kh_oh    = 25.0,    & ! Henry's law constant of OH in pure water (M atm-1) (Klaening et al., 1985)
 142                                    k_oh     = 1.1e9      ! OH reaction rate (mol L-1 s-1) (Ervens and Volkamer, ACP, 2010)
 143
 144       REAL(kind=8)              :: gly_g_atm,  & ! gas-phase concentration in atm
 145                                    f_A1,       & ! fraction of glyoxal in reservoir 1
 146                                    tau1,       & ! characteristical timescale reservoir 1
 147                                    tau2,       & ! characteristical timescale reservoir 2
 148                                    oh_g_atm,   & ! gas-phase OH concentration in atm
 149                                    oh_a,       & ! liquid-phase OH concentration
 150                                    c_tot,      & ! total concentration of dissolved salts (M)
 151                                    Kh_eq,      & ! Henry's law constant at equilibrium (M atm-1)
 152                                    gly_ptot_eq,& ! total glyoxal concentration (reservoirs 1 and 2) at equilibrium
 153                                    gly_r1_eq,  & ! glyoxal concentration (reservoir 1) at equilibrium
 154                                    gly_r2_eq,  & ! glyoxal concentration (reservoir 2) at equilibrium
 155                                    anh4,       & ! ammonium-ion activity (constrained)
 156                                    kII, kI,    & ! second and first order dark rate constants
 157                                    omega         ! mean molecular velocity of glyoxal (cm/s)
 158
 159       ! tendencies
 160       REAL(kind=8)              :: dg_r1,      & ! from gas-phase to reservoir 1
 161                                    dr1_r2,     & ! reservoir 1 to reservoir 2 (or vice versa)
 162                                    dr1_nh4,    & ! reservoir 1 to nh4
 163                                    dr1_oh,     & ! reservoir 1 to oh
 164                                    dg_sfc        ! gas-phase to surface uptake
 165
 166       REAL(kind=8)              :: accloss,    & ! acc. loss
 167                                    scaling       ! tendency scaling
 168
 169       dg_r1   = 0.0
 170       dr1_r2  = 0.0
 171       dr1_nh4 = 0.0
 172       dr1_oh  = 0.0
 173       dg_sfc  = 0.0
 174
 175       ! convert gas-phase glyoxal concentration
 176       ! -------------------------------------------------------------
 177
 178       gly_g_atm = y(igly_g) / conv1a ! mole / mole
 179       gly_g_atm = gly_g_atm * p_atm  ! atm
 180
 181       ! with A2/A1 = Kolig, and Kolig is 1
 182       f_A1 = 0.5
 183       tau1 = 2.5e2 ! s
 184       tau2 = 5.5e3 ! s
 185       IF ( y(ic_as) + y(ic_an) .GT. 12.0 ) THEN
 186         ! with A2/A1 = Kolig, and Kolig is 0.5
 187         f_A1 = 0.6667
 188         tau1 = 4.4e4 ! s
 189         tau2 = 4.7e4 ! s
 190       ENDIF
 191
 192       ! Kampff et al., ES&T, 2013, submitted: kinetic limitation of
 193       ! salting-in for SO4 concentrations > 12 M.
 194       c_tot = MIN( 12.0, y(ic_as) + y(ic_an) )
 195
 196       ! effective Henry's law constant including salting-in effect (Kampff et al.)
 197       ! at equilibrium
 198       ! derived from eqn. 3 in Kampff et al., -0.24 is "salting-in" constant
 199       Kh_eq = Kh_water / 10**(-0.24D0 * c_tot) ! mol L-1 atm-1 == mol kg-1 atm-1
 200
 201       ! gly_g_atm in atm, Kh_eq in mol kg-1 atm-1, water in kg m-3 air
 202       ! total glyoxal concentration (reservoirs 1 and 2) at equilibrium
 203       gly_ptot_eq = gly_g_atm * Kh_eq * y(iwater) * 1e9 ! nmol / m3 air
 204
 205       gly_r1_eq   = gly_ptot_eq *        f_A1  ! in reservoir 1 at equilibrium
 206       gly_r2_eq   = gly_ptot_eq * (1.0 - f_A1) ! in reservoir 2 at equilibrium
 207
 208       ! process tendencies in nmol m-3 s-1
 209       ! from gas-phase to reservoir 1
 210       ! >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 211       IF (.NOT. lfast_tau1) THEN
 212       ! <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 213         dg_r1       = (1.0/tau1) * (gly_r1_eq - y(igly_r1))
 214       ! >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 215       ENDIF
 216       ! <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 217       ! reservoir 1 to reservoir 2
 218       dr1_r2      = (1.0/tau2) * (gly_r2_eq - y(igly_r2))
 219
 220       ! OH reaction (Ervens and Volkamer, Atm. Chem. Phys., 2010)
 221       ! -------------------------------------------------------------
 222
 223       oh_g_atm = y(ioh_g) / conv1a                  ! now its in mole / mole
 224       oh_g_atm = oh_g_atm * p_atm                   ! partial pressure (atm)
 225       oh_a     = oh_g_atm * Kh_oh * y(iwater) * 1e9 ! nmol / m3 air
 226
 227       dr1_oh   = k_oh * y(igly_r1) * oh_a           ! nmole m-3 s-1
 228
 229       ! dark pathway (Noziere et al., J. Phys. Chem., 2009)
 230       ! -------------------------------------------------------------
 231
 232       ! second-order rate constant (mol-1 kg s-1):
 233       anh4 = MAX( 0.0, MIN( 4.0, y(ia_nh4)) ) ! restrict to measured range
 234       kII   = 2.e-10 * exp(1.5 * anh4) * exp(2.5  * y(iph))
 235
 236       ! dark process uses reservoir 1
 237       kI = kII * y(igly_r1) / y(iwater) * 1e-9 ! s-1
 238
 239       dr1_nh4 = kI * y(igly_r1) ! nmol m-3 s-1
 240
 241       ! surface uptake (Ervens and Volkamer, ACP, 2010)
 242       ! -------------------------------------------------------------
 243
 244       ! mean molecular velocity of glyoxal (cm/s)
 245       omega = 1.455e4 * sqrt(t_k / 58.0D0)
 246
 247       ! first order uptake, Fuchs and Sutugin, 1971
 248       ! dCg = 1/4 * gamma * A * |v_mol| * Cg * dt,
 249       ! gamma downscaled to 1.0e-3 according to Waxman et al., 2013
 250       dg_sfc = 0.25D0 * 1.e-3 * y(iarea) * omega * y(igly_g)
 251
 252       ! Numerical integration
 253       ! -------------------------------------------------------------
 254
 255       ! check for undershoots, avoid negative concentrations
 256       ! while ensuring we don't loose mass
 257       IF ( y(igly_g) < eps ) THEN
 258         dg_r1  = 0.0
 259         dg_sfc = 0.0
 260       ELSE
 261         accloss = (dg_r1 + dg_sfc) * dt
 262         IF ( ( y(igly_g) - accloss ) < eps ) THEN
 263           scaling = y(igly_g) / (accloss + eps)
 264           dg_r1  = dg_r1  * scaling
 265           dg_sfc = dg_sfc  * scaling
 266         ENDIF
 267       ENDIF
 268
 269       IF ( y(igly_r1) < eps ) THEN
 270         dr1_r2  = 0.0
 271         dr1_nh4 = 0.0
 272         dr1_oh  = 0.0
 273       ELSE
 274         accloss = (-dg_r1 + dr1_r2 + dr1_nh4 + dr1_oh) * dt
 275         IF ( ( y(igly_r1) - accloss ) < eps) THEN
 276           scaling = y(igly_r1) / (accloss + eps)
 277           dr1_r2  = dr1_r2  * scaling
 278           dr1_nh4 = dr1_nh4 * scaling
 279           dr1_oh  = dr1_oh  * scaling
 280         ENDIF
 281       ENDIF
 282
 283       IF ( y(igly_r2) < eps ) THEN
 284         dr1_r2  = MAX( 0.0, dr1_r2 )
 285       ELSE
 286         accloss = -dr1_r2 * dt
 287         IF ( ( y(igly_r2) - accloss ) < eps ) THEN
 288           scaling = y(igly_r2) / (accloss + eps)
 289           dr1_r2  = dr1_r2  * scaling
 290         ENDIF
 291       ENDIF
 292
 293       ! sum tendencies
 294       dydx(igly_g)   = -dg_r1                          -dg_sfc
 295       dydx(igly_r1)  =  dg_r1 -dr1_r2 -dr1_nh4 -dr1_oh
 296       dydx(igly_r2)  =         dr1_r2
 297       dydx(igly_nh4) =                +dr1_nh4
 298       dydx(igly_oh)  =                          dr1_oh
 299       dydx(igly_sfc) =                                  dg_sfc
 300
 301       END SUBROUTINE glysoa_complex_derivs
 302
 303       SUBROUTINE glysoa_complex(dtchem)
 304
 305       USE module_data_mosaic_therm, ONLY : jaerosolstate, all_liquid, mixed, &
 306                                            jtotal, jliquid, nbin_a, &
 307                                            area_wet_a, gas, water_a, aer, mc, &
 308                                            ph, a_nh4, c_as, c_an, &
 309                                            igly, iho, &
 310                                            iglysoa_r1_a, iglysoa_r2_a, &
 311                                            iglysoa_oh_a, &
 312                                            iglysoa_nh4_a, iglysoa_sfc_a, &
 313                                            iso4_a, ino3_a, inh4_a, jc_h
 314
 315       ! >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 316       USE module_data_mosaic_therm, ONLY: conv1a,  & ! converts q/mol(air) to nq/m^3 (q = mol or g)
 317                                           p_atm      ! pressure (atm)
 318       ! <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 319
 320       REAL(kind=8), INTENT(IN)  :: dtchem
 321
 322       REAL(kind=8)              :: A, conv, y(nspecs), yout(nspecs), &
 323                                    dydx(nspecs), gly_g
 324
 325       INTEGER                   :: i, ii, nbin_proc, bin_proc(nbin_a)
 326
 327       ! >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 328       REAL(kind=8), PARAMETER   :: Kh_water = 4.19e5     ! effective Henry's law constant of glyoxal in pure water (M atm-1) (Ip et al., GRL, 2011)
 329
 330       REAL(kind=8)              :: gly_g_atm,  & ! gas-phase concentration in atm
 331                                    f_A1,       & ! fraction of glyoxal in reservoir 1
 332                                    c_tot,      & ! total concentration of dissolved salts (M)
 333                                    Kh_eq,      & ! Henry's law constant at equilibrium (M atm-1)
 334                                    gly_ptot_eq,& ! total glyoxal concentration (reservoirs 1 and 2) at equilibrium
 335                                    gly_r1_eq,  & ! glyoxal concentration (reservoir 1) at equilibrium
 336                                    deltagly      ! delta to bring r1 in equilibrium
 337       ! <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 338
 339
 340       ! liquid / mixed phase bins available?
 341       ! ---------------------------------------------------------------
 342       nbin_proc   = 0
 343       bin_proc(:) = -1 ! see which bins are either mixed,
 344                        ! or all_liquid (only these are to be processed)
 345       DO i = 1, nbin_a
 346         IF (jaerosolstate(i) == all_liquid .OR. &
 347             jaerosolstate(i) == mixed) THEN
 348           nbin_proc = nbin_proc + 1
 349           bin_proc(nbin_proc) = i
 350         ENDIF
 351       ENDDO
 352       IF (nbin_proc == 0) RETURN
 353
 354       ! aerosol surface area available?
 355       ! ---------------------------------------------------------------
 356
 357       A     =     0.0 ! total aerosol surface area (cm^2 / cm^3)
 358       DO i = 1, nbin_proc
 359         ii = bin_proc(i)
 360         A = A + area_wet_a(ii)
 361       ENDDO
 362       IF (A <= 0) RETURN
 363
 364       ! clean diagnostic arrays
 365       ! ---------------------------------------------------------------
 366
 367       ph(:) = -9999.0 ! aerosol pH
 368       a_nh4(:) = 0.0  ! ammonium ion activity (M, mol/m^3)
 369       c_as(:) = 0.0   ! ammonium sulfate concentration (M, mol/kg)
 370       c_an(:) = 0.0   ! ammonium nitrate concentration (M, mol/kg)
 371
 372       ! get gas-phase
 373       ! ---------------------------------------------------------------
 374
 375       ! gly_g will be re-used for all bins
 376       gly_g     = gas(igly)         ! nmol / m3
 377
 378       DO i = 1, nbin_proc
 379
 380         ii = bin_proc(i)
 381
 382         ! load concentrations array
 383         ! -------------------------------------------------------------
 384
 385         conv        = 1.e-9 / water_a(ii) ! nmol/m^3 (air) -> mol/kg (water)
 386
 387         y(:)        = 0.0
 388
 389         ! nmol/m^3
 390         y(igly_g)   = gly_g
 391         y(igly_r1)  = aer(iglysoa_r1_a,jtotal,ii)
 392         y(igly_r2)  = aer(iglysoa_r2_a,jtotal,ii)
 393         y(igly_nh4) = aer(iglysoa_nh4_a,jtotal,ii)
 394         y(igly_sfc) = aer(iglysoa_sfc_a,jtotal,ii)
 395         y(igly_oh)  = aer(iglysoa_oh_a,jtotal,ii)
 396
 397         y(ic_as)    = aer(iso4_a,jliquid,ii) * conv ! assume we can only form (NH4)2SO4
 398         y(ic_an)    = aer(ino3_a,jliquid,ii) * conv ! assume we can only form NH4NO3
 399         y(ia_nh4)   = aer(inh4_a,jliquid,ii) * conv ! set activity == concentration
 400
 401         y(iph)      = MIN(14.0D0, MAX(0.0D0, -log10(mc(jc_h,ii)) ))
 402
 403         y(iwater)   = water_a(ii) ! kg m-3
 404
 405         y(ioh_g)    = gas(iho) ! nmol m-3
 406
 407         y(iarea)    = area_wet_a(ii) ! cm^2 / cm^3
 408
 409         ! >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 410         IF (lfast_tau1) THEN
 411           ! instantaneous Henry's law equilibrium
 412           ! for fast partitioning (in case characteristic timescale
 413           ! of Kampf et al. is only artefact...)
 414           gly_g_atm = y(igly_g) / conv1a ! mole / mole
 415           gly_g_atm = gly_g_atm * p_atm  ! atm
 416
 417           ! with A2/A1 = Kolig, and Kolig is 1
 418           f_A1 = 0.5
 419           IF ( y(ic_as) + y(ic_an) .GT. 12.0 ) THEN
 420             f_A1 = 0.6667
 421           ENDIF
 422
 423           ! Kampff et al., ES&T, 2013, submitted: kinetic limitation of
 424           ! salting-in for SO4 concentrations > 12 M.
 425           c_tot = MIN( 12.0, y(ic_as) + y(ic_an) )
 426
 427           ! effective Henry's law constant including salting-in effect (Kampff et al.)
 428           ! at equilibrium
 429           ! derived from eqn. 3 in Kampff et al., -0.24 is "salting-in" constant
 430           Kh_eq = Kh_water / 10**(-0.24D0 * c_tot) ! mol L-1 atm-1 == mol kg-1 atm-1
 431
 432           ! gly_g_atm in atm, Kh_eq in mol kg-1 atm-1, water in kg m-3 air
 433           ! total glyoxal concentration (reservoirs 1 and 2) at equilibrium
 434           gly_ptot_eq = gly_g_atm * Kh_eq * y(iwater) * 1e9 ! nmol / m3 air
 435
 436           gly_r1_eq   = gly_ptot_eq *        f_A1  ! in reservoir 1 at equilibrium
 437
 438           deltagly    = gly_r1_eq - y(igly_r1)
 439
 440           y(igly_g)   = y(igly_g) - deltagly
 441           y(igly_r1)  = y(igly_r1) + deltagly
 442
 443           IF (y(igly_g) < 0.0 .OR. y(igly_r1) < 0.0) THEN
 444             WRITE(*,*) "THIS IS NOT RIGHT: ",y(igly_g), y(igly_r1),deltagly
 445           ENDIF
 446         ENDIF
 447         ! <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 448
 449         ! integrate system
 450         ! -------------------------------------------------------------
 451         CALL glysoa_complex_derivs(1._8, y, dtchem, dydx)
 452         CALL rk4(y, dydx, nspecs, 1._8, dtchem, yout, glysoa_complex_derivs)
 453
 454         ! update transported fields
 455         aer(iglysoa_r1_a,jtotal,ii)  = yout(igly_r1)
 456         aer(iglysoa_r2_a,jtotal,ii)  = yout(igly_r2)
 457         aer(iglysoa_nh4_a,jtotal,ii) = yout(igly_nh4)
 458         aer(iglysoa_sfc_a,jtotal,ii) = yout(igly_sfc)
 459         aer(iglysoa_oh_a,jtotal,ii)  = yout(igly_oh)
 460
 461         ! do not put gas-phase glyoxal back yet, as we will use it
 462         ! for the next bin as well...
 463         gly_g                        = yout(igly_g)
 464
 465         ! save diagnostics
 466         c_as(ii)  = yout(ic_as)
 467         c_an(ii)  = yout(ic_an)
 468         ph(ii)    = yout(iph)
 469         a_nh4(ii) = yout(ia_nh4)
 470
 471       ENDDO
 472
 473       ! update gas-phase reservoir, after all bins have been treated.
 474       gas(igly) = gly_g
 475
 476       END SUBROUTINE glysoa_complex
 477
 478
 479       END MODULE module_mosaic_gly
 480