1 /* 2 * Copyright 2018 Google Inc. 3 * 4 * Use of this source code is governed by a BSD-style license that can be 5 * found in the LICENSE file. 6 */ 7 8 #ifndef SkCubicSolver_DEFINED 9 #define SkCubicSolver_DEFINED 10 11 #include "include/core/SkTypes.h" 12 #include "include/private/SkFloatingPoint.h" 13 14 //#define CUBICMAP_TRACK_MAX_ERROR 15 16 namespace SK_OPTS_NS { 17 eval_poly(float t,float b)18 static float eval_poly(float t, float b) { 19 return b; 20 } 21 22 template <typename... Rest> eval_poly(float t,float m,float b,Rest...rest)23 static float eval_poly(float t, float m, float b, Rest... rest) { 24 return eval_poly(t, sk_fmaf(m,t,b), rest...); 25 } 26 cubic_solver(float A,float B,float C,float D)27 inline float cubic_solver(float A, float B, float C, float D) { 28 #ifdef CUBICMAP_TRACK_MAX_ERROR 29 static int max_iters = 0; 30 #endif 31 32 #ifdef SK_DEBUG 33 auto valid = [](float t) { 34 return t >= 0 && t <= 1; 35 }; 36 #endif 37 38 auto guess_nice_cubic_root = [](float a, float b, float c, float d) { 39 return -d; 40 }; 41 float t = guess_nice_cubic_root(A, B, C, D); 42 43 int iters = 0; 44 const int MAX_ITERS = 8; 45 for (; iters < MAX_ITERS; ++iters) { 46 SkASSERT(valid(t)); 47 float f = eval_poly(t, A,B,C,D); // f = At^3 + Bt^2 + Ct + D 48 if (sk_float_abs(f) <= 0.00005f) { 49 break; 50 } 51 float fp = eval_poly(t, 3*A, 2*B, C); // f' = 3At^2 + 2Bt + C 52 float fpp = eval_poly(t, 3*A+3*A, 2*B); // f'' = 6At + 2B 53 54 float numer = 2 * fp * f; 55 float denom = sk_fmaf(2*fp, fp, -(f*fpp)); 56 57 t -= numer / denom; 58 } 59 60 #ifdef CUBICMAP_TRACK_MAX_ERROR 61 if (max_iters < iters) { 62 max_iters = iters; 63 SkDebugf("max_iters %d\n", max_iters); 64 } 65 #endif 66 SkASSERT(valid(t)); 67 return t; 68 } 69 70 } // namespace SK_OPTS_NS 71 #endif 72