A rounding step in AERIUS sub-receptor averaging and the verdict at 0.00499999 mol/ha/yr

Work in progress (web edition of 30 October 2026). This is a working draft under review, published here so that it can be read and checked. It may change, and its findings are not final. It is version 0.3 of StikstofNu research paper D. The English text is the only edition; a Dutch edition follows (Nederlandse versie volgt). Where the draft says it is not published, that describes it before this web edition.

What differs from the draft. References to files of the research materials are shortened to the file or folder name; those materials are not published. Nothing else is changed except what the change list of this edition names; in this paper the rate of Paper 1 is no longer quoted per hexagon, because its denominator is simulated draws (Section 1). Corrections and questions: info@alexandermooij.com.

What the code does, what a Monte Carlo on rescaled OPS output for illustrative sources shows, and what the flip rate does and does not mean. Abstract: a rounding step in AERIUS sub-receptor averaging can move the NOx + NH3 total of a hexagon across 0.00499999 mol/ha/yr, in both directions. The per-key rounding is confirmed in the source and in the compiled jar; the steps after it (sector-group sums, proposed minus reference) are not modelled. The magnitude is conditional and weaker. The table of verdict changes by exact total (Table 3) is less dependent on the window than one rate, but it is still conditional on the synthetic fields (one illustrative source type, an NH3 share uniform from 0 to 1, one OPS build), so we do not call it a general result. We give no single rate: the same flips give 1.58 per 10,000 draws when exact totals are drawn between 0.004 and 0.006 mol/ha/yr and 12.4 per 10,000 when they are drawn between 0.0049 and 0.0051 (Table 4). For the window 0.004 to 0.006, five reproducible regenerated runs give 1.35 to 1.76 per 10,000 (the smallest and largest run, not an interval; pooled 1,236 of 7,800,000, 1.585), and the original committed run gave 225 in 1,560,000 (1.44). Intervals for the pooled rate are in Section 5.4. The rate is conditional on a synthetic distribution and is not a share of real hexagons. Working draft v0.3, 26 October 2026. English edition.

Independence. This is independent desk research by StikstofNu. It is not affiliated with AERIUS, RIVM, the Ministry, a province or a consultant, and none of them reviewed it. Nothing here replaces an official AERIUS Calculator result. The paper makes no policy point: it asks only what holds up under scrutiny. Evidence labels: [V] verified by opening the source; [RUN] computed from committed files or regenerated inputs; [INF] inference from code; [U] unverified.

The answer first

