Skip to content

Claude/optimize sensor resolve 1 xpjc #38

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
HanSur94 merged 12 commits into from
Mar 19, 2026
Merged

Claude/optimize sensor resolve 1 xpjc #38

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
12 commits
Select commit Hold shift + click to select a range
26b7e56
Optimize toStepFunction and resolve segment boundary collection
claude Mar 19, 2026
568f260
Add to_step_function_mex for C-accelerated step-function conversion
claude Mar 19, 2026
9e24ad2
Add SIMD acceleration to to_step_function_mex
claude Mar 19, 2026
1620016
Add tests for toStepFunction and to_step_function_mex
claude Mar 19, 2026
6642fb9
fix: extract toStepFunction to standalone private file
HanSur94 Mar 19, 2026
50ac054
fix: ensure install() triggers build for new to_step_function_mex
HanSur94 Mar 19, 2026
76c6727
bench: add 500M-point resolve stress test with frequent state changes
HanSur94 Mar 19, 2026
ecaba36
fix: probe MEX files by path instead of exist() in stress benchmark
HanSur94 Mar 19, 2026
ac4dbc2
bench: add FastSense plot with thresholds to stress benchmark
HanSur94 Mar 19, 2026
db5c3cc
bench: add JIT warmup to eliminate first-run overhead
HanSur94 Mar 19, 2026
763306b
feat: add JIT warmup to install() for instant first-call performance
HanSur94 Mar 19, 2026
d3356a5
fix: use hidden figure for JIT warmup render
HanSur94 Mar 19, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
284 changes: 284 additions & 0 deletions benchmarks/benchmark_resolve_stress.m
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,284 @@
function benchmark_resolve_stress()
%BENCHMARK_RESOLVE_STRESS Stress test for Sensor.resolve() at 500M points.
% Creates a realistic sensor with 500M datapoints and 4 thresholds,
% each governed by different condition rules across 2 state channels
% that change frequently (~5000 transitions each).
%
% This exercises the full resolve pipeline:
% - Segment boundary collection with ~10K unique boundaries
% - Composite state evaluation at each boundary
% - Condition matching & segment index mapping
% - Batch violation detection (MEX or vectorized)
% - Threshold merge & step-function conversion
%
% Run:
% >> benchmark_resolve_stress

repo_root = fileparts(fileparts(mfilename('fullpath')));
addpath(repo_root);
install();

N = 500e6;

fprintf('\n');
fprintf('============================================================\n');
fprintf(' Sensor.resolve() Stress Test — %.0fM datapoints\n', N / 1e6);
fprintf('============================================================\n');
fprintf(' Building sensor data... ');

tBuild = tic;

% ------------------------------------------------------------------
% 1. Sensor data: 500M points of realistic process signal
% Base signal with multiple frequency components + noise + drift
% ------------------------------------------------------------------
x = linspace(0, 10000, N);

% Layered signal: slow drift + medium oscillation + fast noise
y = 50 ...
+ 10 * sin(2*pi*x / 3000) ... % slow thermal drift
+ 8 * sin(2*pi*x / 200) ... % process oscillation
+ 5 * randn(1, N); % measurement noise

% Inject deliberate excursions to guarantee violations in every regime
% High spikes in running states
y(round(N*0.05):round(N*0.05)+2000) = 95 + 5*randn(1, 2001);
y(round(N*0.35):round(N*0.35)+3000) = 90 + 5*randn(1, 3001);
y(round(N*0.65):round(N*0.65)+1500) = 85 + 3*randn(1, 1501);
% Low dips in evacuated states
y(round(N*0.15):round(N*0.15)+2500) = 10 + 3*randn(1, 2501);
y(round(N*0.55):round(N*0.55)+1000) = 5 + 2*randn(1, 1001);
y(round(N*0.85):round(N*0.85)+2000) = 15 + 4*randn(1, 2001);

fprintf('%.1f s\n', toc(tBuild));

% ------------------------------------------------------------------
% 2. State channels with many transitions (~5000 each)
% This creates ~10K segment boundaries — the core stress factor
% for the resolve pipeline beyond raw data size.
% ------------------------------------------------------------------
fprintf(' Building state channels... ');
tState = tic;

% Machine state: cycles through 0 (idle) → 1 (running) → 2 (evacuated)
% with irregular interval lengths to create realistic patterns
nMachineTransitions = 5000;
rng(42);
machineIntervals = 1.5 + 1.5 * rand(1, nMachineTransitions);
machineX = [0, cumsum(machineIntervals)];
% Scale to fit data range
machineX = machineX / machineX(end) * x(end);
% Cycle through states 0→1→2→1→0→1→2→...
statePattern = [0 1 2 1];
machineY = zeros(1, numel(machineX));
for i = 1:numel(machineX)
machineY(i) = statePattern(mod(i-1, numel(statePattern)) + 1);
end

