Skip to content

Deep audit: DEDX_AUTO returns garbage/unusable results for 174 materials, element-boundary off-by-one, 2 API-reachable memory-safety bugs, plus API/doc inconsistencies — with a plan of action #149

Description

@grzanka

Deep audit of the library sources, the public headers, and the user-facing docs at d5bf0cd. The findings below are all reproduced against a build of main, most of them with a self-contained C program pasted inline. Where an existing issue already covers a finding I reference it instead of restating it.

Method

  • full read of src/*.c, include/*.h, python/, tools/, docs/, index.rst, README.md
  • build with -Wall -Wextra, and a second build with -fsanitize=address,undefined
  • an exhaustive sweep of 101 957 (program, ion, material) triples — every combination that dedx_get_material_list_for_ion() advertises, for every program except DEDX_ESTAR — loading each one and probing the energy range

The sweep alone turns up three systematic problems:

check result
advertised combos that fail to load 407
advertised combos where dedx_get_min/max_energy() is rejected by the loaded dataset 1 568
DEDX_AUTO configs that load with err == DEDX_OK but return ENERGY_OUT_OF_RANGE at every energy 174 (plus 996 more partly broken)

ctest is green (32/32) under ASan+UBSan for the whole suite, so none of this is caught today.


A. Correctness — wrong or unusable numbers

A1 🔴 load_compound() mixes energy grids across tiers, and sums uninitialised heap

src/dedx.c:783-851. load_compound() calls find_data() once per constituent element, reusing one energy[] buffer and then doing the Bragg sum with:

for (j = 0; j < compound_data[0].length; j++) {      /* dedx.c:842 — length of the FIRST element */
    data.data[j] = 0.0;
    for (i = 0; i < length; i++)
        data.data[j] += weight[i] * compound_data[i].data[j];
}
...
return load_data(ws, &data, energy, config, err);    /* energy[] = whatever the LAST element left */

Two independent assumptions are unchecked:

  1. all constituents share one energy grid — the loop length comes from element 0, but the x-axis comes from element n−1;
  2. every stopping_data is filled up to that lengthload_bethe_2() (dedx.c:853-896) writes only data->data[0..121] and never zeroes the rest of the 150-float array, and compound_data comes from a plain malloc().

Both are violated by DEDX_AUTO, which is precisely the mode designed to mix tiers: tabulated elements come back on the 133-point PSTAR-family grid, Bethe-fallback elements on the 122-point Bethe grid (dedx_bethe_energy is exactly the first 122 points of dedx_pstar_energy, so lengths differ but values silently overlap).

A1a — dead configs. When the first constituent is tabulated (133 pts) and the last is Bethe-only, read_embedded_energy_data() has memset energy[122..149] = 0, so the spline knots end at x = 0:

=== BLOOD_ICRP (id 118) === composition Z order: 1 6 7 8 11 12 14 15 16 17 19 20 26 30
load OK (err=0), n=133 bragg_used=1
grid tail: x[119]=900 x[120]=950 x[121]=1000 x[122]=0 ... x[132]=0
check_energy_bounds() uses low=0.001 high=0  -> every query fails
dedx_get_stp(100 MeV) = 0 err=101

dedx_load_config() returns 0 / DEDX_OK, and then every dedx_get_stp() call fails. For protons this hits 46 of 278 advertised materials, including TISSUE_SOFT_ICRP, BLOOD_ICRP, BRAIN_ICRP, LUNG_ICRP, SKIN_ICRP, TESTES_ICRP, CALCIUMCARBONATE, CALCIUMOXIDE, GLASS_PLATE, SARAN, FREON_*, SODIUMNITRATE, TRICHLOROETHYLENE. Across ions 1–18: 174 dead, 996 partly broken, 3 834 ok.

A1b — uninitialised heap in the result. When the last constituent is tabulated the grid survives, but data.data[122..132] is still Σ wᵢ · (uninitialised bytes). The same call then returns different physics depending on what was on the heap:

