Skip to content

Commit babbb17

Browse files
committed
Merge branch 'oksanaguba/eamxx/fixer-clean' (PR #7553)
Adds (dycore) energy fixer. The infrastructure can be used for another process (and will be used for IEFLX), but the current impl is for the dycore. It is identical to the EAM impl, where pressure adjustment and dycore energy leaks are fixed as one global change to temperature for each column/level as the same temperature tendency. There is a check if fixer debug output is on ( enable_energy_fixer_debug_info ) , which computes global energy after the fixer to confirm that it does match. This debug output is used in a test that has postproc *py script to confirm that the energy is fixed up to a tolerance. For full discussion see #7276 . Climo results ( model vs model, with fixer and without from the same branch) https://web.lcrc.anl.gov/public/e3sm/diagnostic_output/ac.onguba/theta/ fixerJuly24-2025/latlon-pdf-vs-default/viewer/lat_lon/index.html . [nonBFB] for any cime runs with eamxx and dycore. Potentially climate-changing, too.
2 parents f5e6bd7 + 0cbd887 commit babbb17

File tree

17 files changed

+673
-83
lines changed

17 files changed

+673
-83
lines changed

cime_config/tests.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -346,7 +346,7 @@
346346
"ERS_Ld5.TL319_oQU240wLI_ais8to30.MPAS_LISIO_JRA1p5.mpaso-ocn_glcshelf",
347347
"SMS_P12x2.ne4pg2_oQU480.WCYCL1850NS.allactive-mach_mods",
348348
"ERS_Ln9.ne4pg2_ne4pg2.F2010-MMF1.eam-mmf_crmout",
349-
"SMS_Lh4.ne4_ne4.F2010-SCREAMv1.eamxx-output-preset-1",
349+
"SMS_Lh4.ne4_ne4.F2010-SCREAMv1.eamxx-output-preset-1--eamxx-fixer_debug_output",
350350
"SMS_Lh4.ne4pg2_ne4pg2.F2010-SCREAMv1.eamxx-output-preset-1--eamxx-prod",
351351
)
352352
},

components/eamxx/cime_config/namelist_defaults_eamxx.xml

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -193,6 +193,10 @@ be lost if SCREAM_HACK_XML is not enabled.
193193
<moisture>moist</moisture>
194194
<!-- Frequency in physics steps to output a global hash over the dycore's
195195
in-fields. <= 0 disables hashing. -->
196+
<enable_column_conservation_checks>false</enable_column_conservation_checks>
197+
<enable_energy_fixer type="logical" doc="Turn on energy fixer for dycore (fixes dycore, pressure adjustment, PD coupling)">true</enable_energy_fixer>
198+
<!-- true means 2 more calls to AllReduce after dycore -->
199+
<enable_energy_fixer_debug_info type="logical" doc="Prints debug info that energy fixer indeed fixed energy (cost of an MPI reduce)">false</enable_energy_fixer_debug_info>
196200
<bfb_hash type="integer">18</bfb_hash>
197201
</homme>
198202

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,2 @@
1+
$CIMEROOT/../components/eamxx/scripts/atmchange homme::enable_energy_fixer_debug_info=true -b
2+
./xmlchange POSTRUN_SCRIPT="$CIMEROOT/../components/eamxx/scripts/check-fixer-output"

components/eamxx/cmake/tpls/CsmShare.cmake

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -35,6 +35,8 @@ macro (CreateCsmShareTarget)
3535
${SCREAM_BASE_DIR}/../../share/util/shr_orb_mod.F90
3636
${SCREAM_BASE_DIR}/../../share/util/shr_strconvert_mod.F90
3737
${SCREAM_BASE_DIR}/../../share/util/shr_sys_mod.F90
38+
${SCREAM_BASE_DIR}/../../share/util/shr_reprosum_mod.F90
39+
${SCREAM_BASE_DIR}/../../share/util/shr_reprosumx86.c
3840
)
3941
# Process genf90 template files. This adds a custom command (and hence target) for each f90 source
4042
# that needs to be built from the genf90 template files listed in GENF90_SOURCE.
@@ -57,7 +59,8 @@ macro (CreateCsmShareTarget)
5759
$<$<AND:$<COMPILE_LANGUAGE:Fortran>,$<CXX_COMPILER_ID:GNU>>:CPRGNU>
5860
$<$<AND:$<COMPILE_LANGUAGE:Fortran>,$<CXX_COMPILER_ID:Intel>>:CPRINTEL>
5961
$<$<AND:$<COMPILE_LANGUAGE:Fortran>,$<CXX_COMPILER_ID:Clang>>:CPRCRAY>
60-
$<$<AND:$<COMPILE_LANGUAGE:Fortran>,$<CXX_COMPILER_ID:clang>>:CPRCRAY>)
62+
$<$<AND:$<COMPILE_LANGUAGE:Fortran>,$<CXX_COMPILER_ID:clang>>:CPRCRAY>
63+
EAMXX_STANDALONE)
6164

