fix SkJumper radial gradient precision

rcp(rsqrt(x)) doesn't have enough precision when x is a coordinate.
(It's fine when x is a color, like in the softlight blend mode.)

Adds a GM to test this.  It used to look quite ugly.

Change-Id: Icec295c2e2f50ae7a5e3e33c62270f632a58f65c
Reviewed-on: https://skia-review.googlesource.com/16914
Reviewed-by: Herb Derby <herb@google.com>
Commit-Queue: Mike Klein <mtklein@chromium.org>

gm/radial_gradient_precision.cpp

@@ -0,0 +1,22 @@
+/*
+ * Copyright 2017 Google Inc.
+ *
+ * Use of this source code is governed by a BSD-style license that can be
+ * found in the LICENSE file.
+ */
+
+#include "gm.h"
+#include "SkGradientShader.h"
+
+// All we're looking for here is that we see a smooth gradient.
+DEF_SIMPLE_GM(radial_gradient_precision, canvas, 200, 200) {
+    SkPoint  center   = {100000, 100000};
+    SkScalar radius   = 40;
+    SkColor  colors[] = {SK_ColorBLACK, SK_ColorWHITE};
+
+    SkPaint p;
+    p.setShader(SkGradientShader::MakeRadial(center, radius,
+                                             colors, nullptr, SK_ARRAY_COUNT(colors),
+                                             SkShader::kRepeat_TileMode));
+    canvas->drawPaint(p);
+}

gn/gm.gni

@@ -235,6 +235,7 @@ gm_sources = [
   "$_gm/poly2poly.cpp",
   "$_gm/polygons.cpp",
   "$_gm/quadpaths.cpp",
+  "$_gm/radial_gradient_precision.cpp",
   "$_gm/readpixels.cpp",
   "$_gm/recordopts.cpp",
   "$_gm/rectangletexture.cpp",

src/jumper/SkJumper_stages.cpp

@@ -1147,7 +1147,7 @@ STAGE(xy_to_unit_angle) {
 STAGE(xy_to_radius) {
     F X2 = r * r,
       Y2 = g * g;
-    r = rcp(rsqrt(X2 + Y2));
+    r = sqrt_(X2 + Y2);
 }
 
 STAGE(save_xy) {

src/jumper/SkJumper_vectors.h

@@ -34,6 +34,7 @@
     SI F   floor_(F v)          { return floorf(v); }
     SI F   rcp   (F v)          { return 1.0f / v; }
     SI F   rsqrt (F v)          { return 1.0f / sqrtf(v); }
+    SI F    sqrt_(F v)          { return sqrtf(v); }
     SI U32 round (F v, F scale) { return (uint32_t)lrintf(v*scale); }
     SI U16 pack(U32 v)          { return (U16)v; }
     SI U8  pack(U16 v)          { return  (U8)v; }
@@ -94,6 +95,7 @@
     SI F   floor_(F v)                           { return vrndmq_f32(v);           }
     SI F   rcp   (F v) { auto e = vrecpeq_f32 (v); return vrecpsq_f32 (v,e  ) * e; }
     SI F   rsqrt (F v) { auto e = vrsqrteq_f32(v); return vrsqrtsq_f32(v,e*e) * e; }
+    SI F    sqrt_(F v)                           { return vsqrtq_f32(v); }
     SI U32 round (F v, F scale)                  { return vcvtnq_u32_f32(v*scale); }
     SI U16 pack(U32 v)                           { return __builtin_convertvector(v, U16); }
     SI U8  pack(U16 v)                           { return __builtin_convertvector(v,  U8); }
@@ -158,6 +160,13 @@
     SI U16 pack(U32 v)                         { return __builtin_convertvector(v, U16); }
     SI U8  pack(U16 v)                         { return __builtin_convertvector(v,  U8); }
 
+    SI F sqrt_(F v) {
+        auto e = vrsqrte_f32(v);  // Estimate and two refinement steps for e = rsqrt(v).
+        e *= vrsqrts_f32(v,e*e);
+        e *= vrsqrts_f32(v,e*e);
+        return v*e;               // sqrt(v) == v*rsqrt(v).
+    }
+
     SI F if_then_else(I32 c, F t, F e) { return vbsl_f32((U32)c,t,e); }
 
     SI F floor_(F v) {
@@ -236,6 +245,7 @@
     SI F   floor_(F v)          { return _mm256_floor_ps(v);    }
     SI F   rcp   (F v)          { return _mm256_rcp_ps  (v);    }
     SI F   rsqrt (F v)          { return _mm256_rsqrt_ps(v);    }
+    SI F    sqrt_(F v)          { return _mm256_sqrt_ps (v);    }
     SI U32 round (F v, F scale) { return _mm256_cvtps_epi32(v*scale); }
 
     SI U16 pack(U32 v) {
@@ -438,8 +448,9 @@
     SI F   min(F a, F b)       { return _mm_min_ps(a,b);    }
     SI F   max(F a, F b)       { return _mm_max_ps(a,b);    }
     SI F   abs_(F v)           { return _mm_and_ps(v, 0-v); }
-    SI F   rcp  (F v)          { return _mm_rcp_ps  (v);    }
-    SI F   rsqrt(F v)          { return _mm_rsqrt_ps(v);    }
+    SI F   rcp   (F v)         { return _mm_rcp_ps  (v);    }
+    SI F   rsqrt (F v)         { return _mm_rsqrt_ps(v);    }
+    SI F    sqrt_(F v)         { return _mm_sqrt_ps (v);    }
     SI U32 round(F v, F scale) { return _mm_cvtps_epi32(v*scale); }
 
     SI U16 pack(U32 v) {