Use of newtonian relaxation within the moist model - #305
Draft
sit23 wants to merge 59 commits into
Draft
Conversation
…thon. This enables comparison between annual-average insolation, and e.g. the Frierson specified incoming radiation.
…own_annual_av Merging changes from two machines into top-down-annual-av branch
Adding in README changes
Including README changes
Adding in newtonian cooling as a radiation scheme option.
… mixed-layer. This sets the mixed-layer temperature to the same as the ground temperature, allowing all the surface fluxes to be calculated as normal, but they do not affect the surface temperature, as this stays fixed.
…ns I can run with the monin-obukhov bl scheme only, without also having a rayleigh-drag at the surface.
…ture field, rather than radiative-convective equilibrium. This has the advantage of then feeling justified in using newtonian cooling as a radiaition scheme in a more complex model.
…ective equilibrium calculation, but using radiative equilibrium for equilibrium temperature.
…n runs, whereby the relative humidty is set to 100 percent in the bottom model layer, and the mixed-layer bc is turned off completely. That way there are no worries about dubious surface fluxes.
…tom model layer to something that is not the atmospheric timestep. That way we can allow moist convection to alter the b/l sphum without any changes being immediatly overruled.
…ium. That way we should get a more convectively unstable atmosphere.
…aranteed to be suitably colder than the surface and so have reasonable surface fluxes.
Merging latest master changes to docs and python.
…Means can have seasonal varying teq in rad equil modes.
…op_down_with_moisture Merging other top-down changes.
…of solar constant. Now incoming sw varies with time when orbit is eccentric.
…at ls=0 as perihelion is always at true-anomaly = 0.. Offset of 1.905637 radians is so that perihelion at ls=251 degrees, as on lmd website. Also as in sol-ls.pro program that Luca wrote that is currently in sharecode. Problem now is that solstices etc in wrong place.
…. Consistent with earth values, etc.
Merging in latest master changes
…moisture source in the model. This often makes it crash if you are doing dry planets other than Earth. Therefore added option, as in mars-dev commit 1ca7b41, such that false temperatures can be fed to escomp so that it doesn't crash, then real temperatures are subsequently used.
Merge in latest master changes.
Merging latest 3 master commits
Conflicts: src/atmos_param/hs_forcing/hs_forcing.F90 src/atmos_spectral/driver/solo/idealized_moist_phys.F90
…sed with gfdl astronomy. Useful as can now apply existing mars orbit numbers to this.
…h t-grnd when it was not present.
…h that this model will assume no short-wave absorption, which we definitely need for Titan. Aborting this for now.
…ure_2026 # Conflicts: # postprocessing/plevel_interpolation/scripts/run_plevel.py # src/atmos_param/hs_forcing/hs_forcing.F90 # src/atmos_spectral/driver/solo/atmosphere.F90 # src/atmos_spectral/driver/solo/idealized_moist_phys.F90 # src/atmos_spectral/driver/solo/mixed_layer.F90 # src/coupler/surface_flux.F90 # src/extra/python/scripts/remove_certain_restart_and_data_files.py # src/shared/astronomy/astronomy.f90 # src/shared/constants/constants.F90
…ostics Both were unconditional and used real (unclamped) temperatures, so on a dry/cold planet run with use_actual_surface_temperatures=.false. (e.g. radiative_eq_mars) they could crash lookup_es_2d with a saturation vapor pressure table overflow, even though neither diagnostic (rh_2m, e_sat) was requested: - surface_flux.F90: the 2m/10m diagnostic block's escomp(temp_2m, ...) now only runs when use_actual_surface_temperatures is true, matching how the other escomp calls in this subroutine already treat that flag. - idealized_moist_phys.F90: escomp(t_surf, e_sat_out) now only runs when the e_sat diagnostic (id_e_sat) is actually registered, matching the existing guard on the send_data call right after it. Neither call existed in top_down_with_moisture's original file - both were added to master after the fork. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
Isolated to two independent issues, both now fixed: 1. hs_forcing's generic tracer_source_sink (trflux/trsink, defaulting to a nonzero surface flux/sink meant for an optional test tracer) applies to every tracer it's passed, including sphum, unless the field_table registers a 'tracer_sms' override for that tracer - which the grey model's field_table doesn't. radiative_eq_mars routes "radiation" through hs_forcing (do_newtonian_cooling_as_rad), so this was silently adding a small but steady moisture source to an otherwise fully dry setup every timestep. Fixed by zeroing trflux/trsink in the test case's own hs_forcing_nml, since this experiment has no optional tracer of its own. 2. gcm_vert_diff_down/up pass the moisture tendency both as a dedicated dt_q argument and, aliased to the exact same memory, as part of the full dt_tracers(:,:,:,:) tendency array - even though the sphum slice of that array is never written by the generic tracer loop (vert_diff_init sets tracers(sphum)%do_vert_diff=.false. specifically to avoid that). Passing the same actual array to two dummy arguments that can each be defined (both intent(out) in gcm_vert_diff_up) isn't standard-conforming regardless of which ends up actually written, and is exactly the kind of thing that can behave inconsistently across compilers/optimization levels. Routed the moisture tendency through its own scratch array (dt_qg_vert_diff) so it never shares memory with dt_tracers(:,:,:,:) during these two calls. Also softened surface_flux.F90's q_atm_in>0. dry-model sanity check to a small tolerance (pure floating-point noise guard, not a real threshold change) now that q_surf0 is matched to q_atm whenever use_actual_surface_temperatures=.false., closing the one remaining path that could otherwise inject an unbounded surface-driven moisture source into a nominally dry run. Verified against radiative_eq_mars_test_case.py both shortened (3 mars days) and at its full, unmodified length (30 mars days): compiles cleanly, runs to completion with no FATAL errors, and produces NaN-free output in a physically plausible Mars temperature range (t_surf 96-230K, temp 118-281K). Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
An aquaplanet-style world with land everywhere (land_option='all_land', no input mask needed) and bucket hydrology active at every gridpoint, given a generous non-limiting initial depth. Unlike a standard mixed-layer aquaplanet - where the ocean is an effectively infinite moisture source/sink - every gridpoint here has a finite water reservoir, so total column water (atmosphere + bucket) should be exactly conserved in the absence of numerical bugs, making it a useful check on tracer-conserving vertical diffusion. Used to test whether the gcm_vert_diff_down/up argument-aliasing fix in the previous commit changes global moisture conservation: ran 10 simulated days both with and without that fix (same source otherwise) and compared area-weighted global-mean total water day by day. Results were bit-identical between the two builds with this compiler/optimization level (gfortran -O2) - the aliasing violation doesn't appear to be miscompiled here, though it remains worth having fixed as a matter of correctness, since that's exactly the kind of undefined behaviour that can vary by compiler or optimization level. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
…tp branch do j=1,size(t,1) should be do j=1,size(t,2) - t's first dimension is the local per-PE longitude extent, the second is latitude, and h_trop/zfull/teq are indexed (lon,lat) here. With domain decomposition splitting latitude across PEs (e.g. X-AXIS=128, Y-AXIS=16 per PE, as trip_test's default 4-core runs use), this looped j up to 128 instead of 16, writing 8x past the end of h_trop/zfull's second dimension - a heap-corrupting out-of-bounds write. This is a pre-existing bug from top_down_with_moisture's original 2018-2020 source (confirmed via git show against that branch directly), not something introduced by the merge. master's own (simpler, without the pure_rad_equil option) version of this same loop has the correct size(t,2), which is what caught this: master's version was never wrong, top_down_with_moisture's superset version - which the merge conflict resolution correctly kept for the pure_rad_equil support master lacks - carried the bug forward silently, apparently never having been trip-tested against a build that would surface it as a crash rather than quiet memory corruption. Found via a "double free or corruption (out)" crash in top_down_test (one of trip_test's default 19 cases) immediately after hs_forcing_init's heat-capacity spin-up completes. Verified fixed by rerunning top_down_test_case.py's full configured length (3 x 90-day segments) to completion with no crash - previously it crashed within seconds of model start. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
Contributor
Author
|
Trip test results (summarised by Claude) are as follows: Trip test results:
|
| Result | Count | Test cases |
|---|---|---|
| ✅ Bit-identical to master | 13 | bucket_model, frierson, giant_planet, held_suarez, MiMA, realistic_continents_fixed_sst, realistic_continents_variable_qflux, variable_co2_rrtm, barotropic_vort_eq_stirring, shallow_water_stirring, column_test, top_down_test, variable_co2_grey |
| ❌ Pre-existing failure, unrelated to this branch | 1 | axisymmetric |
| ⏭️ Excluded - Socrates not fully built in this environment | 3 | socrates_aquaplanet, socrates_aquaplanet_cloud, ape_aquaplanet |
| ⏭️ Excluded - trip_test harness gap (see below) | 2 | grey_mars, radiative_eq_mars |
Two real bugs found and fixed along the way
1. Out-of-bounds array write in hs_forcing.F90. top_down_test
crashed with double free or corruption (out). Root cause: the
extend_tp branch of top_down_newtonian_damping looped
do j=1,size(t,1) instead of do j=1,size(t,2) when indexing
h_trop/zfull (which are (lon,lat)), overrunning the local per-PE
latitude extent by up to 8x under a typical 4-core decomposition. Confirmed
via git show against top_down_with_moisture's original 2018-2020
source that this typo predates the merge - master's own (simpler)
version of this same loop already has the correct bound, which is what
caught it. Fixed, and verified by running top_down_test's full configured
length (3 x 90-day segments) to completion with no crash.
2. Undefined-behaviour argument aliasing in gcm_vert_diff_down/_up.
The moisture tendency was passed both as its own dedicated dt_q argument
and, aliased to the same memory, as part of the full dt_tracers(:,:,:,:)
array - not standard-conforming even though the generic tracer loop never
actually touches that slice. Routed it through a dedicated scratch array
instead. Verified via three independent isolation tests (with vs without
just this fix) that it produces no observable difference under this
compiler (gfortran -O2) - worth having fixed as a correctness matter
regardless, since this is exactly the kind of bug that can vary silently
by compiler or optimization level.
On the excluded cases
- The three Socrates-based cases can't build in this environment
(missing Socrates source files) - unrelated to this branch. grey_mars/radiative_eq_marsfail on a run-length/timestep mismatch:
this branch'strip_test_functions.pydoesn't have the Mars-specific
run-length override that exists only on a separate branch
(mars_dust_lh), never merged in here. Not a regression.
radiative_eq_mars- the test case that motivated this whole
reconciliation - has separately been verified extensively outside
trip_test: both a shortened smoke test and the full, unmodified 30
mars-day run complete cleanly, with NaN-free output in a physically
plausible range (t_surf96-230K,temp118-281K).
An apparent variable_co2_grey divergence turned out to be a false alarm
from stale cached data (the original run crashed and was resumed several
times against the same default work/data directories) - confirmed clean
with three independent fresh-directory comparisons and cleaned up the
stale cache.
Contributor
Author
|
So we've verified that these top-down modifications don't alter existing results, and they do allow the radiative-eq-mars test case to run, which is great. I would suggest that the mars dust branch gets merged first, and then we move onto this one. |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
One test case I've been keen to have for a while is a version of the model where Newtonian relaxation of temperature is treated like the 'radiation scheme' of the full model, rather than it's own isolated version. This is implemented in this branch. In addition, this allows for the top-down option within hs_forcing to interact with the surface module, and this in turn is the scaffolding allowing the radiative_eq_mars test case to run properly.
Creating this p/r as a draft for now to enable further testing and trip tests etc.