27accef223
Review URL: https://codereview.appspot.com/5576043 git-svn-id: http://skia.googlecode.com/svn/trunk@3087 2bbb7eff-a529-9590-31e7-b0007b416f81
65 lines
1.6 KiB
C++
65 lines
1.6 KiB
C++
#include "CubicUtilities.h"
|
|
#include "DataTypes.h"
|
|
#include "QuadraticUtilities.h"
|
|
|
|
const double PI = 4 * atan(1);
|
|
|
|
static bool is_unit_interval(double x) {
|
|
return x > 0 && x < 1;
|
|
}
|
|
|
|
// from SkGeometry.cpp (and Numeric Solutions, 5.6)
|
|
int cubicRoots(double A, double B, double C, double D, double t[3]) {
|
|
if (approximately_zero(A)) { // we're just a quadratic
|
|
return quadraticRoots(B, C, D, t);
|
|
}
|
|
double a, b, c;
|
|
{
|
|
double invA = 1 / A;
|
|
a = B * invA;
|
|
b = C * invA;
|
|
c = D * invA;
|
|
}
|
|
double a2 = a * a;
|
|
double Q = (a2 - b * 3) / 9;
|
|
double R = (2 * a2 * a - 9 * a * b + 27 * c) / 54;
|
|
double Q3 = Q * Q * Q;
|
|
double R2MinusQ3 = R * R - Q3;
|
|
double adiv3 = a / 3;
|
|
double* roots = t;
|
|
double r;
|
|
|
|
if (R2MinusQ3 < 0) // we have 3 real roots
|
|
{
|
|
double theta = acos(R / sqrt(Q3));
|
|
double neg2RootQ = -2 * sqrt(Q);
|
|
|
|
r = neg2RootQ * cos(theta / 3) - adiv3;
|
|
if (is_unit_interval(r))
|
|
*roots++ = r;
|
|
|
|
r = neg2RootQ * cos((theta + 2 * PI) / 3) - adiv3;
|
|
if (is_unit_interval(r))
|
|
*roots++ = r;
|
|
|
|
r = neg2RootQ * cos((theta - 2 * PI) / 3) - adiv3;
|
|
if (is_unit_interval(r))
|
|
*roots++ = r;
|
|
}
|
|
else // we have 1 real root
|
|
{
|
|
double A = fabs(R) + sqrt(R2MinusQ3);
|
|
A = cube_root(A);
|
|
if (R > 0) {
|
|
A = -A;
|
|
}
|
|
if (A != 0) {
|
|
A += Q / A;
|
|
}
|
|
r = A - adiv3;
|
|
if (is_unit_interval(r))
|
|
*roots++ = r;
|
|
}
|
|
return (int)(roots - t);
|
|
}
|