CLAUDE: peaks compare - dual criterion (abs floor OR 4 float ULP at operand scale)

First rung-2 real case (Andrey's 07/14 export) measured the predicted tier
precisely: tasks + corr_td + corr_pd ALL BIT-EXACT @tol 0 (the in-place
FZ-norm is production-identical); peaks ~14% of elements flip by 1-2 float
ULP at their own magnitude (double-math FMA contraction vs Java, absorbed
except at float cast boundaries) - a flat absolute tol cannot cover 1 ULP
of lambda~40 (3.8e-6) while staying meaningful for dx~0.01. Criterion now:
pass if within --peak-tol ABSOLUTE (near-zero floor, default 1e-6) OR
within 4 float ULP of the operand (measured max 2 ULP; near-zero relative
outliers ~21 ULP are 3e-9 absolute = caught by the floor). Verdict on the
fresh case: ALL PASS, viol=0.
Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>

jna/test_pose_corr_jna.cu

@@ -108,7 +108,8 @@ int main(int argc, char** argv){
   std::string srcdir   = "src"; // cwd = repo root; run_cases.sh passes --src explicitly // By Claude on 07/13/2026
   std::string devrt    = "/usr/local/cuda-12.8/lib64/libcudadevrt.a";
   float tol = 0.0f;
-  float peak_tol = 1e-6f; // Java-double-oracle tier for the peak kernel (header note) // By Claude on 07/13/2026
+  float peak_tol = 1e-6f; // Java-double-oracle tier: ABSOLUTE floor (near-zero components) // By Claude on 07/13/2026
+  #define PEAK_TOL_ULP 4  // ...OR within this many float ULP at operand scale (measured max 2 ULP, 07/14/2026) // By Claude on 07/14/2026
   int list_only = 0;
   for (int i = 1; i < argc; i++){
     if      (!strcmp(argv[i],"--data")  && (i+1<argc)) data_dir = argv[++i];
@@ -288,7 +289,7 @@ int main(int argc, char** argv){
       // match by packed index (dense tile -> (tile<<8)|0xff sum slot), oracle row order out
       std::vector<float> res_pd((size_t)n_pd*pd_size, NAN);
       std::vector<float> res_pk((size_t)n_pd*7, NAN);
-      double pd_maxd = 0, pk_maxd = 0; long pd_nan = 0, pk_nan = 0, pd_missing = 0;
+      double pd_maxd = 0, pk_maxd = 0, pk_maxulp = 0; long pd_nan = 0, pk_nan = 0, pd_missing = 0, pk_viol = 0;
       for (int i = 0; i < n_pd; i++){
         const int packed = (pd_idx[i] << 8) | 0xff; // CORR_NTILE_SHIFT
         std::map<int,int>::const_iterator it = by_index.find(packed);
@@ -310,14 +311,24 @@ int main(int argc, char** argv){
           if (gn || en){ if(gn!=en) pk_nan++; continue; }
           const double d = std::fabs((double)gpk[j]-(double)epk[j]);
           if (d > pk_maxd) pk_maxd = d;
+          // dual criterion (measured 07/14/2026 on the first rung-2 case: ~14% of
+          // elements flip by 1-2 float ULP AT OPERAND SCALE - double-math FMA
+          // contraction vs Java, absorbed except at float cast boundaries; an
+          // absolute tol can't cover 1 ULP of lambda~40 = 3.8e-6 while staying
+          // meaningful for dx~0.01): pass if within peak_tol absolute (near-zero
+          // floor) OR within PEAK_TOL_ULP float ULP of the operand magnitude.
+          const double ulp = (double)(std::nextafterf(std::fabs(epk[j]), INFINITY) - std::fabs(epk[j]));
+          const double du = (ulp > 0) ? d/ulp : 0;
+          if (du > pk_maxulp) pk_maxulp = du;
+          if ((d > peak_tol) && (du > PEAK_TOL_ULP)) pk_viol++;
         }
       }
-      const int pd_fail = (pd_maxd>tol)      || (pd_nan>0) || (pd_missing>0);
-      const int pk_fail = (pk_maxd>peak_tol) || (pk_nan>0) || (pd_missing>0);
+      const int pd_fail = (pd_maxd>tol) || (pd_nan>0) || (pd_missing>0);
+      const int pk_fail = (pk_viol>0)   || (pk_nan>0) || (pd_missing>0);
       printf("it%d: compare 'corr_pd' (%d tiles, in-place FZ-norm): max|diff|=%g NaN-mismatch=%ld missing=%ld -> %s\n",
              k, n_pd, pd_maxd, pd_nan, pd_missing, pd_fail?"FAIL":"PASS");
-      printf("it%d: compare 'peaks' (%d tiles vs the double Java oracle): max|diff|=%g NaN-mismatch=%ld -> %s (peak_tol %g)\n",
-             k, n_pd, pk_maxd, pk_nan, pk_fail?"FAIL":"PASS", peak_tol);
+      printf("it%d: compare 'peaks' (%d tiles vs the double Java oracle): max|diff|=%g max=%.2f ULP viol=%ld NaN-mismatch=%ld -> %s (abs %g | %d ULP)\n",
+             k, n_pd, pk_maxd, pk_maxulp, pk_viol, pk_nan, pk_fail?"FAIL":"PASS", peak_tol, PEAK_TOL_ULP);
       if (pd_fail || pk_fail) fail = 1;
       cudaMalloc((void**)&scratch, res_pd.size()*sizeof(float));
       cudaMemcpy(scratch, res_pd.data(), res_pd.size()*sizeof(float), cudaMemcpyHostToDevice);

src/tests/test_pose_corr.cu

@@ -134,7 +134,8 @@ int main(int argc, char **argv)
 {
 	std::string data_dir = "testdata/pose_corr";
 	float tol       = 0.0f;   // the whole chain is expected bit-exact vs the oracle
-	float peak_tol  = 1e-6f;  // Java-double-oracle tier for the peak kernel (header note)
+	float peak_tol  = 1e-6f;  // Java-double-oracle tier: ABSOLUTE floor (near-zero components)
+	#define PEAK_TOL_ULP 4    // ...OR within this many float ULP at operand scale (measured max 2 ULP, 07/14/2026) // By Claude on 07/14/2026
 	int   list_only = 0;
 	for (int i = 1; i < argc; i++) {
 		if      (!strcmp(argv[i], "--data") && (i + 1 < argc)) data_dir = argv[++i];
@@ -506,8 +507,8 @@ int main(int argc, char **argv)
 			// the checkpoint-1 map; save results in ORACLE row order for the Java round trip
 			std::vector<float> res_pd((size_t) n_pd * pd_size, NAN);
 			std::vector<float> res_pk((size_t) n_pd * 7, NAN);
-			double pd_maxd = 0, pk_maxd = 0;
-			long pd_nan = 0, pk_nan = 0, pd_missing = 0;
+			double pd_maxd = 0, pk_maxd = 0, pk_maxulp = 0;
+			long pd_nan = 0, pk_nan = 0, pd_missing = 0, pk_viol = 0;
 			for (int i = 0; i < n_pd; i++) {
 				const int packed = (pd_idx[i] << CORR_NTILE_SHIFT) | 0xff;
 				const std::map<int, int>::const_iterator it = got_by_index.find(packed);
@@ -529,15 +530,22 @@ int main(int argc, char **argv)
 					if (gn || en) { if (gn != en) pk_nan++; continue; }
 					const double d = std::fabs((double) gpk[j] - (double) epk[j]);
 					if (d > pk_maxd) pk_maxd = d;
+					// dual criterion (see test_pose_corr_jna.cu, measured 07/14/2026):
+					// pass if within peak_tol ABSOLUTE (near-zero floor) OR within
+					// PEAK_TOL_ULP float ULP at operand scale (cast-boundary flips)
+					const double ulp = (double) (std::nextafterf(std::fabs(epk[j]), INFINITY) - std::fabs(epk[j]));
+					const double du = (ulp > 0) ? d / ulp : 0;
+					if (du > pk_maxulp) pk_maxulp = du;
+					if ((d > peak_tol) && (du > PEAK_TOL_ULP)) pk_viol++;
 				}
 			}
-			const int pd_fail = (pd_maxd > tol)      || (pd_nan > 0) || (pd_missing > 0);
-			const int pk_fail = (pk_maxd > peak_tol) || (pk_nan > 0) || (pd_missing > 0);
+			const int pd_fail = (pd_maxd > tol) || (pd_nan > 0) || (pd_missing > 0);
+			const int pk_fail = (pk_viol > 0)   || (pk_nan > 0) || (pd_missing > 0);
 			printf("it%d: compare 'corr_pd' (%d tiles, in-place FZ-norm): max|diff|=%g NaN-mismatch=%ld "
 					"missing=%ld -> %s\n", k, n_pd, pd_maxd, pd_nan, pd_missing, pd_fail ? "FAIL" : "PASS");
-			printf("it%d: compare 'peaks' (%d tiles vs the double Java oracle): max|diff|=%g "
-					"NaN-mismatch=%ld -> %s (peak_tol %g)\n",
-					k, n_pd, pk_maxd, pk_nan, pk_fail ? "FAIL" : "PASS", peak_tol);
+			printf("it%d: compare 'peaks' (%d tiles vs the double Java oracle): max|diff|=%g max=%.2f ULP "
+					"viol=%ld NaN-mismatch=%ld -> %s (abs %g | %d ULP)\n",
+					k, n_pd, pk_maxd, pk_maxulp, pk_viol, pk_nan, pk_fail ? "FAIL" : "PASS", peak_tol, PEAK_TOL_ULP);
 			if (pd_fail || pk_fail) fail = 1;
 			checkCudaErrors(cudaMemcpy(gpu_result_scratch.dev(), res_pd.data(),
 					res_pd.size() * sizeof(float), cudaMemcpyHostToDevice));