1. The mechanism is real and is in the code (V, RUN). When a hexagon lies within 300 m of a source, AERIUS 2026.0.0 replaces it by 96 or 384 sub-points and takes the mean of their OPS results. The mean is computed as BigDecimal.valueOf(sum).divide(BigDecimal.valueOf(n), RoundingMode.HALF_UP) (the code names the two variables value and sumCount). That quotient keeps the decimal scale of Double.toString(sum). When the double sum happens to have a short decimal form (3 to 6 decimals), the mean is rounded to that many decimals, which is far coarser than the 0.00499999 bound. The defect is also present in the 2025.3 and 2025 releases.
2. Short decimal sums, and the rounding they cause, are common (RUN). In 4,680 real hexagon results from the production OPS binary, 25.3 % of the sums have at most 6 decimals, and 6.5 % of the means are off by more than 1e-4 relative (maximum 1.14 %). None of the 4,680 changed a verdict, and none lay within 0.0000088 of the bound; that run is too small to say anything about the flip rate (Section 5.1).
3. Both directions exist, but only for the decision variable AERIUS tests (RUN, V). For a single substance, a sub-receptor mean can be wrongly counted (true mean below the bound, rounded up to 0.005) and, in a few rare cases, wrongly lost. AERIUS tests the NOx + NH3 total, which is the double sum of two separately rounded means. When the two sums have different decimal scales, the total can land on either side of the bound. We show a worked case in each direction (Section 5.3), and an exhaustive enumeration (Section 5.2) shows that, in the scale pairs and numbers of points enumerated, same-scale pairs flip only upward and mixed-scale pairs flip both ways.
4. The flip probability by exact total depends less on the window than one rate, but it is conditional on the synthetic fields (RUN). In five reproducible regenerated runs (7,800,000 draws) no verdict changed for an exact NOx + NH3 total below 0.0045. From 0.0045 up to the bound the share of draws that are counted although they are not above the bound is 8.7e-5 (0.0045 to 0.00495), 8.7e-4 (to 0.004995), 1.0e-2 (to 0.0049995) and 1.1e-1 (to 0.00499999). Above the bound the share of draws that are lost is 1.1e-3 up to 0.0051 and 5.0e-5 from 0.0051 to 0.006 (Table 3). A single rate per 10,000 depends on the window of exact totals it is averaged over. The same 1,236 flips give 1.58 per 10,000 over 0.004 to 0.006, 3.12 over 0.0045 to 0.0055, 4.48 over 0.0045 to 0.0051, 6.82 over 0.0048 to 0.0052 and 12.4 over 0.0049 to 0.0051 (Table 4). For the window 0.004 to 0.006 (log-uniform, NH3 share uniform) the five runs give 1.35 to 1.76 per 10,000, which is the smallest and largest run and not an interval. Pooled, it is 1,236 of 7,800,000 (1.585 per 10,000), with a Poisson 95 % interval of 1.50 to 1.68 and cluster bootstrap intervals of 1.41 to 1.91 at the widest (Section 5.4). The original committed run gave 225 in 1,560,000 (1.44). It cannot be reproduced bit for bit, because its input dump was not committed and its line order depended on thread timing (Section 7). A sixth regenerated run (unsorted dump) has no committed input and is not pooled.
5. What the figure is not (U). It is a conditional rate: the probability of a verdict change given that a hexagon's exact NOx + NH3 total lies in the chosen window, drawn log-uniformly, with an NH3 share drawn uniformly, on 390 hexagons of one synthetic source type. It is not the share of real permit hexagons that are misjudged. The real distribution of hexagon totals near 0.005 is not known to us, and the real run cannot test the rate: only 7 of 14,040 real NOx and NH3 pairs have an exact total between 0.004 and 0.006, so about 0.001 flips would be expected there and the zero observed says nothing (Section 5.1). We therefore cannot state how many real projects are affected.
6. Severity is not rated. A flip needs a hexagon whose exact total is close to the bound: none occurred below 0.0045, and in the narrowest band below the bound (0.0049995 to 0.00499999) about one draw in ten flips. Within a band, the effect is the difference between a hexagon that counts and one that does not. We do not rate what that means for real permits. The bound itself is a policy and tool parameter; no source we found prescribes a rounding rule for it (Section 2.3).

1. The claim, and where it comes from

The claim as put to this paper: a rounding defect in sub-receptor averaging can flip the verdict near the 0.00499999 mol/ha/yr bound, in both directions. The earlier wording gave a rate of about 1.4 per 10,000 and stated it per hexagon; its denominator is simulated draws, not hexagons (the original committed run gave 225 in 1,560,000 draws, 1.44 per 10,000 draws; Sections 5.4 and 6).

It is established in the repository in four places [V]:

• Paper 1, Finding F9 and the answer to RQ5: in a Monte Carlo on rescaled real OPS fields, the verdict on the NOx + NH3 total changes in both directions, at a rate of about 1.4 per 10,000 for exact totals between 0.004 and 0.006 mol/ha/y, stated there per hexagon (the denominator is simulated draws: 225 in 1,560,000) (P1-aerius-deconstructed.md, Section 9, F9 and Answer to RQ5).
• The working note f9-and-version-drift.md, Section 2.6, with the table by exact total and the committed output sum_near_summary.txt.
• The evidence register of the supplement, entries E54 and E56 (P1-supplement.md).
• The harness and Monte Carlo sources in f9/ (F9Harness.java, F9SumMonteCarlo.java).

What Paper 1 says, exactly. The line numbers refer to P1-aerius-deconstructed.md as it stands on main, which is version 4 of 22 October 2026. The public English edition is version 4.1 of 24 October 2026; version 4.1 corrects other statements and carries both wordings below unchanged (line 13 and, in version 4.1, line 585). The abstract (line 13) says that the defect changes the verdict on the NOx + NH3 total near the bound, in both directions, at a rate of about 1.4 per 10,000 stated per hexagon, after "in a Monte Carlo on real OPS sub-point fields with illustrative sources". Finding F9 and the answer to RQ5 (line 568) give the same rate, stated per hexagon, "with a total between 0.004 and 0.006 mol/ha/y, in both directions". Only the second wording names the window of exact totals; neither gives an interval, the number of runs or the dependence on the window. The 1.4 is the summed Monte Carlo of the committed run (225 in 1,560,000 is 1.44).

The history matters. Version 3 of Paper 1 claimed "2 in 10,000" from a single-substance Monte Carlo and called the flips one-directional. The interim review (interim-review-r1-v3.md, must-fix 2) pointed out that AERIUS decides on the NOx + NH3 total, not on one substance. The summed case was then run, which gave 1.44 per 10,000 draws (225 in 1,560,000), and showed that the flips go both ways. This paper checks that chain again from the files and from a fresh run.

As worded, the claim reads like a share of hexagons. This paper tests a narrower statement: the probability of a verdict change per draw, given a window of exact totals (Sections 5.4 and 6). That probability depends strongly on the window.

2. Mechanism

2.1 The code path

All code references are to the AERIUS monorepo at fb8812f ("Release 2026.0.0", 2026-10-06) unless noted [V].

• Sub-receptors. A receptor hexagon whose centre is within 150 m of a source with emission is replaced by 384 sub-points, within 300 m by 96. The counts are those of levels 4 and 3 of HexagonEquilateralTrianglesMap, asserted in source/commons-aerius/aerius-commons-util/src/test/java/nl/overheid/aerius/receptor/HexagonEquilateralTrianglesMapTest.java:87-93; the distances of 150 m and 300 m, the first level (3) and the 20 m rule are in source/commons-aerius/aerius-commons-domain/src/main/java/nl/overheid/aerius/shared/domain/theme/own2000/OwN2000SubReceptorConstants.java:29-49. Points more than 20 m from the source are kept (code), and according to the Handboek only those over relevant habitat. The hexagon result is the arithmetic mean of the kept points. This is documented in the user handbook (Handboek 2026, Section 3.6, 384 and 96 sub-receptors and the arithmetic mean, pdf p. 33 [V]; see also legal-and-documentation.md, row d).
• The division. AbstractSubReceptorAggregator.SumResult.setSum (source/commons-aerius/aerius-commons-util/src/main/java/nl/overheid/aerius/receptor/AbstractSubReceptorAggregator.java:203-209) computes, per result key, BigDecimal.valueOf(value).divide(BigDecimal.valueOf(sumCount), RoundingMode.HALF_UP).doubleValue(), where value is the double sum of the sub-point results accumulated in point order by add and setResults (:175-201) and sumCount the number of points that had a result.
• Why this rounds. BigDecimal.valueOf(double) uses Double.toString, which gives the shortest decimal that round-trips (JDK 19 and later; production runs JDK 21). BigDecimal.divide(BigDecimal, RoundingMode) returns a quotient with the scale of the dividend. So the mean is rounded to as many decimals as the sum has in its shortest form. A sum such as 0.4799 has 4 decimals, so its mean over 96 points is rounded to 4 decimals. A sum with 17 significant digits, which is the common case, gets about 16 to 18 decimals and no visible rounding.
• The decision variable. The NOx + NH3 deposition is formed on read as the double sum of the NOx and NH3 values (source/commons-aerius/aerius-commons-shared/src/main/java/nl/overheid/aerius/shared/domain/result/EmissionResults.java:88-103). It is also what is stored: when the results are inserted, determineEmissionResultKeys adds the NOx + NH3 key for every point that has a NOx or NH3 deposition and its value is read through that sum (source/commons-aerius/aerius-commons-repository/src/main/java/nl/overheid/aerius/db/calculator/CalculationRepository.java:562-566 and :576-583). Each of the two values has been rounded by its own division, with its own scale.
• The test. A hexagon counts when ABS(calculation_results.result) > v_permit_lower_bound (source/database/src/main/sql/calculator_common/calculator/10-public/03-functions/02-project_calculation.sql:341), and for the national data the bound is 0.00499999 (source/database/src/data/sql/nl-latest/permit_lower_bounds/load.sql:4, identical in nl-stable). The bound sits 1e-8 below the displayed label 0.005. (The integration-test fixture also carries a second row, policy 0.5 / 0.499999, which is not in nl-latest or nl-stable; the mechanism is the same at any bound.)
• Age. The same division is in release 2025.3 (71a6015, line 209 of the same file) and release 2025 (7bf619a, line 197) [V, GitLab raw files read on 2026-10-25]. No commit exists after fb8812f on any ref (GitLab API, 2026-10-25 17:55 Europe/Amsterdam), so the pin is still the current public release.

2.2 The mechanism in numbers

RoundingDemo.java calls the same BigDecimal operations on hand-picked sums (JDK 21; output in RoundingDemo.out.txt) [RUN].

Table 1. Single-substance means, n = 96.

A one-decimal sum can erase a mean as large as 0.05. A three- to six-decimal sum moves the mean to the next multiple of 0.001 to 0.000001. Only the means near the bound change a verdict.

2.3 Where the bound comes from, and what "wrongly" means

Quotations below were checked against the source files by d_quotes_check.py (d_quotes.json); the source files are not committed [V, fetched 2026-10-26].

• The tool is prescribed. The Raad van State refers to articles 4.15, 6.15 and 7.197h of the Omgevingsregeling for the mandatory use of AERIUS Calculator ("Zie de artikelen 4.15, 6.15 en 7.197h van de Omgevingsregeling.", voorlichting W11.25.00063/IV, 23 May 2025, footnote 19, pdf p. 7). We did not open the consolidated text of those articles; the article numbers rest on that footnote [U for the article text].
• The bound is a policy and tool parameter. The same voorlichting says: "uitgegaan van een ondergrens van 0,005 (afgerond 0,01) mol/ha/jaar. Deze grens is bewust zo laag mogelijk gehouden" (pdf p. 7). It reports that the Commissie Hordijk (2020) noted that "hoewel de onzekerheid in de berekening veel hoger is dan 0,005 mol/ha/jaar, deze waarde voor een beleidstoepassing desondanks nodig is" (pdf p. 8), and that the expert opinion of Petersen finds the 0.005 bound has "geen empirische basis" (pdf p. 10; "heeft volgens het expert-oordeel geen empirische basis"). The government described it on 18 March 2025 as "Op dit moment geldt een nagenoeg nulgrens (0,005 mol/ha/jaar)." (Rijksoverheid.nl news item). Backes and Petersen (21 January 2026, p. 1) describe current practice: Calculator "rondt momenteel af op een honderste (0,01) mol en vereenzelvigt (interne) berekeningsresultaten van minder dan 0,005 mol/ha/jaar met 0 mol/ha/jaar".
• The Handboek. Handboek 2026 (pdf p. 23) calls it "de voor de vergunning geldende rekengrens van 0,005 mol/ha/jaar" and says that a results GML always holds the raw results, "zodat de resultaten van combinatieberekeningen die Calculator uitvoert reproduceerbaar blijven". On pdf p. 112, results below the bound are shown as '-'.
• What we did not find. We found no statute, regulation or document that states the value 0.00499999 or a rounding rule for internal sub-receptor means. The value 0.00499999 is in the code and the database load script (Section 2.1). Whether the bound is set in a legal text is a question for the author (Section 8).
• What "wrongly" means here. "Counted wrongly" and "lost" are judged against the verdict on the exact arithmetic mean of the printed OPS sub-point values, which is the quantity the Handboek describes (arithmetic mean over the sub-receptors, pdf p. 33). It is not a breach of a stated rounding rule, because none was found. That the Handboek keeps raw results unrounded for reproducibility (pdf p. 23) suggests, but does not show, that rounding of internal sums to 3 to 6 decimals is not intended [INF].

3. Data and computation behind the figure

All of this is from f9/ and f9-and-version-drift.md, Section 1 [V for the description, RUN for the reproduction in Section 4].

• Harness. F9Harness.java drives the real AERIUS 2026.0.0 OPS worker classes (OPSSubReceptorCreator, OPSProcessRunner, the BRN, RCP and control-file writers, ResultDataFileReader, OPSSubReceptorAggregator) with the production OPS package 5.3.1.0-5 (amd64 binary, sha256 prefix 992a8170a3e6b897). Per hexagon it records the AERIUS value (output of the real aggregator), a re-implementation of the division, the exact decimal mean of the printed values, and the plain double mean.
• Sources and hexagons. Six 2 x 2 km windows in Natura 2000 areas, 17 random point-source positions, 390 hexagons within 300 m of a source (291 at 96 points and 99 at 384 points), habitat geometry from the public open-data WFS. Point source, height 4 m, spread 2 m, no heat content, NOx 1 to 500 kg/yr and NH3 0.2 to 100 kg/yr. These characteristics are illustrative, not the closed production characteristics (closed dbdata).
• Real run. 204 OPS runs gave 4,680 hexagon results (main.csv).
• Monte Carlo. The 780 printed sub-point fields (390 hexagons times NOx 100 kg and NH3 5 kg) are rescaled: per sample a total mean T is drawn log-uniformly from 0.004 to 0.006 and an NH3 share f uniformly from 0 to 1; the NOx field is scaled to mean T(1 - f) and the NH3 field to T f; every value is re-rounded to the four significant digits OPS prints; each substance is summed in double in point order and divided as AERIUS does; the two means are added. The exact total is the exact decimal sum divided by n. The shape of each field is fixed and only its level moves, while the digit pattern, which drives the effect, is regenerated per sample. 4,000 samples per hexagon give 1,560,000 samples per run.
• Not modelled. Totals over several sector groups (each its own rounded mean), proposed minus reference, the float8 to numeric conversion when results are stored, and production source characteristics.

4. Reproduction

The seeds of the new runs (1 to 4, plus the original 20261022 on the sorted and on the unsorted dump) were fixed before any of them was run. No new hypothesis was registered: these are replications of a computation that already existed.

5. Results

5.1 Real OPS output (4,680 hexagon results) [RUN]

• 1,183 of 4,680 sums (25.3 %) have at most 6 decimals. The rest have 11 to 18.
• The relative error of the AERIUS mean against the exact mean has median 3.5e-15, 99th percentile 6.8e-4 and maximum 1.14 % (sum 2.622 over 96 points gives 0.027 against 0.0273125). 6.5 % of the hexagons are off by more than 1e-4 relative.
• 499 hexagons have an exact mean between 0.0025 and 0.01; the largest absolute error among them is 3.1e-5.
• 0 of 4,680 verdicts differ at 0.00499999. The closest exact mean to the bound is 0.0050088 (a 17-decimal sum), so none lay in a band where a flip is possible. This run therefore shows the size of the rounding, but not a flip.
• The real run cannot test the flip rate. The real run tested single-substance verdicts, but AERIUS decides on the NOx + NH3 total. Pairing the real NOx and NH3 results of the same hexagon and sub-point count gives 14,040 pairs; only 7 have an exact total between 0.004 and 0.006 and 1 between 0.0045 and 0.0051. At the pooled rate of the Monte Carlo (1.585 per 10,000) the expected number of flips among the 7 pairs is about 0.001 (0.0004 when the band rates of Table 3 are used). Observing none is therefore uninformative. The same holds for the statement that no real result lay within 0.0000088 of the bound.

5.2 Exhaustive enumeration on the decimal grid [RUN]

No random numbers and no distribution assumption. A sum with s decimals is k x 10^-s with k not divisible by 10; AERIUS rounds its mean over n points HALF_UP to s decimals; the verdict is compared with the exact mean. We enumerate every such sum whose exact mean lies in 0.004 to 0.006, for n = 96 and n = 384 (d_grid_enumeration.json).

• One substance. For s = 3 to 7 only upward changes occur: the zone of exact means that are counted although they are not above the bound is [0.005 - 0.5 x 10^-s, 0.00499999] (observed minima 0.0045, 0.00495, 0.004995, 0.0049995 and 0.00499995 for s = 3 to 7). For s = 1 the mean is rounded to a multiple of 0.1 and is erased: for n = 96, 2 one-decimal sums have an exact mean in range and 1 of them is lost (exact mean 0.0052); for n = 384 it is 3 of 7 (exact means 0.0055 to 0.0060). (An earlier version of d_grid_enumeration.py admitted one or two grid points outside 0.004 to 0.006; it now tests the window with exact integers, which gives 2 and 1 for n = 96 and 7 and 3 for n = 384, not 4 and 2 and 9 and 4.) For s = 2 nothing changes within these n. For s = 8 and above, a thin zone just above the bound is lost (for s = 8, exact means from 0.0049999901 to 0.0049999949, which round to 0.00499999 and are not above it). That zone is 0.5 x 10^-s wide, so it is negligible for the usual 17-decimal sums.
• Two substances. With both sums at the same scale (3 and 3, or 4 and 4), in the scale pairs and numbers of points enumerated, the total can only be counted wrongly. With different scales (3 and 4, 4 and 5, 3 and 5) both directions occur. The enumeration finds 18,690 counted-wrongly and 14,808 lost pairs for scales 3 and 4 at n = 96, with exact totals in 0.0049 to 0.0051. These are counts of grid points, not rates.

5.3 A worked case in each direction

From the arithmetic (Java, RoundingDemo), n = 96:

Table 2.

From the committed Monte Carlo output (sum_near_flips.csv), draws on rescaled real OPS fields, not real projects [RUN, V]: a hexagon with n = 96 and an NH3 sum of 4 decimals had NOx mean 0.00275611458333333 and NH3 mean 0.0022; the AERIUS total 0.004956 is not counted, the exact total 0.0050040 is above the bound. A hexagon with n = 384 and sums of 6 and 16 decimals had NOx mean 0.002869 and NH3 mean 0.0021311; the AERIUS total 0.0050001 is counted, the exact total 0.0049999 is not.

5.4 The NOx + NH3 Monte Carlo [RUN]

