GEOS-Chem "fullchem" Gas-Phase Mechanism
This is an implementation of the GEOS-Chem "fullchem" gas-phase mechanism in Julia. The original version is used in the GEOS-Chem global 3-D atmospheric chemistry transport model. The version here is adapted from GEOS-Chem version 14.1.1. This mechanism is the result of many journal articles which are cited in API documentation for the GEOSChemGasPhase type.
This implementation is a work in progress. In particular, it does not yet include heterogeneous chemistry.
System overview
First, let's initialize the model and we can also look at the first few ODE equations of the reaction network:
using GasChem, EarthSciMLBase
using DifferentialEquations, ModelingToolkit
using DynamicQuantities, Plots
using ModelingToolkit:t
tspan = (0.0, 360.0)
gc = GEOSChemGasPhase()
equations(gc)[1:5]\[ \begin{align} \frac{\mathrm{d} \mathtt{A3O2}\left( t \right)}{\mathrm{d}t} &= 0.06 \mathtt{j\_69} \mathtt{MEK}\left( t \right) + 0.07 \mathtt{j\_70} \mathtt{RCHO}\left( t \right) + 0.05 \mathtt{j\_98} \mathtt{R4N2}\left( t \right) + 0.36 \mathtt{arrhenius\_121.k}\left( t \right) \mathtt{OH}\left( t \right) \mathtt{RA3P}\left( t \right) - \mathtt{rate\_RO2NO\_a2\_53.k}\left( t \right) \mathtt{NO}\left( t \right) \mathtt{A3O2}\left( t \right) + 0.07 \mathtt{NO}\left( t \right) \mathtt{RCO3}\left( t \right) \mathtt{arrhenius\_70.k}\left( t \right) - \mathtt{NO}\left( t \right) \mathtt{A3O2}\left( t \right) \mathtt{rate\_RO2NO\_b2\_52.k}\left( t \right) + 0.05 \mathtt{NO}\left( t \right) \mathtt{rate\_RO2NO\_b2\_56.k}\left( t \right) \mathtt{R4O2}\left( t \right) + 0.07 \mathtt{RCO3}\left( t \right) \mathtt{arrhenius\_148.k}\left( t \right) \mathtt{MCO3}\left( t \right) + 0.03 \mathtt{RCO3}\left( t \right) \mathtt{arrhenius\_101.k}\left( t \right) \mathtt{HO2}\left( t \right) + 0.07 \mathtt{RCO3}\left( t \right) \mathtt{arrhenius\_117.k}\left( t \right) \mathtt{MO2}\left( t \right) - \mathtt{arrhenius\_115.k}\left( t \right) \mathtt{A3O2}\left( t \right) \mathtt{MCO3}\left( t \right) + \mathtt{arrhenius\_173.k}\left( t \right) \mathtt{Br}\left( t \right) \mathtt{C3H8}\left( t \right) + \mathtt{Cl}\left( t \right) \mathtt{C3H8}\left( t \right) \mathtt{arrhenius\_287.k}\left( t \right) - \mathtt{A3O2}\left( t \right) \mathtt{arrhenius\_113.k}\left( t \right) \mathtt{MCO3}\left( t \right) - \mathtt{A3O2}\left( t \right) \mathtt{constant\_k\_73.k}\left( t \right) \mathtt{MO2}\left( t \right) - \mathtt{A3O2}\left( t \right) \mathtt{HO2}\left( t \right) \mathtt{rate\_RO2HO2\_99.k}\left( t \right) + \mathtt{OH}\left( t \right) \mathtt{C3H8}\left( t \right) \mathtt{tbranch\_51.k}\left( t \right) + 0.05 \mathtt{MCO3}\left( t \right) \mathtt{arrhenius\_132.k}\left( t \right) \mathtt{R4O2}\left( t \right) + 0.03 \mathtt{constant\_k\_85.k}\left( t \right) \mathtt{MO2}\left( t \right) \mathtt{R4O2}\left( t \right) \\ \frac{\mathrm{d} \mathtt{ACET}\left( t \right)}{\mathrm{d}t} &= \mathtt{j\_136} \mathtt{IPRNO3}\left( t \right) - \mathtt{j\_76} \mathtt{ACET}\left( t \right) - \mathtt{j\_77} \mathtt{ACET}\left( t \right) + \mathtt{j\_82} \mathtt{RB3P}\left( t \right) + 0.32 \mathtt{j\_98} \mathtt{R4N2}\left( t \right) + 0.11 \mathtt{PIO2}\left( t \right) \mathtt{NO}\left( t \right) \mathtt{constant\_k\_341.k}\left( t \right) - \mathtt{arrhenius\_171.k}\left( t \right) \mathtt{Br}\left( t \right) \mathtt{ACET}\left( t \right) + \mathtt{NO}\left( t \right) \mathtt{rate\_RO2NO\_a2\_61.k}\left( t \right) \mathtt{B3O2}\left( t \right) + 0.32 \mathtt{NO}\left( t \right) \mathtt{rate\_RO2NO\_b2\_56.k}\left( t \right) \mathtt{R4O2}\left( t \right) + 0.75 \mathtt{constant\_k\_89.k}\left( t \right) \mathtt{B3O2}\left( t \right) \mathtt{MO2}\left( t \right) - \mathtt{Cl}\left( t \right) \mathtt{arrhenius\_288.k}\left( t \right) \mathtt{ACET}\left( t \right) + 0.791 \mathtt{arrhenius\_122.k}\left( t \right) \mathtt{OH}\left( t \right) \mathtt{RB3P}\left( t \right) - \mathtt{acet\_oh\_72.k}\left( t \right) \mathtt{OH}\left( t \right) \mathtt{ACET}\left( t \right) + \mathtt{arrhenius\_143.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{B3O2}\left( t \right) + \mathtt{OH}\left( t \right) \mathtt{IPRNO3}\left( t \right) \mathtt{arrhenius\_382.k}\left( t \right) + 0.32 \mathtt{MCO3}\left( t \right) \mathtt{arrhenius\_132.k}\left( t \right) \mathtt{R4O2}\left( t \right) + \mathtt{MCO3}\left( t \right) \mathtt{B3O2}\left( t \right) \mathtt{arrhenius\_135.k}\left( t \right) + 0.16 \mathtt{constant\_k\_85.k}\left( t \right) \mathtt{MO2}\left( t \right) \mathtt{R4O2}\left( t \right) \\ \frac{\mathrm{d} \mathtt{ACTA}\left( t \right)}{\mathrm{d}t} &= - \mathtt{constant\_k\_285.k}\left( t \right) \mathtt{Cl}\left( t \right) \mathtt{ACTA}\left( t \right) + \mathtt{H2O}\left( t \right) \mathtt{constant\_k\_338.k}\left( t \right) \mathtt{CH3CHOO}\left( t \right) + 0.3 \mathtt{AROMP5}\left( t \right) \mathtt{O3}\left( t \right) \mathtt{constant\_k\_642.k}\left( t \right) + 0.15 \mathtt{AROMP5}\left( t \right) \mathtt{constant\_k\_641.k}\left( t \right) \mathtt{OH}\left( t \right) + \mathtt{arrhenius\_141.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{R4N1}\left( t \right) + \mathtt{arrhenius\_115.k}\left( t \right) \mathtt{A3O2}\left( t \right) \mathtt{MCO3}\left( t \right) + \mathtt{ATO2}\left( t \right) \mathtt{arrhenius\_139.k}\left( t \right) \mathtt{MCO3}\left( t \right) + \mathtt{arrhenius\_147.k}\left( t \right) \mathtt{OTHRO2}\left( t \right) \mathtt{MCO3}\left( t \right) + \mathtt{arrhenius\_138.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{R4O2}\left( t \right) + \mathtt{arrhenius\_146.k}\left( t \right) \mathtt{ETO2}\left( t \right) \mathtt{MCO3}\left( t \right) - \mathtt{arrhenius\_66.k}\left( t \right) \mathtt{OH}\left( t \right) \mathtt{ACTA}\left( t \right) + \mathtt{KO2}\left( t \right) \mathtt{arrhenius\_140.k}\left( t \right) \mathtt{MCO3}\left( t \right) + \mathtt{arrhenius\_143.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{B3O2}\left( t \right) + \mathtt{arrhenius\_116.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{PO2}\left( t \right) + \mathtt{arrhenius\_142.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{PRN1}\left( t \right) + 0.5 \mathtt{OH}\left( t \right) \mathtt{HAC}\left( t \right) \mathtt{rate\_HACOH\_b\_112.k}\left( t \right) + 0.13 \mathtt{MCO3}\left( t \right) \mathtt{arrhenius\_586.k}\left( t \right) \mathtt{HO2}\left( t \right) + \mathtt{MCO3}\left( t \right) \mathtt{MO2}\left( t \right) \mathtt{arrhenius\_131.k}\left( t \right) \\ \frac{\mathrm{d} \mathtt{ALD2}\left( t \right)}{\mathrm{d}t} &= \mathtt{j\_135} \mathtt{ETNO3}\left( t \right) - \mathtt{j\_61} \mathtt{ALD2}\left( t \right) - \mathtt{j\_62} \mathtt{ALD2}\left( t \right) + \mathtt{j\_80} \mathtt{ETP}\left( t \right) + \mathtt{j\_84} \mathtt{PP}\left( t \right) + \mathtt{j\_85} \mathtt{RP}\left( t \right) + 0.32 \mathtt{j\_98} \mathtt{R4N2}\left( t \right) + 0.5 \mathtt{PRPE}\left( t \right) \mathtt{O3}\left( t \right) \mathtt{arrhenius\_103.k}\left( t \right) + \mathtt{arrhenius\_114.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{PO2}\left( t \right) + 0.15 \mathtt{rate\_RO2HO2\_78.k}\left( t \right) \mathtt{KO2}\left( t \right) \mathtt{HO2}\left( t \right) + \mathtt{SO2}\left( t \right) \mathtt{CH3CHOO}\left( t \right) \mathtt{constant\_k\_336.k}\left( t \right) + \mathtt{arrhenius\_381.k}\left( t \right) \mathtt{ETNO3}\left( t \right) \mathtt{OH}\left( t \right) + \mathtt{H2O}\left( t \right) \mathtt{constant\_k\_337.k}\left( t \right) \mathtt{CH3CHOO}\left( t \right) - \mathtt{arrhenius\_41.k}\left( t \right) \mathtt{OH}\left( t \right) \mathtt{ALD2}\left( t \right) + 0.93 \mathtt{arrhenius\_60.k}\left( t \right) \mathtt{NO}\left( t \right) \mathtt{KO2}\left( t \right) + \mathtt{EOH}\left( t \right) \mathtt{Cl}\left( t \right) \mathtt{constant\_k\_284.k}\left( t \right) + \mathtt{EOH}\left( t \right) \mathtt{OH}\left( t \right) \mathtt{constant\_k\_91.k}\left( t \right) + \mathtt{NO}\left( t \right) \mathtt{constant\_k\_334.k}\left( t \right) \mathtt{CH3CHOO}\left( t \right) + \mathtt{NO}\left( t \right) \mathtt{OTHRO2}\left( t \right) \mathtt{arrhenius\_49.k}\left( t \right) + \mathtt{NO}\left( t \right) \mathtt{ETO2}\left( t \right) \mathtt{rate\_RO2NO\_b2\_47.k}\left( t \right) + \mathtt{NO}\left( t \right) \mathtt{arrhenius\_54.k}\left( t \right) \mathtt{PO2}\left( t \right) + \mathtt{NO}\left( t \right) \mathtt{arrhenius\_63.k}\left( t \right) \mathtt{PRN1}\left( t \right) + 0.86 \mathtt{NO}\left( t \right) \mathtt{arrhenius\_58.k}\left( t \right) \mathtt{R4N1}\left( t \right) + 0.32 \mathtt{NO}\left( t \right) \mathtt{rate\_RO2NO\_b2\_56.k}\left( t \right) \mathtt{R4O2}\left( t \right) + 0.5 \mathtt{constant\_k\_90.k}\left( t \right) \mathtt{PRN1}\left( t \right) \mathtt{MO2}\left( t \right) + 2 \left( \mathtt{ETO2}\left( t \right) \right)^{2} \mathtt{constant\_k\_93.k}\left( t \right) + 0.5 \mathtt{constant\_k\_88.k}\left( t \right) \mathtt{KO2}\left( t \right) \mathtt{MO2}\left( t \right) + \mathtt{NO2}\left( t \right) \mathtt{constant\_k\_335.k}\left( t \right) \mathtt{CH3CHOO}\left( t \right) + \mathtt{arrhenius\_137.k}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{PRN1}\left( t \right) + \mathtt{Cl}\left( t \right) \mathtt{OTHRO2}\left( t \right) \mathtt{constant\_k\_282.k}\left( t \right) + \mathtt{Cl}\left( t \right) \mathtt{ETO2}\left( t \right) \mathtt{constant\_k\_281.k}\left( t \right) - \mathtt{Br}\left( t \right) \mathtt{arrhenius\_170.k}\left( t \right) \mathtt{ALD2}\left( t \right) + \mathtt{arrhenius\_147.k}\left( t \right) \mathtt{OTHRO2}\left( t \right) \mathtt{MCO3}\left( t \right) + 0.75 \mathtt{constant\_k\_83.k}\left( t \right) \mathtt{OTHRO2}\left( t \right) \mathtt{MO2}\left( t \right) + 2 \left( \mathtt{OTHRO2}\left( t \right) \right)^{2} \mathtt{constant\_k\_94.k}\left( t \right) + \left( \mathtt{OTHRO2}\left( t \right) \right)^{2} \mathtt{constant\_k\_96.k}\left( t \right) + \mathtt{OTHRO2}\left( t \right) \mathtt{arrhenius\_145.k}\left( t \right) \mathtt{MCO3}\left( t \right) + 0.64 \mathtt{arrhenius\_120.k}\left( t \right) \mathtt{ETP}\left( t \right) \mathtt{OH}\left( t \right) + \mathtt{arrhenius\_134.k}\left( t \right) \mathtt{KO2}\left( t \right) \mathtt{MCO3}\left( t \right) + \mathtt{arrhenius\_146.k}\left( t \right) \mathtt{ETO2}\left( t \right) \mathtt{MCO3}\left( t \right) + \left( \mathtt{ETO2}\left( t \right) \right)^{2} \mathtt{constant\_k\_95.k}\left( t \right) + \mathtt{ETO2}\left( t \right) \mathtt{MCO3}\left( t \right) \mathtt{arrhenius\_144.k}\left( t \right) + 0.75 \mathtt{ETO2}\left( t \right) \mathtt{constant\_k\_82.k}\left( t \right) \mathtt{MO2}\left( t \right) - \mathtt{arrhenius\_42.k}\left( t \right) \mathtt{ALD2}\left( t \right) \mathtt{NO3}\left( t \right) + \mathtt{CH3CHOO}\left( t \right) \mathtt{constant\_k\_333.k}\left( t \right) \mathtt{CO}\left( t \right) + 0.38 \mathtt{constant\_k\_86.k}\left( t \right) \mathtt{MO2}\left( t \right) \mathtt{R4N1}\left( t \right) + 0.32 \mathtt{MCO3}\left( t \right) \mathtt{arrhenius\_132.k}\left( t \right) \mathtt{R4O2}\left( t \right) + 0.75 \mathtt{MCO3}\left( t \right) \mathtt{arrhenius\_136.k}\left( t \right) \mathtt{R4N1}\left( t \right) + 0.16 \mathtt{constant\_k\_85.k}\left( t \right) \mathtt{MO2}\left( t \right) \mathtt{R4O2}\left( t \right) + 0.5 \mathtt{MO2}\left( t \right) \mathtt{constant\_k\_74.k}\left( t \right) \mathtt{PO2}\left( t \right) \\ \frac{\mathrm{d} \mathtt{ALK4}\left( t \right)}{\mathrm{d}t} &= - \mathtt{ALK4}\left( t \right) \mathtt{constant\_k\_290.k}\left( t \right) \mathtt{Cl}\left( t \right) - \mathtt{ALK4}\left( t \right) \mathtt{arrhenius\_55.k}\left( t \right) \mathtt{OH}\left( t \right) - \mathtt{ALK4}\left( t \right) \mathtt{arrhenius\_64.k}\left( t \right) \mathtt{NO3}\left( t \right) \end{align} \]
You can see that each reaction has a rate constant; rate constants are specified at the end of the list of equations:
equations(gc)[end-3:end]\[ \begin{align} \mathtt{arrhenius\_647.k}\left( t \right) &= \left( \frac{\mathtt{arrhenius\_647.K\_300}}{T} \right)^{\mathtt{arrhenius\_647.b0}} \mathtt{arrhenius\_647.a0} e^{\frac{\mathtt{arrhenius\_647.c0}}{T}} \\ \mathtt{arrhenius\_648.k}\left( t \right) &= \left( \frac{\mathtt{arrhenius\_648.K\_300}}{T} \right)^{\mathtt{arrhenius\_648.b0}} \mathtt{arrhenius\_648.a0} e^{\frac{\mathtt{arrhenius\_648.c0}}{T}} \\ \mathtt{arrhenius\_649.k}\left( t \right) &= \left( \frac{\mathtt{arrhenius\_649.K\_300}}{T} \right)^{\mathtt{arrhenius\_649.b0}} \mathtt{arrhenius\_649.a0} e^{\frac{\mathtt{arrhenius\_649.c0}}{T}} \\ \mathtt{constant\_k\_650.k}\left( t \right) &= \mathtt{constant\_k\_650.c} \end{align} \]
Simulation
Now, let's run a simulation and plot the results:
sys = structural_simplify(gc)
vals = ModelingToolkit.get_defaults(sys)
for k in setdiff(unknowns(sys),keys(vals))
vals[k] = 0 # Set variables with no default to zero.
end
prob = ODEProblem(sys, vals, tspan, vals)
sol = solve(prob, AutoTsit5(Rosenbrock23()))
plot(sol, legend = :outertopright, xlabel = "Time (s)",
ylabel = "Concentration (ppb)")