6265

6366
if (${CMAKE_SYSTEM} MATCHES "Linux")
Lines changed: 156 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,156 @@
1+
#!/usr/bin/env python3
2+
3+
#"""
4+
#aaa
5+
#"""
6+
7+
import sys, re, glob, pathlib, argparse, gzip
8+
9+
from utils import run_cmd_no_fail, expect, GoodFormatter
10+
11+
###############################################################################
12+
def parse_command_line(args, description):
13+
###############################################################################
14+
parser = argparse.ArgumentParser(
15+
usage="""\n{0} <CASE_DIR> [<param>=<val>] ...
16+
OR
17+
{0} --help
18+
19+
\033[1mEXAMPLES:\033[0m
20+
\033[1;32m# Run hash checker on /my/case/dir \033[0m
21+
> {0} /my/case/dir
22+
""".format(pathlib.Path(args[0]).name),
23+
description=description,
24+
formatter_class=GoodFormatter
25+
)
26+
27+
parser.add_argument(
28+
"case_dir",
29+
help="The test case you want to check"
30+
)
31+
32+
return parser.parse_args(args[1:])
33+
34+
###############################################################################
35+
def readall(fn):
36+
###############################################################################
37+
with open(fn,'r') as f:
38+
txt = f.read()
39+
return txt
40+
41+
###############################################################################
42+
def greptxt(pattern, txt):
43+
###############################################################################
44+
return re.findall('(?:' + pattern + ').*', txt, flags=re.MULTILINE)
45+
46+
###############################################################################
47+
def grep(pattern, fn):
48+
###############################################################################
49+
txt = readall(fn)
50+
return greptxt(pattern, txt)
51+
52+
###############################################################################
53+
def get_log_glob_from_atm_modelio(case_dir):
54+
###############################################################################
55+
filename = case_dir / 'CaseDocs' / 'atm_modelio.nml'
56+
ln = grep('diro = ', filename)[0]
57+
run_dir = pathlib.Path(ln.split()[2].split('"')[1])
58+
ln = grep('logfile = ', filename)[0]
59+
atm_log_fn = ln.split()[2].split('"')[1]
60+
return str(run_dir / '**' / f'atm.log.*')
61+
62+
###############################################################################
63+
64+
# EAMxx:: energy fixer: T tend added to each physics midlevel 0.000222
65+
# EAMxx:: energy fixer: total energy before fix 32486509719.486938
66+
# EAMxx:: energy fixer: rel energy error after fix -7.63290383866757e-18
67+
68+
def get_fixer_lines(fn,start_from_line):
69+
###############################################################################
70+
fixer_lines = []
71+
error_vals = []
72+
73+
lines = []
74+
with gzip.open(fn,'rt') as file:
75+
start_line_found = False
76+
for line in file:
77+
if start_line_found:
78+
lines.append(line)
79+
elif start_from_line in line:
80+
start_line_found = True
81+
82+
i = 0
83+
while i < len(lines):
84+
line = lines[i]
85+
i = i+1
86+
# eamxx hash line has the form "eamxx hash> date=YYYY-MM-DD-XXXXX (STRING), naccum=INT
87+
# The INT at the end says how many of the following line contain hashes for this proc-step
88+
if "EAMxx:: energy fixer: rel energy error" in line:
89+
fixer_lines.append(line)
90+
lline = (line.strip('\n'))
91+
errr = float(lline.split()[8])
92+
error_vals.append(errr)
93+
94+
return fixer_lines,error_vals
95+
96+
###############################################################################
97+
def get_model_start_of_step_lines (atm_log):
98+
###############################################################################
99+
lines = []
100+
with gzip.open(atm_log,'rt') as file:
101+
for line in file:
102+
if "model beg-of-step time" in line:
103+
lines.append(line)
104+
return lines
105+
106+
###############################################################################
107+
def check_fixer_output(case_dir):
108+
###############################################################################
109+
110+
####################### TOLERANCE
111+
TOL = 1e-12
112+
113+
case_dir_p = pathlib.Path(case_dir)
114+
expect(case_dir_p.is_dir(), f"{case_dir} is not a dir")
115+
116+
# Look for the two atm.log files.
117+
glob_pat = get_log_glob_from_atm_modelio(case_dir_p)
118+
atm_fns = glob.glob(glob_pat, recursive=True)
119+
if len(atm_fns) == 0:
120+
print('Could not find atm.log files with glob string {}'.format(glob_pat))
121+
return False
122+
if len(atm_fns) == 1:
123+
print('Found output file {}'.format(atm_fns[0]))
124+
125+
start_line = get_model_start_of_step_lines(atm_fns[0])[0]
126+
127+
print (start_line)
128+
129+
# Extract hash lines, along with their timestamps, but ignore anything
130+
# before the line $start_line
131+
f = atm_fns[0]
132+
fixer_lines,errs = get_fixer_lines(f,start_line)
133+
print('Array of rel. energy errors after energy fixer:',errs)
134+
errsa = [abs(ele) for ele in errs]
135+
mmax = max(errs)
136+
print('Abs. max. for the rel. energy errors:',mmax)
137+
138+
if mmax < TOL:
139+
print('SUCCESS')
140+
return True
141+
else:
142+
print('FAIL, abs. max is less than tolerance, which is ', TOL)
143+
return False
144+
145+
146+
147+
###############################################################################
148+
def _main_func(description):
149+
###############################################################################
150+
success = check_fixer_output(**vars(parse_command_line(sys.argv, description)))
151+
sys.exit(0 if success else 1)
152+
153+
###############################################################################
154+
155+
if (__name__ == "__main__"):
156+
_main_func(__doc__)

components/eamxx/src/control/atmosphere_driver.cpp

Lines changed: 3 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -12,7 +12,7 @@
1212
#include "share/util/eamxx_timing.hpp"
1313
#include "share/util/eamxx_utils.hpp"
1414
#include "share/io/eamxx_io_utils.hpp"
15-
#include "share/property_checks/mass_and_energy_column_conservation_check.hpp"
15+
#include "share/property_checks/mass_and_energy_conservation_check.hpp"
1616

1717
#include <ekat_assert.hpp>
1818
#include <ekat_string_utils.hpp>
@@ -420,13 +420,12 @@ void AtmosphereDriver::setup_column_conservation_checks ()
420420
{
421421
// Query m_atm_process_group if any process enables the conservation check,
422422
// and if not, return before creating and passing the check.
423-
if (not m_atm_process_group->are_column_conservation_checks_enabled()) {
423+
if (not m_atm_process_group->are_conservation_checks_enabled()) {
424424
return;
425425
}
426426

427427
auto phys_grid = m_grids_manager->get_grid("physics");
428428
const auto phys_grid_name = phys_grid->name();
429-
430429
// Get fields needed to run the mass and energy conservation checks. Require that
431430
// all fields exist.
432431
EKAT_REQUIRE_MSG (
@@ -467,7 +466,7 @@ void AtmosphereDriver::setup_column_conservation_checks ()
467466
const auto heat_flux = m_field_mgr->get_field("heat_flux", phys_grid_name);
468467

469468
auto conservation_check =
470-
std::make_shared<MassAndEnergyColumnConservationCheck>(phys_grid,
469+
std::make_shared<MassAndEnergyConservationCheck>(m_atm_comm,phys_grid,
471470
mass_error_tol, energy_error_tol,
472471
pseudo_density, ps, phis,
473472
horiz_winds, T_mid, qv,

components/eamxx/src/dynamics/homme/eamxx_homme_process_interface.cpp

Lines changed: 40 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -252,7 +252,18 @@ void HommeDynamics::set_grids (const std::shared_ptr<const GridsManager> grids_m
252252

253253
// Create separate remapper for initial_conditions
254254
m_ic_remapper = grids_manager->create_remapper(m_cgll_grid,m_dyn_grid);
255-
}
255+
256+
// Layout for 2D (1d horiz X 1d vertical) variable
257+
const auto& grid_name = m_phys_grid->name();
258+
// Boundary flux fields for energy and mass conservation checks
259+
if (has_energy_fixer()) {
260+
add_field<Computed>("vapor_flux", pg_scalar2d, kg/(m2*s), grid_name);
261+
add_field<Computed>("water_flux", pg_scalar2d, m/s, grid_name);
262+
add_field<Computed>("ice_flux", pg_scalar2d, m/s, grid_name);
263+
add_field<Computed>("heat_flux", pg_scalar2d, W/m2, grid_name);
264+
}
265+
266+
}//set_grids
256267

257268
size_t HommeDynamics::requested_buffer_size_in_bytes() const
258269
{
@@ -466,7 +477,8 @@ void HommeDynamics::initialize_impl (const RunType run_type)
466477

467478
// Initialize Rayleigh friction variables
468479
rayleigh_friction_init();
469-
}
480+
481+
}//initialize_impl
470482

471483
void HommeDynamics::run_impl (const double dt)
472484
{
@@ -744,11 +756,21 @@ void HommeDynamics::homme_post_process (const double dt) {
744756
// Store T at end of the dyn timestep (to back out tendencies later)
745757
T_prev(ilev) = T_val;
746758
});
747-
});
759+
}); //op()
748760

749761
// Apply Rayleigh friction to update temperature and horiz_winds
750762
rayleigh_friction_apply(dt);
751-
}
763+
764+
if (has_energy_fixer()) {
765+
766+
get_field_out("vapor_flux").deep_copy(0);
767+
get_field_out("ice_flux").deep_copy(0);
768+
get_field_out("water_flux").deep_copy(0);
769+
get_field_out("heat_flux").deep_copy(0);
770+
771+
}; //if fixer
772+
773+
}//homme_post_proc
752774

753775
void HommeDynamics::
754776
create_helper_field (const std::string& name,
@@ -1045,6 +1067,9 @@ void HommeDynamics::restart_homme_state () {
10451067

10461068
// Erase also qv_prev_phys (if we created it).
10471069
m_helper_fields.erase("qv_prev_phys");
1070+
1071+
// Update the time stamp of the fields we inited in here (to avoid triggering invalid output in IO)
1072+
get_field_out("pseudo_density",pgn).get_header().get_tracking().update_time_stamp(start_of_step_ts());
10481073
}
10491074

10501075
void HommeDynamics::initialize_homme_state () {
@@ -1221,6 +1246,17 @@ void HommeDynamics::initialize_homme_state () {
12211246

12221247
// Can clean up the IC remapper now.
12231248
m_ic_remapper = nullptr;
1249+
1250+
// Update the time stamp of the fields we inited in here (to avoid triggering invalid output in IO)
1251+
get_field_out("pseudo_density",rgn).get_header().get_tracking().update_time_stamp(start_of_step_ts());
1252+
get_internal_field("v_dyn").get_header().get_tracking().update_time_stamp(start_of_step_ts());
1253+
get_internal_field("dp3d_dyn").get_header().get_tracking().update_time_stamp(start_of_step_ts());
1254+
get_internal_field("ps_dyn").get_header().get_tracking().update_time_stamp(start_of_step_ts());
1255+
get_internal_field("phis_dyn").get_header().get_tracking().update_time_stamp(start_of_step_ts());
1256+
get_internal_field("vtheta_dp_dyn").get_header().get_tracking().update_time_stamp(start_of_step_ts());
1257+
if (not params.theta_hydrostatic_mode) {
1258+
get_internal_field("w_int_dyn").get_header().get_tracking().update_time_stamp(start_of_step_ts());
1259+
}
12241260
}
12251261
// =========================================================================================
12261262
void HommeDynamics::

components/eamxx/src/share/CMakeLists.txt

Lines changed: 7 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -34,14 +34,15 @@ set(SHARE_SRC
3434
property_checks/property_check.cpp
3535
property_checks/field_nan_check.cpp
3636
property_checks/field_within_interval_check.cpp
37-
property_checks/mass_and_energy_column_conservation_check.cpp
37+
property_checks/mass_and_energy_conservation_check.cpp
3838
util/eamxx_data_interpolation.cpp
3939
util/eamxx_fv_phys_rrtmgp_active_gases_workaround.cpp
4040
util/eamxx_time_interpolation.cpp
4141
util/eamxx_time_stamp.cpp
4242
util/eamxx_timing.cpp
4343
util/eamxx_utils.cpp
4444
util/eamxx_bfbhash.cpp
45+
util/eamxx_repro_sum_mod.F90
4546
)
4647

4748
# Append ETI sources (I didn't do it above for clarity of reading)
@@ -123,6 +124,9 @@ if (EAMXX_ENABLE_EXPERIMENTAL_CODE)
123124
endif()
124125

125126
add_library(scream_share ${SHARE_SRC})
127+
set_target_properties(scream_share PROPERTIES Fortran_MODULE_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/modules)
128+
target_include_directories (scream_share PRIVATE ${CMAKE_CURRENT_BINARY_DIR}/modules)
129+
126130
target_include_directories(scream_share PUBLIC
127131
${SCREAM_SRC_DIR}
128132
${SCREAM_BIN_DIR}/src
@@ -132,6 +136,8 @@ target_include_directories(scream_share PUBLIC
132136
# Used in the data interpolation
133137
target_link_libraries(scream_share PUBLIC stdc++fs)
134138

139+
target_link_libraries(scream_share PUBLIC csm_share)
140+
135141
if (GPTL_PATH)
136142
target_include_directories(scream_share PUBLIC ${GPTL_PATH})
137143
else()

0 commit comments

Comments
 (0)