Five regenerated runs of F9SumMonteCarlo on the sorted dump (seeds 20261022, 1, 2, 3 and 4; 1,560,000 draws each, 7,800,000 in all) are reproducible bit for bit from the committed dump. An independent rerun on 26 October 2026 gave byte-identical output files (hashes in SHA256SUMS-D-regenerated). A sixth run on the unsorted dump has no committed input and is reported separately. The window is the range of exact NOx + NH3 totals from which draws are taken: 0.004 to 0.006, log-uniform. The results below are computed by d_window_uncertainty.py (d_window_uncertainty.json) and cross-checked against the decimal-exact bands of d_mc_analysis.json.

Table 3. Verdict changes by exact total, five reproducible runs (7,800,000 draws).

The band edges are post hoc: they follow the decimal-scale zones of Section 5.2 and were not fixed before the runs. The rates describe these runs and are not tests of fixed boundaries.

This table depends less on the window than a single rate, but it is conditional on the synthetic fields: one illustrative source type, an NH3 share uniform from 0 to 1, 390 hexagons from 17 positions and one OPS build. The rate already differs by sub-point count and by source (Section 5.4), so the band rates are not a general property of AERIUS. Flips occur only for exact totals from 0.0045 upwards: below the bound only upward flips (counted wrongly), above it only losses. The rate rises steeply towards the bound (Table 3): 1.0e-2 from 0.004995 and 1.1e-1, about one draw in ten, in the last band, which is less than 5e-7 wide. The wrong counts do not all sit in that band: 194 of the 677 are in it and 159 are in [0.0045, 0.00495). Draws just outside 0.004 to 0.006 (136 of 7,800,000, after the four-digit re-rounding of the fields) are counted in the first and last band.

Table 4. The same flips, averaged over different windows of exact totals.

No flip is added or removed between the rows; only the denominator changes. A single rate per 10,000 is therefore a property of the window, and this paper does not offer one as a prevalence. The window 0.004 to 0.006 is the one used by the earlier work (Paper 1, Finding F9) and is arbitrary.

Table 5. The headline rate over runs (window 0.004 to 0.006).

• Pool and range. The pool is the five reproducible runs: 1,236 flips in 7,800,000 draws, 1.585 per 10,000. The range 1.35 to 1.76 is the smallest and largest of the five run rates. It is a spread of runs, not a confidence interval of any quantity. The unsorted-dump run (222 flips, 1.42) lies inside it; a pool of all six runs would be (1,458 of 9,360,000, 1.56) and is not used. The committed run (225 in 1,560,000, 1.44) lies inside the range as well but is not in the pool, because its input is not committed.
• Sampling intervals for the pooled rate. Poisson exact 95 % interval for 1,236 events: 1.50 to 1.68 per 10,000 (for the committed 225: 1.26 to 1.64). The five run counts have a standard deviation of 23.3 flips against 15.7 expected from Poisson noise (chi-square 8.8 on 4 degrees of freedom, p = 0.066), a variance ratio of about 2.2. With only five runs the test has little power, so extra scatter between runs is not ruled out and the Poisson interval may be too narrow; the cluster bootstrap intervals below are the ones to rely on. Cluster bootstrap 95 % intervals (10,000 resamples, seed 7), treating the fixed 390 hexagons as a sample: resampling hexagons 1.47 to 1.70, resampling the 17 sources 1.44 to 1.78, resampling the 6 sites 1.41 to 1.91. These cover the choice of hexagons, sources and sites under the model. They do not cover the structural uncertainty below.
• Structural uncertainty (not in any interval). The distribution of totals and the window (Table 4); one source type at one height; one OPS build; one aggregation key and no proposed-minus-reference step (Section 7). The rate also differs between hexagons: the 300 hexagons with at most 96 sub-points used (this group includes 384-point hexagons from which points were removed) give 1.85 per 10,000 and the 90 with more give 0.72; per-source rates run from 1.2 to 3.5; 33 of the 390 hexagons never flipped in any of the five runs.
• Both directions are real: 55 % of the pooled flips are counted wrongly and 45 % are lost.
• A flip needs a short component sum. In the five runs 10 of 1,236 flips (0.8 %) had no component with at most 6 decimals, while 70 % of all samples have no short component. 163 to 200 of the 390 hexagons (42 to 51 %) flipped at least once per run.
• Display window. In the committed run, 4 of 1,560,000 samples gave a counted total below the label 0.005 (shown as "-" in the PDF overview and 0,00 at two decimals), and 3 of them had an exact total above the bound. In the five runs, 43 of 7,800,000 samples (5.5e-6), and 26 of those had an exact total above the bound.

5.5 The single-substance variant [RUN]

The version 3 design (one substance, 2,000 samples per field, means from 0.004 to 0.006) gave 302 flips in 1,560,000 samples in the committed run, 297 counted wrongly and 5 lost (1.94 per 10,000). Three regenerated runs gave 278, 270 and 285 flips (pooled 833 of 4,680,000, 1.78 per 10,000), of which 823 were counted wrongly and 10 lost. The plain double mean never changed a verdict (0 in all runs). So for one substance the effect is almost entirely one-directional; the two-directional result needs the two-substance total.

6. What the flip rate means, and what it does not

• It is a conditional rate on a synthetic distribution and a chosen window. The condition is an exact total between 0.004 and 0.006, log-uniform, an NH3 share uniform over 0 to 1, 390 hexagons of one illustrative source type, 96 or 384 points (or fewer where habitat or the 20 m rule removes points). The same flips give 1.58 to 12.4 per 10,000 depending on the window (Table 4). It is not a prevalence among real hexagons.
• It does not say how many real projects are affected. That needs the real distribution of the highest hexagon total near 0.005 in permit calculations, which we do not have. The real run cannot inform it: 7 of 14,040 real pairs have an exact total in the window, and about 0.001 flips would be expected there (Section 5.1).
• The range is a spread of runs. The five reproducible runs span 1.35 to 1.76 per 10,000 and pool to 1.585; the Poisson and cluster bootstrap intervals are in Section 5.4, and none of them covers the structural uncertainty.
• Wording proposed for use. "Rounding before averaging can move the NOx + NH3 total of a hexagon across 0.00499999 mol/ha/yr, in both directions (confirmed in the AERIUS source). In a Monte Carlo on rescaled real OPS sub-point fields (390 hexagons, one illustrative source type) with exact totals drawn between 0.004 and 0.006 mol/ha/yr, the verdict changed in 1.35 to 1.76 of 10,000 draws (smallest and largest of five reproducible runs; pooled 1.585). The rate depends on the window: it is 12.4 per 10,000 when only totals between 0.0049 and 0.0051 are drawn. The original committed run gave 225 in 1,560,000 (1.44 per 10,000). It is not a share of real hexagons."

7. Limits and open items

• Illustrative sources. Heights, spreads and emission levels are not production values; other source types change the within-hexagon gradient and the digit patterns.
• One OPS build. The harness used the amd64 OPS binary. Production uses the aarch64 build, and in Paper 2 the two builds differ in the fourth printed digit for about 0.2 % of hexagons. Whether that also holds within 300 m of a source was not tested: there was no aarch64 run of this harness. Not running it is a choice, not a constraint: the aarch64 OPS binary is available on the box. Whether to run it is a decision for the author (Section 8).
• Sector groups. AERIUS rounds each sector-group mean and sums them; each rounding adds a chance to flip. Not modelled.
• Aggregation key and proposed minus reference. AERIUS aggregates per source key and substance and then sums, and the legal test applies to the project calculation (proposed minus reference). The Monte Carlo models one key. Rounding errors at different scales in the two situations need not cancel. The direction and size of this effect are not bounded here.
• Storage. The float8 to numeric conversion on insert and the subtraction of a reference situation are not modelled.
• Two runs have no committed input. The committed run (225 flips) and the sixth regenerated run (unsorted dump, 222 flips) cannot be rerun. The pool and the range use only the five runs on the committed sorted dump.
• The committed run is not bit-reproducible. The original dump of sub-point fields was not committed, and F9Harness writes it from several threads, so its line order, which sets which random draws a hexagon receives, differed in the regenerated dump. The regenerated dump is committed here, sorted, so later runs are deterministic. The rate itself reproduces within the range above.
• OPS-derived data. f9_dump_sorted.csv holds OPS printed values, like the committed main.csv. Paper 1, Section 7.3, raises the OPS terms of use. This file is withheld from any public release because of those terms, not for privacy reasons. It stays in the private repository, no derived file is published in its place, and it is listed in PUBLIC_RELEASE_WITHHELD.txt. This is not a promise to contact anyone.
• Not a statement about intent or fault. The code does what it says. Whether the rounding is acceptable is a question for the maintainers, and Paper 1 treats it as an ordinary bug report (Finding F9). No contact with the AERIUS maintainers has been made.

8. Open decisions

1. Wording of the rate. Section 6 proposes stating the rate as the range 1.35 to 1.76 per 10,000 (smallest and largest of five reproducible runs; pooled 1.585), with its dependence on the window, and with the committed 225 in 1,560,000 (1.44) as the original figure. The wording is the author's to confirm or revise.
2. Decided, not open: f9_dump_sorted.csv is withheld from public release (OPS terms) and only the summaries are published; no derived file replaces it.
3. Whether to state a legal reading of the 0.005 bound, and which sources to cite for it (Section 2.3 lists what was found and what was not).
4. Whether to contact the AERIUS maintainers before any publication.
5. Whether to run the production-like checks that would turn this from a conditional to a real rate: an aarch64 harness run, and a Connect job with a source within 300 m of habitat. Neither was done, and no Connect job was submitted for this paper.