scMachine = StateChannel('machine');
scMachine.X = machineX;
scMachine.Y = machineY;

% Operating zone: string-valued, toggles between A/B/C
% Different transition frequency than machine state → many unique combos
nZoneTransitions = 4000;
zoneIntervals = 2.0 + 2.0 * rand(1, nZoneTransitions);
zoneX = [0, cumsum(zoneIntervals)];
zoneX = zoneX / zoneX(end) * x(end);
zonePattern = {'A', 'B', 'C'};
zoneY = cell(1, numel(zoneX));
for i = 1:numel(zoneX)
zoneY{i} = zonePattern{mod(i-1, numel(zonePattern)) + 1};
end

scZone = StateChannel('zone');
scZone.X = zoneX;
scZone.Y = zoneY;

fprintf('%.1f s\n', toc(tState));
fprintf(' Machine transitions: %d\n', numel(machineX));
fprintf(' Zone transitions: %d\n', numel(zoneX));

% ------------------------------------------------------------------
% 3. Build sensor with 4 threshold rules (different conditions each)
% ------------------------------------------------------------------
fprintf(' Configuring sensor + thresholds...\n');

s = Sensor('stress_test', 'Name', 'Process Stress Test', 'Units', 'bar');
s.X = x;
s.Y = y;
s.addStateChannel(scMachine);
s.addStateChannel(scZone);

% Rule 1: Upper alarm when running (machine==1), any zone
% → Active in ~25% of timeline, single-condition
s.addThresholdRule(struct('machine', 1), 75, ...
'Direction', 'upper', 'Label', 'HH (running)', ...
'Color', [0.9 0.1 0.1], 'LineStyle', '--');

% Rule 2: Lower alarm when evacuated (machine==2), any zone
% → Active in ~25% of timeline, single-condition
s.addThresholdRule(struct('machine', 2), 25, ...
'Direction', 'lower', 'Label', 'LL (evacuated)', ...
'Color', [0.1 0.1 0.9], 'LineStyle', '--');

% Rule 3: Upper alarm when running AND zone B (two conditions)
% → Active in ~8% of timeline (intersection), stricter limit
s.addThresholdRule(struct('machine', 1, 'zone', 'B'), 65, ...
'Direction', 'upper', 'Label', 'HH (running+B)', ...
'Color', [0.8 0.0 0.8], 'LineStyle', ':');

% Rule 4: Lower alarm when idle AND zone C (two conditions)
% → Active in ~8% of timeline (intersection)
s.addThresholdRule(struct('machine', 0, 'zone', 'C'), 30, ...
'Direction', 'lower', 'Label', 'LL (idle+C)', ...
'Color', [0.0 0.5 0.8], 'LineStyle', ':');

fprintf(' Rules: %d (%d upper, %d lower)\n', ...
numel(s.ThresholdRules), 2, 2);

% ------------------------------------------------------------------
% 4. MEX availability check
% MEX files live in private/ dirs, so exist() can't see them from
% here. Probe the actual file paths instead.
% ------------------------------------------------------------------
ext = mexext();
fsPriv = fullfile(repo_root, 'libs', 'FastSense', 'private');
stPriv = fullfile(repo_root, 'libs', 'SensorThreshold', 'private');
mexProbes = {
'binary_search_mex', fsPriv
'compute_violations_mex', stPriv
'to_step_function_mex', stPriv
};
fprintf('\n MEX status:\n');
for i = 1:size(mexProbes, 1)
mexFile = fullfile(mexProbes{i,2}, [mexProbes{i,1} '.' ext]);
if exist(mexFile, 'file')
fprintf(' %-30s compiled\n', mexProbes{i,1});
else
fprintf(' %-30s NOT compiled (MATLAB fallback)\n', mexProbes{i,1});
end
end

% ------------------------------------------------------------------
% 5. JIT warmup — resolve a tiny sensor to force-compile all code paths
% ------------------------------------------------------------------
fprintf('\n JIT warmup... ');
tWarm = tic;
sw = Sensor('warmup');
sw.X = [0 1 2 3]; sw.Y = [50 60 40 70];
scw1 = StateChannel('machine'); scw1.X = [0 1 2]; scw1.Y = [0 1 2];
scw2 = StateChannel('zone'); scw2.X = [0 1]; scw2.Y = {'A','B'};
sw.addStateChannel(scw1); sw.addStateChannel(scw2);
sw.addThresholdRule(struct('machine', 1), 55, 'Direction', 'upper');
sw.addThresholdRule(struct('machine', 2), 35, 'Direction', 'lower');
sw.addThresholdRule(struct('machine', 1, 'zone', 'B'), 50, 'Direction', 'upper');
sw.addThresholdRule(struct('machine', 0, 'zone', 'A'), 40, 'Direction', 'lower');
sw.resolve();
clear sw scw1 scw2;
fprintf('%.3f s\n', toc(tWarm));

% ------------------------------------------------------------------
% 6. Run resolve() with timing
% ------------------------------------------------------------------
nRuns = 3;
fprintf('\n Resolving (%d runs)...\n', nRuns);

times = zeros(1, nRuns);
for r = 1:nRuns
% Rebuild sensor each run to avoid cache effects
sr = Sensor('stress_test');
sr.X = x; sr.Y = y;
sr.addStateChannel(scMachine);
sr.addStateChannel(scZone);
sr.addThresholdRule(struct('machine', 1), 75, ...
'Direction', 'upper', 'Label', 'HH (running)');
sr.addThresholdRule(struct('machine', 2), 25, ...
'Direction', 'lower', 'Label', 'LL (evacuated)');
sr.addThresholdRule(struct('machine', 1, 'zone', 'B'), 65, ...
'Direction', 'upper', 'Label', 'HH (running+B)');
sr.addThresholdRule(struct('machine', 0, 'zone', 'C'), 30, ...
'Direction', 'lower', 'Label', 'LL (idle+C)');

tic;
sr.resolve();
times(r) = toc;

fprintf(' Run %d: %.3f s\n', r, times(r));
end

% Use the last run's sensor for result inspection
s = sr;

% ------------------------------------------------------------------
% 7. Results summary
% ------------------------------------------------------------------
fprintf('\n============================================================\n');
fprintf(' RESULTS\n');
fprintf('============================================================\n');
fprintf(' Data points: %s\n', format_size(N));
fprintf(' State transitions: %d (machine) + %d (zone)\n', ...
numel(machineX), numel(zoneX));
fprintf(' Threshold rules: %d\n', numel(s.ThresholdRules));
fprintf(' Resolved lines: %d\n', numel(s.ResolvedThresholds));

totalViol = 0;
for v = 1:numel(s.ResolvedViolations)
nv = numel(s.ResolvedViolations(v).X);
totalViol = totalViol + nv;
fprintf(' [%s] %s: %s violations\n', ...
s.ResolvedViolations(v).Direction, ...
s.ResolvedViolations(v).Label, ...
format_size(nv));
end
fprintf(' Total violations: %s\n', format_size(totalViol));
fprintf('\n');
fprintf(' resolve() time:\n');
fprintf(' median: %.3f s\n', median(times));
fprintf(' min: %.3f s\n', min(times));
fprintf(' max: %.3f s\n', max(times));
fprintf(' throughput: %.1f M pts/s\n', (N / 1e6) / median(times));

% Memory estimate (approximate)
memBytes = N * 8 * 2; % X + Y arrays
for v = 1:numel(s.ResolvedViolations)
memBytes = memBytes + numel(s.ResolvedViolations(v).X) * 8 * 2;
end
for t = 1:numel(s.ResolvedThresholds)
memBytes = memBytes + numel(s.ResolvedThresholds(t).X) * 8 * 2;
end
fprintf(' approx memory: %.1f GB (sensor) + %.1f MB (results)\n', ...
N * 8 * 2 / 1e9, (memBytes - N*8*2) / 1e6);

% ------------------------------------------------------------------
% 8. Plot with FastSense
% ------------------------------------------------------------------
fprintf('\n Rendering with FastSense...\n');
tPlot = tic;

fp = FastSense();
fp.addSensor(s, 'ShowThresholds', true);
fp.render();
title('Sensor.resolve() Stress Test — 500M pts, 4 thresholds, ~9K state transitions');
xlabel('Time [s]');
ylabel('Process Value [bar]');

fprintf(' Render time: %.3f s\n', toc(tPlot));

fprintf('\n============================================================\n');
fprintf(' Done.\n');
fprintf('============================================================\n');
end


function s = format_size(N)
if N >= 1e9
s = sprintf('%.1fB', N / 1e9);
elseif N >= 1e6
s = sprintf('%.1fM', N / 1e6);
elseif N >= 1e3
s = sprintf('%.1fK', N / 1e3);
else
s = sprintf('%d', N);
end
end
65 changes: 62 additions & 3 deletions install.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,8 @@ function install()
% 2. Adds examples, benchmarks, and tests to the path
% 3. Compiles MEX accelerators if not yet built (first run only)
% 4. Verifies the installation (first run only)
% 5. JIT warmup — runs a tiny end-to-end workflow to force-compile
% all hot code paths (once per MATLAB session)
%
% Directories added:
% libs/FastSense — core plotting engine
Expand Down Expand Up @@ -60,6 +62,9 @@ function install()
if needs_build(root)
first_run(root);
end

% --- Once per session: JIT warmup ---
jit_warmup();
end

function yes = needs_build(root)
Expand All @@ -69,9 +74,19 @@ function install()
return;
end
mex_dir = fullfile(root, 'libs', 'FastSense', 'private');
probe = fullfile(mex_dir, ['binary_search_mex.' mexext()]);
probe_oct = fullfile(mex_dir, 'binary_search_mex.mex');
yes = exist(probe, 'file') ~= 3 && exist(probe_oct, 'file') ~= 3;
sensor_dir = fullfile(root, 'libs', 'SensorThreshold', 'private');
% Probe a representative MEX from each library — if any are missing,
% trigger build_mex() which will compile only the missing ones.
probes = {
fullfile(mex_dir, ['binary_search_mex.' mexext()])
fullfile(mex_dir, 'binary_search_mex.mex')
fullfile(sensor_dir, ['to_step_function_mex.' mexext()])
fullfile(sensor_dir, 'to_step_function_mex.mex')
};
% Need build if either probe set is missing
core_ok = exist(probes{1}, 'file') == 3 || exist(probes{2}, 'file') == 3;
step_ok = exist(probes{3}, 'file') == 3 || exist(probes{4}, 'file') == 3;
yes = ~core_ok || ~step_ok;
end

function first_run(root)
Expand Down Expand Up @@ -159,3 +174,47 @@ function verify_installation(root)
fprintf(' %d/%d checks passed, %d warnings\n', n_ok, total, n_warn);
end
end


function jit_warmup()
%JIT_WARMUP Force MATLAB's JIT to compile all hot code paths once per session.
% Runs a tiny end-to-end workflow (sensor creation, state channels,
% threshold resolve, downsampling, rendering setup) on trivial data.
% Subsequent calls are no-ops. Adds ~0.1-0.3 s on first install() call
% but eliminates multi-second JIT overhead on real workloads.
persistent warmedUp;
if ~isempty(warmedUp)
return;
end
warmedUp = true;

try
% --- Sensor + StateChannel + resolve pipeline ---
sw = Sensor('__jit_warmup__');
sw.X = [0 1 2 3 4 5];
sw.Y = [50 60 40 70 30 80];

sc1 = StateChannel('s1');
sc1.X = [0 2 4]; sc1.Y = [0 1 2];
sc2 = StateChannel('s2');
sc2.X = [0 3]; sc2.Y = {'A', 'B'};
sw.addStateChannel(sc1);
sw.addStateChannel(sc2);

sw.addThresholdRule(struct('s1', 1), 65, 'Direction', 'upper');
sw.addThresholdRule(struct('s1', 2), 35, 'Direction', 'lower');
sw.addThresholdRule(struct('s1', 1, 's2', 'B'), 55, 'Direction', 'upper');
sw.addThresholdRule(struct('s1', 0, 's2', 'A'), 45, 'Direction', 'lower');
sw.resolve();

% --- FastSense downsample + render setup (hidden figure) ---
fig = figure('Visible', 'off');
ax = axes('Parent', fig);
fp = FastSense('Parent', ax);
fp.addSensor(sw, 'ShowThresholds', true);
fp.render();
close(fig);
catch
% Warmup is best-effort — never block install on failure
end
end
3 changes: 3 additions & 0 deletions libs/FastSense/build_mex.m
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,7 @@ function build_mex()
% lttb_core_mex.c — Largest-Triangle-Three-Buckets kernel
% violation_cull_mex.c — threshold violation culling
% compute_violations_mex.c — batch violation detection for resolve()
% to_step_function_mex.c — SIMD step-function conversion for thresholds
% build_store_mex.c — bulk SQLite writer for DataStore init
% mksqlite.c — SQLite3 MEX interface (bundled sqlite3.c)
%
Expand Down Expand Up @@ -137,6 +138,7 @@ function build_mex()
'compute_violations_mex.c', 'compute_violations_mex', {{}}, {{}}
'resolve_disk_mex.c', 'resolve_disk_mex', {{sqlite3_src}}, {sqlite3_flags}
'build_store_mex.c', 'build_store_mex', {{sqlite3_src}}, {sqlite3_flags}
'to_step_function_mex.c', 'to_step_function_mex', {{}}, {{}}
};

mksqlite_src = fullfile(rootDir, 'mksqlite.c');
Expand Down Expand Up @@ -228,6 +230,7 @@ function build_mex()
copy_mex_to(outDir, sensorPrivDir, 'violation_cull_mex');
copy_mex_to(outDir, sensorPrivDir, 'compute_violations_mex');
copy_mex_to(outDir, sensorPrivDir, 'resolve_disk_mex');
copy_mex_to(outDir, sensorPrivDir, 'to_step_function_mex');
end

function compile_mex(src_file, out_name, outDir, include_flag, opt_flags, compiler, extra_srcs)
Expand Down
Loading
Loading