xref: /aosp_15_r20/external/XNNPACK/src/math/sqrt-u32-scalar-clz-binsearch.c (revision 4bdc94577ba0e567308109d787f7fec7b531ce36)
1 // Copyright 2022 Google LLC
2 //
3 // This source code is licensed under the BSD-style license found in the
4 // LICENSE file in the root directory of this source tree.
5 
6 #include <assert.h>
7 #include <stddef.h>
8 
9 #include <xnnpack/math.h>
10 #include <xnnpack/math-stubs.h>
11 
12 
xnn_math_u32_sqrt__scalar_clz_binsearch(size_t n,const uint32_t * input,uint32_t * output)13 void xnn_math_u32_sqrt__scalar_clz_binsearch(
14     size_t n,
15     const uint32_t* input,
16     uint32_t* output)
17 {
18   assert(n % sizeof(uint32_t) == 0);
19 
20   for (; n != 0; n -= sizeof(uint32_t)) {
21     const uint32_t vx = *input++;
22 
23     // Based on Hacker's Delight, Figure 11-3.
24     uint32_t vb = (UINT32_C(1) << ((33 - math_clz_u32(vx)) / 2)) - 1;
25     uint32_t va = (vb + 3) / 2;
26     do {
27       const uint32_t vm = (va + vb) >> 1;
28       assert(vm <= UINT32_C(65535));
29       if XNN_UNPREDICTABLE(vm * vm > vx) {
30         vb = vm - 1;
31       } else {
32         va = vm + 1;
33       }
34     } while XNN_LIKELY(vb >= va);
35 
36     uint32_t vy = va - 1;
37     // vy is sqrt(vx) rounded down. Do the final rounding up if needed.
38     if XNN_UNPREDICTABLE(va * vy < vx) {
39       vy += 1;
40     }
41 
42     *output++ = vy;
43   }
44 }
45