9. Reproduction

All in rounding_subreceptor/ unless noted. Prerequisites as in f9-and-version-drift.md, Section 2.5 (JDK 21, AERIUS 2026.0.0 commons in the local Maven repository, OPS 5.3.1.0-5 with ops_terrain_properties_20250619).

cd f9 && mvn dependency:build-classpath -Dmdep.outputFile=cp.txt
javac -d classes -cp $(cat cp.txt) src/main/java/nl/overheid/aerius/ops/*.java
java -Xmx10g -Df9.dump=dump.csv -cp classes:$(cat cp.txt) nl.overheid.aerius.ops.F9Harness <opsroot> <rundir> \
  data/sources.csv data/hexes.csv dumprun.csv "NOX:100;NH3:5" 4
sort the dump by (site, pos, receptor, substance)     # as f9_dump_sorted.csv
java -cp classes:$(cat cp.txt) nl.overheid.aerius.ops.F9SumMonteCarlo f9_dump_sorted.csv out.csv 4000 0.004 0.006 <seed>
python3 -I d_mc_analysis.py name=out.csv      # prints JSON
python3 -I d_window_uncertainty.py <repo-root> name=out.csv ...   # five sorted-dump runs; windows, intervals, real-run power
python3 -I d_quotes_check.py <dir-with-sources>              # quotations of Section 2.3
python3 -I d_main_csv_check.py <repo-root>
python3 -I d_grid_enumeration.py
javac RoundingDemo.java && java RoundingDemo

check_numbers.py alone does not catch a mutation of the data behind the registered figures (for example the sums of Tables 1 and 2): that is caught by rerunning d_build_numbers.py, which asserts them, and by SHA256SUMS-D and the manifest. SHA256SUMS-D lists the hashes of the files committed in rounding_subreceptor/ (check with sha256sum -c SHA256SUMS-D in that directory). The regenerated large outputs (about 220 MB each, not committed) are listed in SHA256SUMS-D-regenerated, which is a record and cannot be checked in place. numbers.json registers the figures of this paper against d_numbers.json.

References

Code (all at the pins of Paper 1)

• [C1] AERIUS monorepo, gitlab.com/AERIUS/AERIUS @ fb8812fd74f9b5a548c8772fde3a7e21be238e54 (Release 2026.0.0); 71a601536cab (Release 2025.3); 7bf619a (Release 2025).
• [C2] OPS 5.3.1.0-5 production package, nexus.aerius.nl/repository/resources (anonymous).

Repository

• [P1] P1-aerius-deconstructed.md, Finding F9 and Section 5.8; P1-supplement.md, E54 and E56.
• [N1] f9-and-version-drift.md; interim-review-r1-v3.md; legal-and-documentation.md, row d.
• [T1] f9/ (harness, Monte Carlo, main.csv, sum_near_summary.txt, sum_near_flips.csv).

Documents

• [R2] Afdeling advisering van de Raad van State, voorlichting W11.25.00063/IV, 23 May 2025, pdf pp. 7, 8, 10. https://zoek.officielebekendmakingen.nl/blg-1198305.pdf (read 2026-10-26).

• [R3] Backes, C. and Petersen, A., Reflectie over de (juridische) noodzaak van flankerend beleid bij de introductie van een RKO, 21 January 2026, p. 1 (Tweede Kamer 2026D14542). https://www.tweedekamer.nl/downloads/document?id=2026D14542 (read 2026-10-26).

• [R4] Rijksoverheid.nl, Rekenkundige ondergrens voor juridische toets naar Raad van State, 18 March 2025. https://www.rijksoverheid.nl/actueel/nieuws/2025/03/18/rekenkundige-ondergrens-voor-juridische-toets-naar-raad-van-state (read 2026-10-26).

• [R1] AERIUS, Handboek werken met AERIUS Calculator 2026 (version 2026-09-29), Section 3.6 (sub-receptors), pdf p. 33. https://www.aeriusproducten.nl/site/binaries/site-content/collections/documents/2026/10/06/handboek-werken-met-aerius-calculator-2026/handboek-werken-met-calculator-2026.pdf (read 2026-10-25; quotations of pdf pp. 23 and 112 checked 2026-10-26).