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/base/SkFloatingPoint.h" 13 14 namespace SK_OPTS_NS { 15 eval_poly(float t,float b)16 static float eval_poly(float t, float b) { 17 return b; 18 } 19 20 template <typename... Rest> eval_poly(float t,float m,float b,Rest...rest)21 static float eval_poly(float t, float m, float b, Rest... rest) { 22 return eval_poly(t, sk_fmaf(m,t,b), rest...); 23 } 24 cubic_solver(float A,float B,float C,float D)25 inline float cubic_solver(float A, float B, float C, float D) { 26 27 #ifdef SK_DEBUG 28 auto valid = [](float t) { 29 return t >= 0 && t <= 1; 30 }; 31 #endif 32 33 auto guess_nice_cubic_root = [](float a, float b, float c, float d) { 34 return -d; 35 }; 36 float t = guess_nice_cubic_root(A, B, C, D); 37 38 int iters = 0; 39 const int MAX_ITERS = 8; 40 for (; iters < MAX_ITERS; ++iters) { 41 SkASSERT(valid(t)); 42 float f = eval_poly(t, A,B,C,D); // f = At^3 + Bt^2 + Ct + D 43 if (sk_float_abs(f) <= 0.00005f) { 44 break; 45 } 46 float fp = eval_poly(t, 3*A, 2*B, C); // f' = 3At^2 + 2Bt + C 47 float fpp = eval_poly(t, 3*A+3*A, 2*B); // f'' = 6At + 2B 48 49 float numer = 2 * fp * f; 50 float denom = sk_fmaf(2*fp, fp, -(f*fpp)); 51 52 t -= numer / denom; 53 } 54 55 SkASSERT(valid(t)); 56 return t; 57 } 58 59 } // namespace SK_OPTS_NS 60 #endif 61