/* poison the exact block load_compound() is about to malloc, then query */
static void poison(unsigned char b, int n) {
    size_t sz = sizeof(stopping_data) * (size_t)n;      /* 612 * n */
    unsigned char *p = malloc(sz); memset(p, b, sz); free(p);
}
poison(pat, 3);
dedx_get_simple_stp_for_program(DEDX_AUTO, 1, DEDX_CADMIUM_TUNGSTATE, 2000.0f, &err);
heap fill     stp(2000 MeV)  stp(5000 MeV)   stp(100 MeV)
0x00               0.903953       0.946743        4.36013
0x3c               0.903953       0.946743        4.36013
0x41               0.903953        4.71554        4.36013
0x7f            1.05969e+38    1.05969e+38        4.36013
0xbe               0.787708       0.830498        4.36013

err == DEDX_OK every time. The physically correct value is ≈ 0.9 MeV cm²/g.

A1c — silent range truncation. Reversed ordering (Bethe element first) truncates the table to 122 points, so BORONCARBIDE under DEDX_AUTO silently tops out at 1 000 MeV while dedx_get_max_energy(DEDX_AUTO, 1) still advertises 10 000.

Fix. load_compound() must carry the energy grid per constituent, verify all constituents resolved to the same grid (identical length and identical knots), and fail with a real error otherwise — or resample onto a common grid. Independently, load_bethe_2() should memset(data, 0, sizeof(*data)) like read_embedded_stopping_data() already does, and load_compound() should use min(length_i) rather than length[0].

A2 🔴 Element/compound boundary is off by one: material 99 is a compound, not an element

dedx_elements.h runs DEDX_HYDROGEN = 1 … DEDX_CALIFORNIUM = 98, then DEDX_A150_TISSUE_EQUIVALENT_PLASTIC = 99. The header comment gets it right ("Elemental ions (Z=1–98)"), and dedx_embedded_metadata.c:41 gets it right (if (id < 99)), but the rest of the library uses <= 99 / > 99 as the element/compound test:

site code
src/dedx.c:316 if (material <= 99) return element_supported_for_ion(...)
src/dedx.c:753 if (prog == DEDX_AUTO && target_load > 0 && target_load <= 99) → Bethe fallback
src/dedx.c:858 if (config->target > 99) { *err = COMBINATION_NOT_FOUND; }
src/dedx_validate.c:90 if (config->target > 0 && config->target <= 99) return 0; → skip Bragg decomposition
src/dedx_embedded_metadata.c:101 state == 2 && !is_gas(target) && target <= 99 → ×1.13 solid-state I correction

So A150 tissue-equivalent plastic is evaluated by the Bethe path as elemental einsteinium (TZ = 99, TA = dedx_amu[98] = 252 u), and its mean excitation potential gets the elemental solid-state correction:

PSTAR   p 100 MeV in A150 (tabulated) = 7.3376
BETHE   p 100 MeV in A150             = 5.1799   <-- Z/A of Es-252, not of A150
BETHE   p 100 MeV in ACETONE (id 100) = 7.3905   (decomposed correctly)
BETHE   p 100 MeV in WATER            = 7.3009

loader i_value for A150 (id 99)     = 73.5630 eV  bragg_used=0   <-- 65.1 x 1.13
loader i_value for ACETONE (id 100) = 64.2000 eV  bragg_used=1

Composition data for id 99 exists and is correct (H .1013, C .7755, N .0351, O .0523, F .0174, Ca .0184) — it is simply never used. dedx_get_material_list_for_ion(DEDX_BETHE_EXT00/DEDX_AUTO, …) advertises A150 as available, so the wrong number is served without any signal.

Fix. Introduce one named constant (e.g. DEDX_MAX_ELEMENT_ID 98) and replace every 99 literal. A regression test asserting bragg_used == 1 and Bethe ≈ tabulated within a few % for id 99 would have caught this.

A3 🟠 dedx_get_i_value() hardcodes the gas state, and the compound loader uses it

src/dedx.c:224:

float dedx_get_i_value(int target, int *err) {
    return dedx_internal_get_i_value(target, DEDX_GAS, err);   /* state is never a parameter */
}

DEDX_GAS suppresses the ×1.13 condensed correction, so the public accessor disagrees with what the library itself uses:

dedx_get_i_value(CARBON)               = 70.000 eV
loader-resolved cfg->i_value(CARBON)   = 81.000 eV   (compound_state = DEDX_CONDENSED)

Worse, load_compound() (dedx.c:818, 827) calls the public accessor to seed _temp_i_value for each constituent, so Bragg-decomposed compounds are computed with gas-phase elemental I-values while dedx_internal_evaluate_i_pot() (dedx_validate.c:42, 58) uses config->compound_state. Two code paths, two answers, for the same element.

