Conversation
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## master #230 +/- ##
==========================================
- Coverage 78.49% 78.47% -0.03%
==========================================
Files 56 56
Lines 9755 9805 +50
==========================================
+ Hits 7657 7694 +37
- Misses 2098 2111 +13 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
Corrected typo in README regarding coagulation kernel.
jcurtis2
left a comment
There was a problem hiding this comment.
Mostly looks straightforward to get this merged. Here are my comments:
- As I can't offer suggestions on files you haven't touched: quickly just add the output directory
scenarios/7_drydep/outto the.gitignore. Looks like we've gotten lazy about this for new scenarios but lets at least fix this one here. - You've added this
do_aero_dilutionto the sectional spec file to disable dilution because the modal lacks dilution entirely so that you can compare the two fairly. In general, this is likely unnecessary as we typically disable dilution by having our input files set in a particular way that causes this process to be zero. Dilution depends on (1) background mixing rate and (2) mixing height changes So for (1) we've controlled it by havingrate = 0inaero_back.datso that there is no mixing, but we've keep mixing from changes in box height inheight.dat. Mixing height does impact the dry deposition rates and the results will change if we changeheight.datto lose our boundary layer evolution. I don't know if we consider that important for this scenario at this time. So this feels like a workaround for a deficiency in modal model by introducing something permanent in the sectional model. But this all isn't 100% without impact, its just whether or not that impact matters. - Possibly worth a sentence in the PR description that you change the particle size distribution (larger size, broader distribution) where you discuss the example scenario.
- Not directly related to your changes but since you are making edits in that scenario directory and comparing the modal to sectional, you have the units of
tot_mass_concin the modal processing code to be micrograms per cubic meter.
| call bin_kernel(bin_grid_size(bin_grid), bin_grid%centers, aero_data, & | ||
| run_sect_opt%coag_kernel_type, env_state, k_bin) | ||
| call smooth_bin_kernel(bin_grid_size(bin_grid), k_bin, ck) | ||
| do i = 1,bin_grid_size(bin_grid) | ||
| do j = 1,bin_grid_size(bin_grid) | ||
| ck(i,j) = ck(i,j) * 1d6 ! m^3/s to cm^3/s | ||
| end do | ||
| end do | ||
|
|
||
| ! multiply kernel with constant timestep and logarithmic grid distance | ||
| do i = 1,bin_grid_size(bin_grid) | ||
| do j = 1,bin_grid_size(bin_grid) | ||
| ck(i,j) = ck(i,j) * run_sect_opt%del_t * bin_grid%widths(i) | ||
| end do | ||
| end do |
There was a problem hiding this comment.
Not your mistake, but I think we should take this time to wrap this code in a if (run_sect_opt%do_coagulation) to ensure that this is only called when coagulation is on.
Doing this change, you can fix your sectional spec files that required you to have do_coagulation yes and setting the kernel to zero.
|
|
||
| end do | ||
| else | ||
| return |
There was a problem hiding this comment.
Change this to a die_msg to handle the other loss processes that aren't currently implemented so that if someone turns them on, they fail. Right now, this reads like someone could have these loss processes on and have the simulation run successfully.
The modal model (scenario_update_aero_modes) has a similar issue that I didn't catch in the previous review. It felt more complete when looking at it as it added a check for SCENARIO_LOSS_FUNCTION_CONSTANT (although it only returns) but it doesn't handle the others (there's no final else in that function).
| call spec_file_read_logical(file, 'do_optical', run_sect_opt%do_optical) | ||
| if (run_sect_opt%do_optical) then | ||
| call spec_file_die_msg(527436819, file, & | ||
| "sectional run does not support optical properties calculation") | ||
| end if |
There was a problem hiding this comment.
While neither the optical properties or MOSAIC are supported, move this block inside the logic for checking MOSAIC. The particle model (in run_part.F90) checks for do_mosaic and then if do_mosaic is true, it will check if optical properties are true. It will just make it easier in the future when MOSAIC is supported.
Lines 204 to 212 in 7984020
| !> Updates an array (of size equal to the number of sections) containing | ||
| !> the dry deposition velocity for each bin in a sectional simulation. | ||
| !! | ||
| !! This is a diagnostic helper: it returns deposition velocities (m s^{-1}) | ||
| !! rather than loss rates, by multiplying the rate from | ||
| !! scenario_loss_rate_drydep() back by the mixing layer height. It is not | ||
| !! used by the sectional time-stepping itself, which calls | ||
| !! scenario_binned_loss(). | ||
| subroutine scenario_section_drydep_rates(scenario, bin_grid, aero_data, & | ||
| env_state, rates) | ||
|
|
||
| !> Scenario data. | ||
| type(scenario_t), intent(in) :: scenario | ||
| !> Bin grid. | ||
| type(bin_grid_t), intent(in) :: bin_grid | ||
| !> Aerosol data. | ||
| type(aero_data_t), intent(in) :: aero_data | ||
| !> Environmental state. | ||
| type(env_state_t), intent(in) :: env_state | ||
| !> Deposition velocities for each section/bin (m s^{-1}). | ||
| real(kind=dp), intent(inout) :: rates(:) | ||
|
|
||
| integer :: i_bin | ||
| real(kind=dp) :: density, vol | ||
|
|
||
| call assert_msg(516274839, size(rates) == bin_grid_size(bin_grid), & | ||
| "rates array size must match the number of bins") | ||
|
|
||
| density = aero_data%density(1) | ||
|
|
||
| do i_bin = 1,bin_grid_size(bin_grid) | ||
| vol = aero_data_rad2vol(aero_data, bin_grid%centers(i_bin)) | ||
| rates(i_bin) = scenario_loss_rate_drydep(vol, density, & | ||
| aero_data, env_state, scenario) * env_state%height | ||
| end do | ||
|
|
||
| end subroutine scenario_section_drydep_rates |
There was a problem hiding this comment.
Is this worth keeping/maintaining?
|
|
||
| if (present(scenario)) then | ||
| call drydep_params_input_netcdf(scenario%drydep, ncid) | ||
| end if |
There was a problem hiding this comment.
This is called unconditionally but should be guarded. The drydep_params_output_netcdf is only called when the loss process is dry deposition
if (scenario%loss_function_type == SCENARIO_LOSS_FUNCTION_DRYDEP) then
call drydep_params_output_netcdf(scenario%drydep, ncid)
end if
but drydep_params_input_netcdf would always expect it provided scenario is passed to input_sectional.
| only active process. The sectional case runs with the zero | ||
| coagulation kernel. It also sets do_aero_dilution to no; the modal | ||
| representation has no dilution treatment. |
There was a problem hiding this comment.
Just a note to edit this if changes are made regarding coagulation or dilution.
There was a problem hiding this comment.
You can remove this spec file entirely - this directory isn't doing any dilution.
There was a problem hiding this comment.
Change the spec file and related output processing for the sectional run as the particle run has no dilution.
Co-authored-by: Jeffrey Curtis <jcurtis2@illinois.edu>
Co-authored-by: Jeffrey Curtis <jcurtis2@illinois.edu>
Summary
Adds dry deposition support to the sectional representation. Also corrects the
deposition velocity formula shared by all representations to match
Zhang et al. (2001) and Emerson et al. (2020) as published.
Sectional dry deposition
scenario_binned_loss()applies per-bin loss each timestep via the samescenario_loss_rate()dispatcher already used by the particle-resolvedand modal paths, so the sectional model supports every
loss_function_type(constant, volume, drydep, chamber), not just drydeposition.
Each bin decays as
m(t+dt) = m(t) * exp(-rate*dt)— the exact solutionof
dm/dt = -rate*m. Mass (vol_conc) isprognostic and
num_concis diagnosed, matching the convention already used bysectional coagulation.
output_sectional()/input_sectional()gain ascenarioargument sodry deposition parameters are written to and read from NetCDF output,
mirroring
output_modal()/input_modal().Deposition velocity formula change
The
R_a*R_s*V_sterm comes from the two-layer flux balance(Seinfeld and Pandis, 1998); however, it is not a physically representative term, but
rather a direct artifact of the mathematical manipulation of the multilayer resistance analogy
system of equations with settling acting through all layers. Strictly speaking, neither
Zhang et al. (2001) nor Emerson et al. (2020) use this term. Both follow Slinn (1982)'s
derivation where settling is treated purely in parallel. Since
drydep_params_tis explicitlyparameterized to switch between Z01-and E20 coefficient sets, retaining the cross
term would implement neither parameterization as published.
Applied in all three places:
scenario_loss_rate_drydep(particle/sectional),
scenario_integrated_loss_rate_drydep(modal,analytic), and
dep_vel_integrand(modal, QUADPACK).New spec-file option
do_aero_dilution(sectional runs) skips both background-exchangedilution and the mixing-layer-growth correction when set to
no. Themodal representation has no dilution treatment at all, so this is needed
for the modal/sectional comparison and lets dry deposition be
studied in isolation.
Tests
test_loss_16— sectional dry deposition regression testtest_loss_17— checks sectional against particle-resolvedReference files for
test_loss_13,test_loss_14(QUADPACK), andtest_loss_16are regenerated due to formula update.Example scenario
scenarios/7_drydepnow demonstrates both modal and sectional representations. Thesectional case is added as
drydep_sect.specwith1_run_sect.shand2_process_sect.sh, mirroring the existing modal scripts. It uses thestandard
extract_sectional_*tools rather than a dedicatedpost-processor. A
READMEdocuments both cases.