Validation as OpenModelica's library testing: omc's tube comparison in Julia

README.md

@@ -18,6 +18,12 @@ Each model is tested through up to four pipeline phases:
 
 Results are compared against Dymola-generated reference trajectories from the MAP-LIB reference results repository.
 
+Each listed signal is compared as OpenModelica's library testing compares it (omc's `diffSimulationResults` with
+relTol 3e-3, relTolDiffMinMax 3e-3, rangeDelta 1e-3; ported to Julia as `tube_compare`, `src/tubes.jl`): OM.jl's
+values at every time of the reference must lie in a band around it, 0.1 % of the time span wide in time and 0.3 %
+of the value and of the signal's range high; a reference event is not compared. A value that is not finite fails.
+The report's signal viewer plots the band at 101 evenly spaced times.
+
 ## Latest coverage
 
 <!-- coverage:begin (scripts/publish_coverage.jl) -->

models/models.toml

@@ -17,8 +17,11 @@ reference_source = "MAP-LIB_ReferenceResults (Dymola, branch v3.2.3)"
 #                                                  (default = solver default)
 #   solverReltol   = <float>                       (default: the experiment's Tolerance)
 #   solverAtol     = <float>                       (0.0 = SciML default ~1e-6)
-#                  Note: distinct from atol/reltol above which are validation
-#                  tolerances for trajectory comparison. Use solverReltol /
+#                  Note: distinct from atol/reltol, the tolerances of an inline
+#                  `reference` (final values). A reference trajectory (CSV) is
+#                  compared as OpenModelica's library testing does (tube_compare:
+#                  relTol 3e-3, relTolDiffMinMax 3e-3, rangeDelta 1e-3, 2026-10-28);
+#                  atol/reltol do not apply to it. Use solverReltol /
 #                  solverAtol to loosen the integrator (e.g. near-singular
 #                  Jacobians that force dt below eps under default tolerances).
 #   observedFilter = ["regex", ...]                 ([] = OM.simulate default)

src/OMLibraryTesting.jl

@@ -41,6 +41,7 @@ import OMFrontend
 
 include("types.jl")
 include("registry.jl")
+include("tubes.jl")
 include("comparison.jl")
 include("discovery.jl")
 include("runner.jl")

src/discovery.jl

@@ -33,7 +33,7 @@
     discover_experiments(; library, version, omc_path, filter) -> Vector{ModelSpec}
 
 Use omc to enumerate all experiment models in a Modelica library.
-Returns a ModelSpec for each model with stopTime from the experiment annotation.
+Returns a ModelSpec for each model with stopTime, Tolerance and startTime from the experiment annotation.
 Without omc on the path (a CI runner), the list recorded from omc in `models/`
 (`record_experiments`) is used instead.
 """
@@ -85,12 +85,14 @@ function discover_experiments(; library::String = "Modelica",
         #= The experiment's Tolerance is the solver's relative tolerance, as
            in the tools that produced the references (default 1e-6). =#
         tolerance = length(parts) >= 3 ? something(tryparse(Float64, parts[3]), 0.0) : 0.0
+        #= The experiment's StartTime (lists recorded before 2026-10-28 have none: 0). =#
+        startTime = length(parts) >= 4 ? something(tryparse(Float64, parts[4]), 0.0) : 0.0
         domain = _extract_domain(name)
         key = _name_to_key(name)
         push!(specs, ModelSpec(name, key, domain, stopTime,
                                UNKNOWN, Dict{String, Float64}(),
                                0.01, 3e-3, "", Dict{String, String}(), "",
-                               Set{Phase}(), "", 0.0, "", 0.0, tolerance))
+                               Set{Phase}(), "", 0.0, "", 0.0, tolerance, String[], 0.0, startTime))
     end
     sort!(specs, by = s -> (s.domain, s.name))
     if cache
@@ -103,15 +105,15 @@ function discover_experiments(; library::String = "Modelica",
     return specs
 end
 
-#= The experiment models of a library as omc lists them: one `name|stopTime|tolerance` line each. =#
+#= The experiment models of a library as omc lists them: one `name|stopTime|tolerance|startTime` line each. =#
 function _omc_experiments(library::String, version::String, omc_path::String)::String
     script = """
     loadModel($library, {"$version"});
     names := getClassNames($library, recursive=true, qualified=true);
     for n in names loop
       if isExperiment(n) then
         (startTime, stopTime, tolerance, numberOfIntervals, interval) := getSimulationOptions(n);