Fix. Add a state parameter (dedx_get_i_value_for_state()), or at minimum make load_compound() call dedx_internal_get_i_value(z, config->compound_state, err). The header should state which state the value refers to.

A4 🟠 dedx_get_min_energy() / dedx_get_max_energy() advertise bounds the data rejects

1 568 of the swept combos are rejected at an endpoint the accessors advertise. Breakdown:

prog  name   ion  combos  boundmismatch  advertised
9     ICRU     6     114            110  [8.33333e-05, 833.333]
10    AUTO     1     279            138  [0.001, 10000]
10    AUTO     2     279            130  [0.00025, 250]
10    AUTO   3-18    279         61 each  [0.025, 1000]
10    AUTO     6     279            275  [8.33333e-05, 833.333]

Two causes: the DEDX_CARBON special case (dedx.c:406, 448) only holds for the handful of targets ICRU90_C covers, and DEDX_AUTO's Bethe tier has a different range from the tabulated tier it mirrors. The code comment at dedx.c:400-405 calls these "best-effort hints, not authoritative", but dedx.h:257-269 documents them as "the minimum/maximum valid energy for a program/ion combination" with no caveat — and dedx_get_inverse_csda() / dedx_get_csda() (dedx_tools.c:106-107, 189) use them as hard integration limits. There is also no default: in either switch, so an unknown program silently returns 0.

Fix. Either make the bounds target-aware (dedx_get_energy_range(program, ion, target, *lo, *hi)), or derive them from the loaded dataset (dedx_get_config_energy_range(ws, cfg, …)) and document the program-level ones as advisory. Add a default: returning a sentinel.

A5 🟠 A missing density row blocks tabulated programs — 407 false positives in the availability API

dedx_internal_validate_rho() (dedx_validate.c:19-26) runs unconditionally for every program and propagates DEDX_ERR_TARGET_NOT_FOUND when the target has no embedded density — even though ρ is irrelevant to a table lookup. Exactly one material is affected (FERROUSOXIDE, id 159), but it kills that material for all 407 program/ion pairs:

PSTAR  ion=1 target=159: rc=-1 err=201
       ... with cfg->rho=5.7 : rc=0 err=0    <-- the table was there all along
ICRU73 ion=6 target=159: rc=-1 err=201
       ... with cfg->rho=5.7 : rc=0 err=0

material_id_supported() (dedx.c:306) does not model this, so dedx_get_material_list_for_ion() — documented as "a material appears here if and only if dedx_load_config() is expected to succeed" (dedx.h:222-230) — is wrong for all of them.

Fix. Only require ρ when the program actually needs it (program >= DEDX_DEFAULT, or DEDX_AUTO when it falls through to Bethe). Add the missing ferrous-oxide density row. Add a CI test that loads every advertised combination.

A6 🟠 A zero in the embedded ICRU73 table is served as a valid result and silently downgrades interpolation

Exactly one non-positive value exists in the whole ICRU73 table — ion Z=11 (Na), target 18 (Ar), first energy point:

ICRU73  Na in Ar:  E=0.025 -> stp=0        err=0
                   E=0.03  -> stp=2671     err=0

Two consequences:

  • dedx_get_stp() returns 0.0 with DEDX_OK; dedx_get_csda()/dedx_get_inverse_stp() then hit 1.0/stp (dedx_tools.c:57, 66 guard it to INFINITY, which quietly poisons the quadrature).
  • dedx_internal_calculate_coefficients() (dedx_spline.c:101-105) reacts to any non-positive sample by silently switching the entire table to linear interpolation, discarding the caller's DEDX_INTERPOLATION_LOG_LOG request. loaded_data->interpolation_mode still reports LOG_LOG, so there is no way to detect it.

Fix. Validate the generated tables at build time (tools/dat2c.py should refuse or flag non-positive stopping values), decide what the missing point means (drop it, or extrapolate), and record the effective interpolation mode so the silent downgrade is observable.

A7 🟠 An invalid mstar_mode produces an all-zero table with no error

dedx_internal_calculate_mspaul_coef() (dedx_mpaul.c:134-137) falls through to an empty else { /* illegal mode */ }, leaving FOUT = 0, so output = 0:

mstar_mode default : stp(100) = 262.6096  err=0
mstar_mode 'Z'     : load rc=0 err=0
mstar_mode 'Z'     : stp(100) = 0.0000    err=0

Nothing validates config->mstar_mode against the documented DEDX_MSTAR_MODE_* set.

Fix. Validate mstar_mode in dedx_internal_validate_config(); return a new DEDX_ERR_INVALID_MSTAR_MODE (or reuse DEDX_ERR_ION_NOT_SUPPORTED_MSTAR, currently marked "reserved legacy code") from the mpaul fall-through.

A8 🟡 check_ion() is skipped on the custom-compound path

dedx_load_config() (dedx.c:471-474) dispatches straight to load_compound() when elements_id != NULL, and only load_config_clean() calls check_ion() (dedx.c:681). A carbon ion with DEDX_PSTAR (proton-only) and a custom compound gets DEDX_ERR_COMBINATION_NOT_FOUND (202) instead of DEDX_ERR_ION_NOT_SUPPORTED (207). Move the check into dedx_internal_validate_config().


B. Memory safety

B1 🔴 NULL dereference in load_compound() from a plain public-API call

dedx.c:845 reads weight[i] where weight = config->elements_mass_fraction. dedx_internal_validate_config() only runs dedx_internal_evaluate_compound() (which derives the weights) when program >= 100 or target == 0 (dedx_validate.c:202, 221). A tabulated program + a non-zero target + an element list therefore reaches the sum with weight == NULL:

cfg->program = DEDX_PSTAR;   cfg->ion = DEDX_PROTON;   cfg->target = DEDX_WATER;
cfg->elements_length = 2;
cfg->elements_id = malloc(2 * sizeof(int));  cfg->elements_id[0] = 1;  cfg->elements_id[1] = 8;
/* elements_atoms and elements_mass_fraction left NULL — permitted by the header docs */
dedx_load_config(ws, cfg, &err);
dedx.c:845:35: runtime error: load of null pointer of type 'float'
AddressSanitizer: SEGV on unknown address 0x000000000000
    #0 load_compound  src/dedx.c:845
    #1 dedx_load_config src/dedx.c:472

B2 🔴 Out-of-bounds read for elements_id[i] <= 0

dedx_internal_get_atom_mass() / get_nucleon() (dedx_periodic_table.c:12-26) guard only the upper bound (if (id < 113)), so dedx_amu[id - 1] indexes negatively. elements_id[] is caller-supplied and never range-checked:

custom compound with element Z=0  -> dedx_periodic_table.c:15: index -1 out of bounds for 'float [112]'
custom compound with element Z=-5 -> dedx_periodic_table.c:15: index -6 out of bounds for 'float [112]'

dedx_get_nucleon_number() / dedx_get_atom_mass() in dedx.c:230-246 already guard ion < 1 and even comment on this; the guard belongs in the periodic-table functions themselves.

B3 🟠 dedx_get_ion_list(DEDX_ESTAR) hands out 1001, which overflows dedx_ion_table[120][40]

Distinct from #110 (which is about callers passing an invalid id): here the library's own list is the source of the bad id. dedx_program_available_ions[3] = {1001, -1} — the documented way to enumerate ions. The shipped CLI does exactly that:

