convert poly into double

Author: hulxv

libc/src/math/generic/asinpif16.cpp

@@ -93,22 +93,22 @@ LLVM_LIBC_FUNCTION(float16, asinpif16, (float16 x)) {
   //
   // it's very accurate in the range [0, 0.5] and has a maximum error of
   // 0.0000000000000001 in the range [0, 0.5].
-  static constexpr float POLY_COEFFS[10] = {
-      0.318309886183791f,   // x^1
-      0.0530516476972984f,  // x^3
-      0.0238732414637843f,  // x^5
-      0.0142102627760621f,  // x^7
-      0.00967087327815336f, // x^9
-      0.00712127941391293f, // x^11
-      0.00552355646848375f, // x^13
-      0.00444514782463692f, // x^15
-      0.00367705242846804f, // x^17
-      0.00310721681820837f  // x^19
+  static constexpr double POLY_COEFFS[10] = {
+      0.318309886183791,   // x^1
+      0.0530516476972984,  // x^3
+      0.0238732414637843,  // x^5
+      0.0142102627760621,  // x^7
+      0.00967087327815336, // x^9
+      0.00712127941391293, // x^11
+      0.00552355646848375, // x^13
+      0.00444514782463692, // x^15
+      0.00367705242846804, // x^17
+      0.00310721681820837  // x^19
   };
 
   // polynomial evaluation using horner's method
   // work only for |x| in [0, 0.5]
-  auto asinpi_polyeval = [](float xsq) -> float {
+  auto asinpi_polyeval = [](double xsq) -> double {
     return fputil::polyeval(xsq, POLY_COEFFS[0], POLY_COEFFS[1], POLY_COEFFS[2],
                             POLY_COEFFS[3], POLY_COEFFS[4], POLY_COEFFS[5],
                             POLY_COEFFS[6], POLY_COEFFS[7], POLY_COEFFS[8],
@@ -118,8 +118,9 @@ LLVM_LIBC_FUNCTION(float16, asinpif16, (float16 x)) {
   // if |x| <= 0.5:
   if (LIBC_UNLIKELY(x_abs <= ONE_OVER_TWO)) {
     // Use polynomial approximation of asin(x)/pi in the range [0, 0.5]
-    float16 xsq = x * x;
-    float result = x * asinpi_polyeval(xsq);
+    double xd = static_cast<double>(x);
+    double xsq = xd * xd;
+    double result = xd * asinpi_polyeval(xsq);
     return fputil::cast<float16>(signed_result(result));
   }
 
@@ -150,11 +151,11 @@ LLVM_LIBC_FUNCTION(float16, asinpif16, (float16 x)) {
   //             = 0.5 - 0.5 * x
   //             = multiply_add(-0.5, x, 0.5)
 
-  float u = static_cast<float>(
+  double u = static_cast<double>(
       fputil::multiply_add(-ONE_OVER_TWO, x_abs, ONE_OVER_TWO));
-  float u_sqrt = fputil::sqrt<float>(static_cast<float>(u));
-  float asinpi_sqrt_u = u_sqrt * asinpi_polyeval(u);
-  float result = fputil::multiply_add(-2.0f16, asinpi_sqrt_u, ONE_OVER_TWO);
+  double u_sqrt = fputil::sqrt<double>(static_cast<double>(u));
+  double asinpi_sqrt_u = u_sqrt * asinpi_polyeval(u);
+  double result = fputil::multiply_add(-2.0f16, asinpi_sqrt_u, ONE_OVER_TWO);
 
   return fputil::cast<float16>(signed_result(result));
 }