Fastjx_speedup

[ctessum]

@jialinl6 would you be willing to take a look at this?

The first commit just makes Fast-JX run (about 20x) faster without changing the results it gives.

When I was working on that, however, I noticed that WL says it is effective wavelengths of bins covering 117-850 nm, and then below that we loop through each of the wavelengths and multiply by the difference between one effective wavelength and the next one.

This seems kind of like we're doing an integral, e.g. $\int j dW \approx \sum_i j_i \Delta W_i$ .

If that's correct, then it seems like we might want to be using the width of each bin instead of the difference between the center of one bin and the next one.

So in the second commit I changed it to calculate the bin widths (assuming that the first bins starts at 117 and the last one ends at 850, and then use those bin widths to calculate the J values.

Does that seem correct?

[jialinl6]

Use WL edges

Hi, @ctessum. Thank you for pointing that out! I just found out that I shouldn't do the integral using either effective or boundaries wavelengths.

The basic formula to calculate the photolysis rate is: J = ʃ Φ(λ,T) σ(λ,T) F(θ,λ) dλ, which requires an integral. But the wavelength bins in Fast-JX are very odd. At low wavelengths (i.e. hard UV), they are actually non-continuous, with a single “bin” actually being a weighted combination of several sub-bins. Because I couldn't find the accurate boundaries of the 18 bins, I compromised to use the effective wavelength to do the integral. But it turned out that I shouldn't. In fast-jx's Fortran codes, I found that GEOS-CHEM only uses the effective bins instead of boundary bins to calculate the J rate and it just summed up the actinic_flux * cross-section in each bin.

[1182 line in https://github.com/geoschem/geos-chem/blob/main/GeosCore/fast_jx_mod.F90] QQQT is the interpolated cross-section, and FFF is the actinic flux, where K represents the number in 18 effective wavelength bins, and L represents the layer number in a column.

So in Fast-JX.jl, we don't need to integrate actinic_flux cross-section quantum yield with respect to wavelength. The function should look like this:

function j_mean_NO2(flux,time, lat, Temperature)
    csa = cos_solar_zenith_angle(lat, time)
    actinic_flux = csa*flux
    j = 0
    if Temperature <= 200
        j = sum(actinic_flux .* table_σ_NO2_jx[1] * ϕ_NO2_jx)
    elseif Temperature >= 294
        j = sum(actinic_flux .* table_σ_NO2_jx[2] * ϕ_NO2_jx)
    elseif 200 < Temperature < 294
        k = (Temperature - 200) / (294 - 200)
        σ = (table_σ_NO2_jx[2] - table_σ_NO2_jx[1]) * k + table_σ_NO2_jx[1]
        j = sum(actinic_flux .* σ * ϕ_NO2_jx)
    end

I will change all these functions and update the solar zenith angle function that accounts for longitude.

[ctessum]

Sounds good!

[ctessum]

I'm going to merge the version of this that just speeds up FastJX without making any changes to the result, and then you can start with that version to make the changes you mention above. Thanks!