-        print(typeNameString(n) + "|" + String(stopTime) + "|" + String(tolerance) + "\\n");
+        print(typeNameString(n) + "|" + String(stopTime) + "|" + String(tolerance) + "|" + String(startTime) + "\\n");
       end if;
     end for;
     """
@@ -265,7 +267,7 @@ function merge_overrides!(specs::Vector{ModelSpec},
                                 atol, reltol, referenceFile,
                                 sig_map, issue, skip, solverName, dtmaxVal, initAlgName,
                                 solverAtolVal, solverReltolVal, observedFilterVal,
-                                maxitersVal))
+                                maxitersVal, spec.startTime))
     end
     n_broken = count(s -> s.expected == BROKEN, merged)
     n_skipped = count(s -> !isempty(s.skipPhases), merged)
@@ -316,7 +318,7 @@ function auto_detect_references!(specs::Vector{ModelSpec},
                                      spec.signalMapping, spec.issue, spec.skipPhases,
                                      spec.solver, spec.dtmax, spec.initAlg,
                                      spec.solverAtol, spec.solverReltol, spec.observedFilter,
-                                     spec.maxiters))
+                                     spec.maxiters, spec.startTime))
         else
             push!(result, spec)
         end

src/report.jl

@@ -321,6 +321,7 @@ const _SIGNALS_CSS = raw"""
       .sig-plot .ref { fill: none; stroke: #555; stroke-width: 1.5; }
       .sig-plot .omjl { fill: none; stroke: #0366d6; stroke-width: 1.5; stroke-dasharray: 5 3; }
       .sig-plot .worst { fill: none; stroke: #b31d28; stroke-width: 2; }
+      .sig-plot .worst-line { stroke: #b31d28; stroke-width: 1; stroke-dasharray: 3 2; }
       .sig-plot .cursor { stroke: #999; stroke-dasharray: 2 2; }
       .sig-plot text { font-size: 10px; fill: #666; font-family: monospace; }
       #error-modal-text.values { color: #333; }
@@ -395,9 +396,14 @@ function renderSignals(td, m) {
   var failed = sigs.filter(function(s) { return !s.passed; }).length;
   var head = document.createElement('div');
   head.className = 'sig-head';
+  var banded = sigs.some(function(s) { return s.low; });   // OpenModelica's comparison (older files: atol/reltol)
   head.innerHTML = '<span class="sig-key-ref">— reference</span> &nbsp; <span class="sig-key-omjl">- - OM.jl</span> &nbsp; ' +
-    'grey: the tolerance around the reference (atol ' + fmt(m.atol) + ', reltol ' + fmt(m.reltol) + '); values with 6 significant digits. ';
-  head.appendChild(document.createTextNode(sigs.length + ' signals compared, ' + failed + ' outside the tolerance' +
+    (banded ? 'grey: the band OM.jl must stay in, as OpenModelica\'s library testing compares (0.1 % of the time span wide ' +
+              'in time, 0.3 % of the value and of the signal\'s range high); checked at every time of the reference except ' +
+              'its events; red: where OM.jl is furthest outside it' :
+              'grey: the tolerance around the reference (atol ' + fmt(m.atol) + ', reltol ' + fmt(m.reltol) + ')') +
+    '; values with 6 significant digits. ';
+  head.appendChild(document.createTextNode(sigs.length + ' signals compared, ' + failed + ' failed' +
     (m.validated ? '.' : '; the model did not validate.') + ' Point at a plot to read its values; click a name for all of them.'));
   td.appendChild(head);
   var grid = document.createElement('div');
@@ -413,15 +419,27 @@ function signalCard(s, m) {
   name.href = '#'; name.textContent = s.name; name.title = 'All compared values';
   name.onclick = function() { showValues(s, m); return false; };
   mark.className = s.passed ? 'sig-mark-ok' : 'sig-mark-fail';
-  mark.textContent = s.passed ? ' ✓' : ' ✗ outside the tolerance';
+  mark.textContent = s.passed ? ' ✓' : (s.low ? ' ✗ outside the band' : ' ✗ outside the tolerance');
   title.appendChild(name); title.appendChild(mark);
   stats.className = 'sig-stats';
-  stats.textContent = s.max_abs_err === null ? 'not compared: no OM.jl value' :
-    'max |error| ' + fmt(s.max_abs_err) + ' (relative ' + fmt(s.max_rel_err) + ') at t = ' + fmt(s.worst_time);
+  stats.textContent = signalStats(s);
   readout.className = 'sig-readout';
   card.appendChild(title); card.appendChild(stats); card.appendChild(signalPlot(s, m, readout)); card.appendChild(readout);
   return card;
 }
+function signalStats(s) {
+  if (s.max_abs_err === null) return 'not compared: no OM.jl value';
+  var text = 'max |error| ' + fmt(s.max_abs_err) + ' (relative ' + fmt(s.max_rel_err) + ')';
+  if (s.compared === undefined) {   // older files: where the largest error is (0 when none)
+    if (s.worst_time !== null && s.worst_time !== undefined) text += ' at t = ' + fmt(s.worst_time);
+    return text;
+  }
+  return text + (s.outside ?
+    '; outside the band at ' + s.outside + ' of ' + s.compared + ' times, furthest at t = ' + fmt(s.worst_time) +
+      ' (OM.jl ' + fmt(s.worst_omjl) + ', reference ' + fmt(s.worst_reference) + ')' :
+    '; in the band at all ' + (s.compared - (s.events || 0)) + ' compared times' +
+      (s.events ? ' (' + s.events + ' at events, not compared)' : ''));
+}
 function svgEl(name, attrs) {
   var e = document.createElementNS(SVGNS, name);
   for (var k in attrs) e.setAttribute(k, attrs[k]);
@@ -441,11 +459,15 @@ function signalPlot(s, m, readout) {
   var W = 360, H = 160, L = 62, R = 8, T = 8, B = 20;
   var t = sampleTimes(s, m), ref = s.reference, om = s.omjl, n = t.length, i;
   var tol = ref.map(function(v) { return v === null ? null : m.atol + m.reltol * Math.abs(v); });
+  var high = s.high || ref.map(function(v, k) { return v === null ? null : v + tol[k]; });   // older files: atol/reltol
+  var low = s.low || ref.map(function(v, k) { return v === null ? null : v - tol[k]; });
   var ys = [];
   for (i = 0; i < n; i++) {
-    if (ref[i] !== null) ys.push(ref[i] + tol[i], ref[i] - tol[i]);
+    if (high[i] !== null && low[i] !== null) ys.push(high[i], low[i]);
     if (om[i] !== null) ys.push(om[i]);
   }
+  var worstOmjl = !s.passed && s.worst_omjl !== undefined && s.worst_omjl !== null;   // where it is, between samples too
+  if (worstOmjl) ys.push(s.worst_omjl);
   var y0 = ys.length ? Math.min.apply(null, ys) : 0, y1 = ys.length ? Math.max.apply(null, ys) : 1;
   if (y1 - y0 <= 1e-12 * Math.max(1, Math.abs(y1))) { var c = y0, d = 1e-6 * Math.max(1, Math.abs(c)); y0 = c - d; y1 = c + d; }
   var pad = 0.05 * (y1 - y0); y0 -= pad; y1 += pad;
@@ -455,11 +477,14 @@ function signalPlot(s, m, readout) {
   var svg = svgEl('svg', {viewBox: '0 0 ' + W + ' ' + H, width: W, height: H, 'class': 'sig-plot'});
   svg.appendChild(svgEl('rect', {x: L, y: T, width: W - L - R, height: H - T - B, 'class': 'frame'}));
   var up = [], lo = [];
-  for (i = 0; i < n; i++) if (ref[i] !== null) { up.push(X(t[i]) + ',' + Y(ref[i] + tol[i])); lo.unshift(X(t[i]) + ',' + Y(ref[i] - tol[i])); }
+  for (i = 0; i < n; i++) if (high[i] !== null && low[i] !== null) { up.push(X(t[i]) + ',' + Y(high[i])); lo.unshift(X(t[i]) + ',' + Y(low[i])); }
   if (up.length) svg.appendChild(svgEl('polygon', {points: up.concat(lo).join(' '), 'class': 'band'}));
   svgLine(svg, t, ref, X, Y, 'ref');
   svgLine(svg, t, om, X, Y, 'omjl');
-  if (!s.passed && s.worst_time !== null) {
+  if (worstOmjl && s.worst_time !== null) {
+    svg.appendChild(svgEl('line', {x1: X(s.worst_time), x2: X(s.worst_time), y1: T, y2: H - B, 'class': 'worst-line'}));
+    svg.appendChild(svgEl('circle', {cx: X(s.worst_time), cy: Y(s.worst_omjl), r: 4, 'class': 'worst'}));
+  } else if (!s.passed && s.worst_time !== null && s.worst_time !== undefined) {   // older files: the nearest sample
     var w = Math.round((s.worst_time - t0) / (t1 - t0) * (n - 1));
     w = Math.max(0, Math.min(n - 1, w));
     if (om[w] !== null) svg.appendChild(svgEl('circle', {cx: X(t[w]), cy: Y(om[w]), r: 4, 'class': 'worst'}));
@@ -485,10 +510,11 @@ function signalPlot(s, m, readout) {
 function showValues(s, m) {
   var t = sampleTimes(s, m), lines = [];
   function pad(text) { return (text + '                ').slice(0, 16); }
-  lines.push(pad('t') + pad('reference') + pad('OM.jl') + '|diff|');
+  lines.push(pad('t') + pad('reference') + pad('OM.jl') + pad('|diff|') + (s.low ? pad('band low') + 'band high' : ''));
   for (var i = 0; i < t.length; i++) {
     var r = s.reference[i], o = s.omjl[i];
-    lines.push(pad(fmt(t[i])) + pad(fmt(r)) + pad(fmt(o)) + (r !== null && o !== null ? fmt(Math.abs(o - r)) : '—'));
+    lines.push(pad(fmt(t[i])) + pad(fmt(r)) + pad(fmt(o)) + (r !== null && o !== null ? pad(fmt(Math.abs(o - r))) : pad('—')) +
+               (s.low ? pad(fmt(s.low[i])) + fmt(s.high[i]) : ''));
   }
   document.getElementById('error-modal-title').textContent = s.name + ': the compared values';
   var pre = document.getElementById('error-modal-text');

src/runner.jl

@@ -391,7 +391,8 @@ function _run_model_phases(model_name::String,
                            solver_atol::Float64 = 0.0,
                            solver_reltol::Float64 = 0.0,
                            observed_filter::Vector{String} = String[],
-                           maxiters::Float64 = 0.0)
+                           maxiters::Float64 = 0.0,
+                           start_time::Float64 = 0.0)
     results = Tuple{Int, Bool, Float64, Union{Nothing, String}}[]
     highest = 0
     sol = nothing
@@ -435,7 +436,7 @@ function _run_model_phases(model_name::String,
                 if solver_atol > 0.0
                     _tol = (; _tol..., abstol = solver_atol)
                 end
-                sol = OM.simulate(model_name; msl..., stopTime = stop_time, _sa..., _extra..., _ia..., _of..., _mi..., _tol...)
+                sol = OM.simulate(model_name; msl..., startTime = start_time, stopTime = stop_time, _sa..., _extra..., _ia..., _of..., _mi..., _tol...)
                 #= Terminated is a successful end per Modelica semantics: a model
                    calling terminate() completed its intended scenario. =#
                 if sol.retcode != ReturnCode.Success && sol.retcode != ReturnCode.Terminated
@@ -447,7 +448,7 @@ function _run_model_phases(model_name::String,
                     temp_spec = ModelSpec(
                         model_name, "", "", stop_time, UNKNOWN,
                         reference, atol, reltol, reference_file, signal_mapping, "",
-                        Set{Phase}())
+                        Set{Phase}(), "", 0.0, "", 0.0, 0.0, String[], 0.0, start_time)
                     (passed, comparisons) = validate_against_reference(sol, temp_spec, ref_dir;
                         terminatedAt = sol.retcode == ReturnCode.Terminated ? Float64(sol.t[end]) : nothing)
                     signals = [_signal_tuple(c) for c in comparisons]
@@ -505,11 +506,12 @@ end
 
 #= A SignalComparison's fields, for a worker's result (plain types across processes). =#
 const SIGNAL_TUPLE = Tuple{String, Bool, Float64, Float64, Float64, Float64, Float64,
-                           Vector{Float64}, Vector{Float64}, Vector{Float64}}
+                           Vector{Float64}, Vector{Float64}, Vector{Float64},
+                           Vector{Float64}, Vector{Float64}, Int, Int, Int}
 
 _signal_tuple(c::SignalComparison)::SIGNAL_TUPLE =
     (c.name, c.passed, c.max_abs_err, c.max_rel_err, c.worst_time, c.worst_actual, c.worst_expected,
-     c.times, c.expected, c.actual)
+     c.times, c.expected, c.actual, c.low, c.high, c.n_compared, c.n_events, c.n_outside)
 
 # ── Worker Execution with Timeout ───────────────────────────────────────
 
@@ -593,7 +595,7 @@ function run_on_worker(mgr::WorkerManager, spec::ModelSpec,
         spec.referenceFile, spec.reference, spec.atol, spec.reltol,
         spec.signalMapping, mgr.ref_dir, phase_ints, check_sim_code,
         spec.solver, spec.dtmax, spec.initAlg,
-        spec.solverAtol, spec.solverReltol, spec.observedFilter, spec.maxiters)
+        spec.solverAtol, spec.solverReltol, spec.observedFilter, spec.maxiters, spec.startTime)
 
     #= Wait for the result in a local task and poll that task. isready(future) asks
        the worker, and a worker running code that never yields (a pegged solve or

src/signals.jl

@@ -16,13 +16,19 @@ The signals the validations of `results` compared, as JSON, by model name:
 ```
 {"models": {"<model>": {"validated": true, "atol": 1e-4, "reltol": 3e-3, "time": [t0, t1, n],
                         "signals": [{"name": "J1.w", "passed": true, "max_abs_err": ..., "max_rel_err": ...,
-                                     "worst_time": ..., "reference": [...], "omjl": [...]}]}}}
+                                     "worst_time": ..., "worst_omjl": ..., "worst_reference": ...,
+                                     "compared": 5051, "events": 20, "outside": 0,
+                                     "reference": [...], "omjl": [...], "low": [...], "high": [...]}]}}}
 ```
 
 `time` is the comparison's `n` evenly spaced times from `t0` to `t1` (a signal compared at other times
 has its own `time`); `reference` and `omjl` are the values there, 6 significant digits, `null` where
 there is none (OM.jl's after the comparison stopped: the signal is not in the solution, or a value is
-not finite). Models without compared signals are left out.
+not finite); `low` and `high` the band OM.jl's values must lie in there (left out when the comparison
+stopped before it). The signal was compared as OpenModelica's library testing does (`compare_signal`) at
+`compared` times of the reference, `events` of them at its events (not compared), `outside` outside the
+band: `worst_time` is where OM.jl is furthest outside it (`null` when inside everywhere), `worst_omjl` and
+`worst_reference` the values there (left out without one). Models without compared signals are left out.
 """
 function signals_json(results::Vector{ModelResult})::String
     local io = IOBuffer()
@@ -70,6 +76,13 @@ function _model_signals_json(io::IO, r::ModelResult)::Nothing
         _json_number(io, s.max_rel_err)
         print(io, ",\"worst_time\":")
         _json_number(io, s.worst_time)
+        if isfinite(s.worst_time)
+            print(io, ",\"worst_omjl\":")
+            _json_number(io, s.worst_actual)
+            print(io, ",\"worst_reference\":")
+            _json_number(io, s.worst_expected)
+        end
+        print(io, ",\"compared\":", s.n_compared, ",\"events\":", s.n_events, ",\"outside\":", s.n_outside)
         if !isempty(s.times) && s.times != times
             print(io, ",\"time\":")
             _json_times(io, s.times)
@@ -78,6 +91,12 @@ function _model_signals_json(io::IO, r::ModelResult)::Nothing
         _json_numbers(io, s.expected)
         print(io, ",\"omjl\":")
         _json_numbers(io, s.actual)
+        if !isempty(s.low)
+            print(io, ",\"low\":")
+            _json_numbers(io, s.low)
+            print(io, ",\"high\":")
+            _json_numbers(io, s.high)
+        end
         print(io, "}")
     end
     print(io, "]}")