$ ./getdedx 3 -1 276 100
src/dedx.c:205:26: runtime error: index 1001 out of bounds for type 'char [120][40]'
ESTAR can handle the following ions:
1001: ?6B?R(B?BB9B??B???AB`?Ax??A?u?A+v?A?A4

And getdedx -1 CARBON WATER 100 — a documented invocation (index.rst:188: "Passing -1 for the program, ion, or target slot lists the available values") — is a hard ASan crash, because dedx_get_ion_list(-1) indexes dedx_program_available_ions[-1]:

src/dedx.c:370:43: runtime error: index -1 out of bounds for type 'int [110][20]'
AddressSanitizer: global-buffer-overflow ... READ of size 4 at ...
    #0 main examples/dedx_get.c:75
0x... is located 0 bytes after global variable 'dedx_available_programs'

CI only exercises getdedx -1 1 276 100 (numeric ion), which avoids the lookup. Same unguarded indexing in dedx_get_material_list() (dedx.c:266) and dedx_get_ion_list() (dedx.c:370) — please extend #110 to cover the list accessors and the ESTAR row, and fix examples/dedx_get.c to resolve names only after the program id is known good.

B4 🟡 dedx_get_error_code() is an unbounded strcpy into a caller buffer of undocumented size

dedx.c:111-190. The longest message is 55 characters ("Embedded elemental composition metadata is unavailable."), so callers need ≥ 56 bytes — nowhere documented (dedx.h:113-117 just says "caller-allocated"). Add #define DEDX_ERROR_STRING_MAX 64, document it, and consider a bounded variant.

B5 🟡 Unchecked calloc() in load_bethe_2()

dedx.c:885-886 allocates the Bethe workspace and dereferences it immediately at dedx.c:888 with no NULL check — the only unchecked allocation left in the library.


C. API contracts and error handling

C1 🟠 dedx_get_stp() cannot tell a valid config from a stale or foreign one

dedx.c:486-493 validates only 0 <= cfg_id < ws->active_datasets. It never checks config->loaded, and dedx_config carries no back-pointer to its workspace. Both of these return a plausible wrong number with err == DEDX_OK:

stp(wsA, cfgA, 100) = 7.2861   (water, correct)
stp(wsB, cfgA, 100) = 3.6271   <-- cfgA evaluated against a different workspace, err=0

cfgB->cfg_id=0 loaded=0 -> stp=7.2861 err=0   <-- never loaded; calloc's 0 is a valid id

The second is nasty: calloc()-zeroed cfg_id == 0 is exactly the id of the first loaded dataset, and dedx.h:318 tells users to calloc() the struct.

Fix. Check config->loaded; initialise cfg_id to -1 in a dedx_config_init() helper (or store a workspace tag/generation counter and compare).

C2 🟠 dedx_load_config() is not idempotent and there is no unload

Every successful load increments active_datasets; nothing ever releases a slot. Re-loading the same config on a 1-slot workspace leaves the config unusable:

1st load: rc=0  err=0   cfg_id=0  loaded=1
2nd load: rc=-1 err=203 cfg_id=-1 loaded=0

dedx_get_inverse_stp() (dedx_tools.c:137) calls dedx_load_config() unconditionally, so it works exactly once per workspace slot, while dedx_get_csda() (dedx_tools.c:181) guards on config->loaded:

1st: E=66.0138 err=0
2nd: E=-1      err=203

Fix. Make dedx_load_config() reuse the existing slot when config->loaded and the workspace matches; add dedx_unload_config(); make the two tools functions agree.

C3 🟠 dedx_get_csda() / dedx_get_inverse_stp() read *err as an input

dedx_tools.c:135, 140, 179, 183 do if (*err != 0) return -1; before ever writing *err, and neither sets DEDX_OK on success. With an uninitialised int err this is UB; with a reused variable it is a silent wrong answer:

dedx_get_csda with stale err=101 -> -1        (err=101)   <-- correct answer is 7.72
dedx_get_csda with err=0         -> 7.72118   (err=0)

Fix. *err = DEDX_OK; at entry of every function taking an int *err, as dedx_internal_validate_config() already does (dedx_validate.c:190).

C4 🟡 convert_units(): the second switch's default discards the accumulated factor

dedx_tools.c:228-230 assigns rather than validating, so an unknown target unit silently yields an identity conversion instead of an error; and a material with no density gives conversion_rate = 0 (or a divide-by-zero when converting to mass stopping power), with new_values[] written before the error is returned:

MeVcm2/g -> unit 99                        : rc=0   out={1.0000, 2.0000}   (silently 1:1)
MeVcm2/g -> MeV/cm with unknown material   : rc=201 out={0.0000, 0.0000}

Reject unknown units, bail before writing on error. (Naming is already tracked in #112.)

C5 🟡 DEDX_ERR_INCONSISTENT_COMPOUND (211) has no message

It is the only code in dedx_error.h missing from the dedx_get_error_code() switch:

dedx_get_error_code(211) -> "No such error code."
dedx_get_error_code(212) -> "Interpolation mode is not supported."

A table-driven mapping plus a test that walks every DEDX_ERR_* macro would keep the two files in sync permanently.

C6 🟡 config->i_value is silently ignored by tabulated programs

dedx.h:332 documents i_value as an optional input defaulting to "tabulated ICRU values", with no hint that it only affects Bethe-type programs:

PSTAR water, I = default -> 7.2861
PSTAR water, I = 200 eV  -> 7.2861   (identical)

Either reject a caller-supplied i_value for programs that cannot honour it, or document the restriction. (Related: #6.)

C7 🟡 find_data() ignores the energy-read error

dedx.c:765 calls dedx_internal_read_energy_data() without checking *err; on the MSTAR path dedx_internal_convert_energy_to_mstar() then resets *err = DEDX_OK (dedx_mstar.c:52), erasing it, and proceeds over an untouched stack buffer. Not reachable today (embedded data is always present) but it is the same class of bug as A1.


D. Documentation

D1 🔴 The list accessors are documented as 0-terminated; they are -1-terminated

dedx.h:202-211 and dedx.h:251-255:

"Return a null-terminated list … Pointer to a static array terminated by 0" — for dedx_get_program_list(), dedx_get_material_list() and dedx_get_ion_list().

dedx_get_program_list() = 1 2 3 4 5 6 7 9 10 100 101 -1

A consumer following the header runs off the end of every one of these arrays. The implementation comments (dedx.c:257, 262, 363), dedx_wrappers.h:45-48, and dedx.h:231-233 all correctly say -1, so the three @return lines are simply stale. This should be fixed first — it is a one-line doc change that prevents an OOB read in user code.

D2 🟠 Wrapper return contracts do not match the implementation

function header says actually returns
dedx_get_stp_table_size() "Number of data points, or 0 if not supported" (dedx_wrappers.h:85) -1 on failure
dedx_fill_default_energy_stp_table() "Number of points filled, or negative error code" (dedx_wrappers.h:96) 0 on success

python/libdedx/_api.py:97-108 codes against the actual behaviour, so the docs are the odd one out — but get_default_table() also has an unreachable if n == 0 branch because the C side never returns 0.

D3 🟠 DEDX_AUTO is absent from all user-facing documentation

grep -n AUTO README.md index.rst examples/README.md docs/zero hits. The program tables in README.md:12-22 and index.rst:52-90 list eight programs and omit both DEDX_AUTO and DEDX_DEFAULT, even though dedx.h:52-63 gives DEDX_AUTO a 12-line rationale and it is the mode with the most severe bug in this report (A1). Tracked as #145 — raising priority given A1.

D4 🟡 Undocumented buffer contracts

  • dedx_get_composition() (dedx.h:155-161) takes float composition[][2] with no stated capacity; internal callers hard-code [20][2] and the largest real composition is 14 rows. Publish DEDX_MAX_COMPOSITION_ELEMENTS and take a capacity argument.
  • dedx_fill_program_list() / dedx_fill_material_list() / dedx_fill_ion_list() (dedx_wrappers.h:18-35) have no capacity parameter at all — the caller must guess. dedx_fill_material_list() needs up to 281 ints.
  • dedx_get_error_code() — see B4.

D5 🟡 dedx_get_min/max_energy() documented as authoritative but implemented as advisory — see A4.


E. Architecture and build hygiene

E1 🟠 No warning flags anywhere in the build

Neither CMakeLists.txt nor src/CMakeLists.txt sets -Wall/-Wextra//W4. Building by hand surfaces 8 live -Wsign-compare warnings (dedx.c:79, 842, 887; dedx_validate.c:56, 74, 120, 252, 270) — all int i against unsigned int loop bounds, several of them in the very loops involved in A1. Recommend -Wall -Wextra -Wshadow -Wconversion on the library target, -Werror in CI.

E2 🟠 CI has no sanitizer coverage

ci.yml runs Valgrind on one binary (test_bethe_ext00). Valgrind cannot see global-buffer-overflows on the static tables (B3), and memcheck stays quiet on A1b because the uninitialised bytes never reach a branch. An ASan+UBSan ctest job would have caught B1, B2 and B3 immediately, and a "load every advertised combination" test would have caught A1, A2 and A5.

E3 🟡 The Bethe evaluator is duplicated

evaluate_bethe_model() (dedx_bethe.c:249-315) and evaluate_bethe_model_LEext() (dedx_bethe.c:75-143) are ~60 identical lines apart from the Lindhard-Scharff prologue. Any physics change has to be made twice. Extract the shared core. (Related: #95.) The three while (1) golden-section loops (dedx_bethe.c:154, 185, 219) also have no iteration cap and the int *err threaded through them is never read.

E4 🟡 Embedded metadata has no generator — provenance is lost

tools/dat2c.py regenerates only the stopping-power tables. src/data/embedded/dedx_composition.h (compositions) and the density / I-value / gas-state / effective-charge rows in dedx_metadata.h have no generator: grep -rn "gas_states\|effective_charge\|composition" tools/ → no hits, even though data/raw/gas_states.dat and data/raw/effective_charge.dat are retained "as the remaining raw metadata inputs" (data/README.md:31-33). They are now hand-maintained C blobs, which is why A5's missing ferrous-oxide density cannot be fixed at the source. Add a generator + a build-time consistency check (every id in dedx_material_table has a density row; every compound has a composition; no non-positive stopping values — see A6).

E5 🟡 Dead and internally inconsistent embedded data

src/data/embedded/dedx_estar.h, dedx_icru90_e.h, dedx_icru90_pos.h are generated but included by no translation unit. dedx_estar.h is also self-inconsistent — dedx_estar_energy[133] vs dedx_estar_stp[1][76][132] — so wiring it into the generic dedx_embedded_program_data struct (which derives energy_len from the energy array) would compute strides of 133 over a 132-wide table and read out of bounds. Fix the generator or the data before #5 (add ESTAR tables) is picked up. (Adjacent: #126.)

E6 🟡 Availability matrices are hand-maintained

Already acknowledged in dedx_program_const.h:63-65 and tracked by #138. Worth noting that dedx_program_available_materials[3] (ESTAR) is the only row that includes material id 0 ("(N/A)"), and that the ESTAR ion row is the source of B3.

E7 — already tracked


Plan of action

Sequenced so that each phase leaves the tree green and the next phase cheaper.

Phase 0 — stop the bleeding (docs + guards only, no behaviour change)

  1. D1 — fix the three @return lines to say -1. One-line change, prevents OOB reads in downstream code. Ship immediately.
  2. B2, B4, B5 — add the lower-bound guard in dedx_periodic_table.c, DEDX_ERROR_STRING_MAX, and the calloc NULL check.
  3. B1 — reject elements_id != NULL without weights in dedx_internal_validate_config() (DEDX_ERR_INCONSISTENT_COMPOUND).
  4. C5 — table-driven error strings + a test that walks every DEDX_ERR_*.

Phase 1 — make the failures visible

  1. Add -Wall -Wextra -Werror to the library target and clear the 8 sign-compare warnings (E1).
  2. Add an ASan+UBSan ctest job to ci.yml, and extend Valgrind to the whole suite (E2).
  3. Add tests/test_availability_exhaustive.c: for every combination dedx_get_material_list_for_ion() advertises, assert dedx_load_config() succeeds, the advertised energy bounds are accepted, and every sampled value is finite and > 0. This is the regression net for A1/A4/A5/A6 — it currently reports 407 + 1 568 + 174 failures.

Phase 2 — the correctness fixes

  1. A1 (highest severity). Per-constituent energy grids in load_compound(); reject or resample mismatched grids; memset in load_bethe_2(); use min(length_i).
  2. A2. DEDX_MAX_ELEMENT_ID 98; replace every <= 99 / > 99; regression test on id 99.
  3. A5. Require ρ only for programs that use it; add the ferrous-oxide density row; teach material_id_supported() about it.
  4. A3. State-aware I-value accessor; make load_compound() use compound_state.
  5. A7, A8. Validate mstar_mode; move check_ion() into dedx_internal_validate_config().
  6. A6. Build-time data validation in tools/dat2c.py; decide the fate of the ICRU73 Na→Ar zero; expose the effective interpolation mode.

Phase 3 — API contracts

  1. C1, C2. dedx_config_init() with cfg_id = -1; check config->loaded in dedx_get_stp(); workspace tagging; dedx_unload_config(); idempotent dedx_load_config().
  2. C3. *err = DEDX_OK at entry everywhere.
  3. A4. Target-aware / dataset-derived energy-range API; default: in both switches.
  4. C4, C6. convert_units() validation; decide and document i_value semantics.
  5. B3. Bounds-check the list accessors, fix the ESTAR ion row, fix examples/dedx_get.c.

Phase 4 — architecture and data

  1. E4. Generator + build-time consistency check for compositions and metadata.
  2. E3. De-duplicate the Bethe evaluator; cap the golden-section loops.
  3. E5. Remove or fix the unused ESTAR/ICRU90-e/ICRU90-pos headers.
  4. D2, D3, D4. Wrapper return contracts, DEDX_AUTO/DEDX_DEFAULT documentation (docs: document DEDX_DEFAULT and DEDX_AUTO programs, and document available materials #145), buffer-capacity constants.

Suggested split: A1 and A2 each deserve their own bug issue and their own PR — they change numbers users may already depend on and need release notes. Everything in Phase 0 can go in one small PR. I am happy to open the sub-issues and start on Phase 0 + Phase 1 if that ordering looks right.


Reproducers (self-contained, build against a normal libdedx.a)
/* A1a — DEDX_AUTO loads OK, then fails at every energy.  Needs -Isrc for the internal header. */
#include <dedx.h>
#include <stdio.h>
#include <stdlib.h>
#include "dedx_lookup_data.h"

int main(void) {
    int err = 0;
    dedx_workspace *ws = dedx_allocate_workspace(1, &err);
    dedx_config *cfg = calloc(1, sizeof(dedx_config));
    cfg->program = DEDX_AUTO; cfg->ion = 1; cfg->target = DEDX_TISSUE_SOFT_ICRP;
    err = 0;
    printf("load rc=%d err=%d\n", dedx_load_config(ws, cfg, &err), err);
    dedx_internal_lookup_data *d = ws->loaded_data[cfg->cfg_id];
    printf("n=%d  x[121]=%g  x[122]=%g  x[%d]=%g\n", d->n, d->base[121].x, d->base[122].x,
           d->n - 1, d->base[d->n - 1].x);
    int e = 0;
    float v = dedx_get_stp(ws, cfg, 100.0f, &e);
    printf("stp(100 MeV) = %g err=%d\n", v, e);
    return 0;
}
/* load rc=0 err=0 ; n=133 x[121]=1000 x[122]=0 x[132]=0 ; stp(100 MeV) = 0 err=101 */
/* A2 — material id 99 (A150 plastic) evaluated as elemental einsteinium */
#include <dedx.h>
#include <dedx_wrappers.h>
#include <stdio.h>
int main(void) {
    int err = 0;
    printf("name(99)=%s  name(98)=%s\n", dedx_get_material_name(99), dedx_get_material_name(98));
    printf("PSTAR p 100 MeV in A150 = %.4f\n",
           dedx_get_simple_stp_for_program(DEDX_PSTAR, 1, 99, 100.0f, &err));
    printf("BETHE p 100 MeV in A150 = %.4f   <-- Z=99, A=252\n",
           dedx_get_simple_stp_for_program(DEDX_BETHE_EXT00, 1, 99, 100.0f, &err));
    printf("BETHE p 100 MeV in ACETONE(100) = %.4f\n",
           dedx_get_simple_stp_for_program(DEDX_BETHE_EXT00, 1, 100, 100.0f, &err));
    return 0;
}
/* A150_TISSUE_EQUIVALENT_PLASTIC  CALIFORNIUM ; 7.3376 ; 5.1799 ; 7.3905 */
/* B1 — NULL deref;  B2 — negative index.  Build with -fsanitize=address,undefined */
#include <dedx.h>
#include <stdlib.h>
int main(void) {
    int err = 0;
    dedx_workspace *ws = dedx_allocate_workspace(1, &err);
    dedx_config *cfg = calloc(1, sizeof(dedx_config));
    cfg->program = DEDX_PSTAR; cfg->ion = 1; cfg->target = DEDX_WATER;
    cfg->elements_length = 2;
    cfg->elements_id = malloc(2 * sizeof(int));
    cfg->elements_id[0] = 1; cfg->elements_id[1] = 8;   /* no weights -> SEGV at dedx.c:845 */
    return dedx_load_config(ws, cfg, &err);
}

The exhaustive sweep (407 / 1 568 / 174) is a short program over
dedx_get_program_list()dedx_get_ion_list()dedx_get_material_list_for_ion();
I will contribute it as tests/test_availability_exhaustive.c in Phase 1.

Metadata

Metadata

Assignees

No one assigned

    Labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions