kernels/custom: grid_sampler_2d fp16 — accumulate in fp32

Match the precision of the portable kernel (after pytorch#19117)
and avoid fp16 catastrophic cancellation on weight computation. The NEON
half variant previously did interpolation weight computation and FMA
accumulation in fp16 via vmul_f16 / vfma_f16; this change loads fp16,
promotes to float32x4 via vcvt_f32_f16, does the four-corner FMA chain in
fp32, and casts back to fp16 on store.

Speed impact: two vcvt per 4-channel group — single-cycle on modern ARM,
unmeasurable at op level in a full-model benchmark (3.5 ms for a typical
call shape, unchanged).

Precision impact: max_abs vs an fp32-then-down-cast reference drops from
~0.1 to 0 on the shapes the polycam depth model uses.

Author: jgibson2

kernels/custom/neon_grid_sampler_2d.cpp

@@ -146,6 +146,8 @@ inline void bilinear_sample_all_channels(
 }
 
 // Process one output spatial location for all channels using NEON, fp16 variant.
+// Math happens in fp32 (matches the portable kernel's precision); loads/stores
+// are fp16. vcvt_f32_f16 / vcvt_f16_f32 are single-cycle on modern ARM.
 inline void bilinear_sample_all_channels_fp16(
     const __fp16* input_n,
     __fp16* output_n,
@@ -161,10 +163,10 @@ inline void bilinear_sample_all_channels_fp16(
     float fx = gx - static_cast<float>(x0);
     float fy = gy - static_cast<float>(y0);
 
-    float16x4_t vw_tl = vdup_n_f16(static_cast<__fp16>((1.0f - fx) * (1.0f - fy)));
-    float16x4_t vw_tr = vdup_n_f16(static_cast<__fp16>(fx * (1.0f - fy)));
-    float16x4_t vw_bl = vdup_n_f16(static_cast<__fp16>((1.0f - fx) * fy));
-    float16x4_t vw_br = vdup_n_f16(static_cast<__fp16>(fx * fy));
+    float32x4_t vw_tl = vdupq_n_f32((1.0f - fx) * (1.0f - fy));
+    float32x4_t vw_tr = vdupq_n_f32(fx * (1.0f - fy));
+    float32x4_t vw_bl = vdupq_n_f32((1.0f - fx) * fy);
+    float32x4_t vw_br = vdupq_n_f32(fx * fy);
 
     bool tl_valid = (static_cast<unsigned>(x0) < static_cast<unsigned>(W_in) &&
                      static_cast<unsigned>(y0) < static_cast<unsigned>(H_in));
@@ -214,18 +216,18 @@ inline void bilinear_sample_all_channels_fp16(
             br[2] = p2[off_br]; br[3] = p3[off_br];
         }
 
-        float16x4_t v_tl = vld1_f16(tl);
-        float16x4_t v_tr = vld1_f16(tr);
-        float16x4_t v_bl = vld1_f16(bl);
-        float16x4_t v_br = vld1_f16(br);
+        float32x4_t v_tl = vcvt_f32_f16(vld1_f16(tl));
+        float32x4_t v_tr = vcvt_f32_f16(vld1_f16(tr));
+        float32x4_t v_bl = vcvt_f32_f16(vld1_f16(bl));
+        float32x4_t v_br = vcvt_f32_f16(vld1_f16(br));
 
-        float16x4_t result = vmul_f16(vw_tl, v_tl);
-        result = vfma_f16(result, vw_tr, v_tr);
-        result = vfma_f16(result, vw_bl, v_bl);
-        result = vfma_f16(result, vw_br, v_br);
+        float32x4_t result = vmulq_f32(vw_tl, v_tl);
+        result = vfmaq_f32(result, vw_tr, v_tr);
+        result = vfmaq_f32(result, vw_bl, v_bl);
+        result = vfmaq_f32(result, vw_br, v_br);
 
         __fp16 res[4];
-        vst1_f16(res, result);
+        vst1_f16(res, vcvt_f16_f32(result));
         output_n[(c + 0) * spatial_out + out_off] = res[0];
         output_n[(c + 1) * spatial_out + out_off] = res[1];
         output_n[(c + 2) * spatial_out + out_off] = res[2];