Skip to content

Use the real mesh geometry for plain Elmer mesh grids - #193

Merged
BenjaminRodenberg merged 5 commits into
mainfrom
elmer-mesh-real-coordinates
Oct 26, 2026
Merged

BenjaminRodenberg merged 5 commits into
mainfrom
elmer-mesh-real-coordinates

Conversation

@BenjaminRodenberg

@BenjaminRodenberg BenjaminRodenberg commented Oct 20, 2026 •

Copy link
Copy Markdown
Collaborator

Grids built from --elmer-mesh alone (no --netcdf-dem-mesh / --netcdf-dem-mesh-unstructured) filled the column geometry with placeholder values from the initial implementation: every column at 75 N 40 W (lon 320), a random elevation between 0 and 100 m, and flat surfaces. This PR replaces them with the geometry of the mesh, one commit each:

  1. lat/lon: cell centers of the mesh, shared with the DEM branches (_set_cell_center_coordinates). Longitudes are now -180..180, as for all other grid types. 320 and -40 give the same hour angle modulo 360 deg; only the shading azimuth tests h < 0, and Elmer grids cannot shade.
  2. z: mean of the node z coordinates of each cell, which is where reader.py writes DEM elevations. A mesh without elevation data gives 0 m.
  3. slopes: plane through the nodes of each triangle, rotated from projected x/y to east/north (meridian convergence, scale factor). slope_beta/slope_gamma use the formulas of the MATLAB branch, moved unchanged into _slope_angles.
  4. Elmer/Ice fallback: when coupled to Elmer/Ice, slopes are reset to flat columns with a warning instead of aborting the run.

Where results change, and where they do not

Measured with --dump-reference --random-seed 0 and tools/compare_snapshots.py, default run (1 Jan 1979, 8 x 3 h), per commit against its parent.

Configuration lat/lon z slopes
--matlab-mesh, --netcdf-mesh unchanged unchanged bit-identical (refactor only)
--elmer-mesh + --netcdf-dem-mesh(-unstructured) unchanged (already used cell centers) unchanged unchanged (still flat)
--elmer-mesh alone, mesh without elevation (greenland_mesh_v0) changes: Tsurf up to 0.06 K changes: Tsurf +1.0 K on average (random 0-100 m -> 0 m) unchanged (flat)
same, DEM-amended mesh (reader.py) changes changes strongly: Tsurf -7.6 K on average, up to 20 K changes: Tsurf > 0.1 K in 1 % of columns, up to 4 K
CI pseudo-coupling step (greenland_mesh_v0, --couple-to-icon-atmo --couple-to-icon-land --couple-to-elmer-ice --fake-coupling) bit-identical bit-identical bit-identical

Notes:

  • The CI step does not change because ICON-Land (fake) provides the surface energy balance, through which lat/lon act; EBFM's own energy balance only goes to the ebm_* diagnostics there. A standalone --elmer-mesh run on the same mesh does change.
  • The January runs are mostly in polar night, so the lat/lon and slope effects are much smaller than in summer.
  • The output file of plain --elmer-mesh runs now carries the real lon/lat and z.

Behaviour with Elmer/Ice coupling

When coupled to Elmer/Ice (--couple-to-elmer-ice), the columns stay flat and a warning is logged, since the slopes would not follow the surface elevation received from Elmer/Ice (it does not send dhdx/dhdy yet). Such runs therefore get the same flat columns as before. ElmerIce.validate_grid is renamed to ElmerIce.prepare_grid because it now adapts the grid; shading and 2-D elevation fields are still rejected.

Not covered

  • The --netcdf-dem-mesh(-unstructured) branches still use flat slopes; their elevations are interpolated to cell centers only, so there are no node elevations to fit a plane through.
  • slope_gamma has the wrong sign where slope_y == 0 exactly; this is pre-existing, moved unchanged, and fixed separately.

Generated with Claude Code

@BenjaminRodenberg BenjaminRodenberg self-assigned this Oct 20, 2026
@BenjaminRodenberg BenjaminRodenberg added the bug Something isn't working label Oct 20, 2026
@BenjaminRodenberg
BenjaminRodenberg force-pushed the elmer-mesh-real-coordinates branch from 6239bab to e3d753c Compare October 20, 2026 15:19
@BenjaminRodenberg BenjaminRodenberg changed the title Place Elmer mesh columns at their cell centers Use the real mesh geometry for plain Elmer mesh grids Oct 20, 2026
@BenjaminRodenberg BenjaminRodenberg added this to the v0.9.0 milestone Oct 20, 2026
@BenjaminRodenberg
BenjaminRodenberg force-pushed the elmer-mesh-real-coordinates branch 2 times, most recently from 3eaa683 to e053416 Compare October 26, 2026 15:38
Plain --elmer-mesh grids (no DEM) put every column at lat 75, lon 320
(75 N 40 W), placeholder values from the initial implementation. Use
the mesh cell centers instead, as the --netcdf-dem-mesh branch already
did; both branches now share _set_cell_center_coordinates.

Longitudes are now -180..180 instead of 0..360. Both give the same
hour angle modulo 360 deg, so cos/sin(h) are unchanged; only the
shading azimuth's h < 0 test differs, and Elmer grids cannot shade.

On greenland_mesh_v0 (default 1-day January run) Tsurf changes by up
to 0.06 K. The CI pseudo-coupling step is bit-identical because its
surface energy balance comes from (fake) ICON-Land.

Assisted-by: Claude Code:claude-opus-5-5
Plain --elmer-mesh grids (no DEM file) filled grid["z"] with random
values between 0 and 100 m. Use the mean of the node z coordinates of
each cell instead, which is where reader.py stores DEM elevations and
what README.md describes for this input. A mesh without elevation data
now gives 0 m everywhere.

Assisted-by: Claude Code:claude-opus-5-5
Plain --elmer-mesh grids treated every column as flat. Use the plane
through the nodes of each triangle instead. Its gradient in projected
x/y is rotated to local east/north with the meridian convergence of
the mesh CRS and scaled by the projection scale factor, since
slope_x/slope_y are east/north gradients and EPSG:3413/3031 x/y are
not. Mesh now keeps its CRS for this.

slope_beta/slope_gamma use the same formulas as the MATLAB branch
(arctan of the slope, slope_azimuth).

A mesh without elevation data, such as greenland_mesh_v0 in CI, still
gives exactly zero slopes. With elevation data, --couple-to-elmer-ice
now aborts at startup in ElmerIce.validate_grid.

Assisted-by: Claude Code:claude-opus-5-5
Elmer mesh grids now carry real slopes, which ElmerIce rejected since
they cannot follow the surface elevation it sends back (no dhdx/dhdy
yet). That turned every Elmer/Ice run on a mesh with elevation data
into an AssertionError. Zero the slopes with a warning instead, which
is what these runs used before the mesh slopes were introduced.

Shading and 2-D elevation fields are still rejected. validate_grid is
renamed to prepare_grid since it now modifies the grid.

Assisted-by: Claude Code:claude-opus-5-5
@BenjaminRodenberg
BenjaminRodenberg force-pushed the elmer-mesh-real-coordinates branch from e053416 to 9f82868 Compare October 26, 2026 16:26
@BenjaminRodenberg
BenjaminRodenberg merged commit e9381e7 into main Oct 26, 2026
34 checks passed
@BenjaminRodenberg
BenjaminRodenberg deleted the elmer-mesh-real-coordinates branch October 26, 2026 17:28
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

bug Something isn't working

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant