benchmark pdist2

## Benchmark script for pdist2: original vs blocked
## Outputs CSV with N, t_orig, t_blocked, speedup

clear all;

## Add paths for implementations
addpath("/tmp");
addpath("inst");

## Create renamed copies to avoid function name conflicts
copyfile("/tmp/pdist2_orig.m", "/tmp/pdist2_original.m");
copyfile("/tmp/pdist2_blocked.m", "/tmp/pdist2_new.m");

## Modify function names in copied files
fid = fopen("/tmp/pdist2_original.m", "r");
content = char(fread(fid)');
fclose(fid);
content = strrep(content, "function [D, I] = pdist2 (", "function [D, I] = pdist2_original (");
fid = fopen("/tmp/pdist2_original.m", "w");
fwrite(fid, content);
fclose(fid);

fid = fopen("/tmp/pdist2_new.m", "r");
content = char(fread(fid)');
fclose(fid);
content = strrep(content, "function [D, I] = pdist2 (", "function [D, I] = pdist2_new (");
fid = fopen("/tmp/pdist2_new.m", "w");
fwrite(fid, content);
fclose(fid);

## Clear function cache and reload
clear functions;
rehash;

## Test N values: small values below threshold (500), then larger values above
N_values = [5, 10, 20, 30, 50, 75, 100, 150, 200, 250, 300, 350, 400, 450, 500, 600, 700, 800, 1000, 1500, 2000, 3000, 4000, 5000];
P = 50;  # dimensions
n_runs = 5;  # number of timed runs
n_warmup = 3;  # warmup runs

## Open CSV file
fid = fopen("benchmark_pdist2_results.csv", "w");
fprintf(fid, "N,t_orig,t_blocked,speedup\n");

printf("Benchmarking pdist2: original vs blocked (threshold=500)\n");
printf("=========================================================\n");
printf("%-6s  %-12s  %-12s  %-8s\n", "N", "t_orig", "t_blocked", "speedup");
printf("------  ------------  ------------  --------\n");

for n_idx = 1:length(N_values)
  N = N_values(n_idx);
  X = rand(N, P);
  Y = rand(N, P);  # square case M=N

  ## Warmup runs for original
  for w = 1:n_warmup
    D = pdist2_original(X, Y);
  endfor

  ## Timed runs for original
  times_orig = zeros(1, n_runs);
  for r = 1:n_runs
    tic;
    D = pdist2_original(X, Y);
    times_orig(r) = toc;
  endfor
  t_orig = mean(times_orig);

  ## Warmup runs for blocked
  for w = 1:n_warmup
    D = pdist2_new(X, Y);
  endfor

  ## Timed runs for blocked
  times_blocked = zeros(1, n_runs);
  for r = 1:n_runs
    tic;
    D = pdist2_new(X, Y);
    times_blocked(r) = toc;
  endfor
  t_blocked = mean(times_blocked);

  speedup = t_orig / t_blocked;

  printf("%-6d  %-12.6f  %-12.6f  %-8.4f\n", N, t_orig, t_blocked, speedup);
  fprintf(fid, "%d,%.6f,%.6f,%.4f\n", N, t_orig, t_blocked, speedup);
endfor

fclose(fid);
printf("\nResults saved to benchmark_pdist2_results.csv\n");

## Cleanup
delete("/tmp/pdist2_original.m");
delete("/tmp/pdist2_new.m");