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