| OLD | NEW |
| (Empty) |
| 1 #include "CubicUtilities.h" | |
| 2 #include "Intersection_Tests.h" | |
| 3 #include "QuadraticUtilities.h" | |
| 4 #include "QuarticRoot.h" | |
| 5 | |
| 6 double mulA[] = {-3, -1, 1, 3}; | |
| 7 size_t mulACount = sizeof(mulA) / sizeof(mulA[0]); | |
| 8 double rootB[] = {-9, -6, -3, -1, 0, 1, 3, 6, 9}; | |
| 9 size_t rootBCount = sizeof(rootB) / sizeof(rootB[0]); | |
| 10 double rootC[] = {-8, -6, -2, -1, 0, 1, 2, 6, 8}; | |
| 11 size_t rootCCount = sizeof(rootC) / sizeof(rootC[0]); | |
| 12 double rootD[] = {-7, -4, -1, 0, 1, 2, 5}; | |
| 13 size_t rootDCount = sizeof(rootD) / sizeof(rootD[0]); | |
| 14 double rootE[] = {-5, -1, 0, 1, 7}; | |
| 15 size_t rootECount = sizeof(rootE) / sizeof(rootE[0]); | |
| 16 | |
| 17 | |
| 18 static void quadraticTest(bool limit) { | |
| 19 // (x - a)(x - b) == x^2 - (a + b)x + ab | |
| 20 for (size_t aIndex = 0; aIndex < mulACount; ++aIndex) { | |
| 21 for (size_t bIndex = 0; bIndex < rootBCount; ++bIndex) { | |
| 22 for (size_t cIndex = 0; cIndex < rootCCount; ++cIndex) { | |
| 23 const double A = mulA[aIndex]; | |
| 24 double B = rootB[bIndex]; | |
| 25 double C = rootC[cIndex]; | |
| 26 if (limit) { | |
| 27 B = (B - 6) / 12; | |
| 28 C = (C - 6) / 12; | |
| 29 } | |
| 30 const double b = A * (B + C); | |
| 31 const double c = A * B * C; | |
| 32 double roots[2]; | |
| 33 const int rootCount = limit ? quadraticRootsValidT(A, b, c, root
s) | |
| 34 : quadraticRootsReal(A, b, c, roots); | |
| 35 int expected; | |
| 36 if (limit) { | |
| 37 expected = B <= 0 && B >= -1; | |
| 38 expected += B != C && C <= 0 && C >= -1; | |
| 39 } else { | |
| 40 expected = 1 + (B != C); | |
| 41 } | |
| 42 SkASSERT(rootCount == expected); | |
| 43 if (!rootCount) { | |
| 44 continue; | |
| 45 } | |
| 46 SkASSERT(approximately_equal(roots[0], -B) | |
| 47 || approximately_equal(roots[0], -C)); | |
| 48 if (expected > 1) { | |
| 49 SkASSERT(!approximately_equal(roots[0], roots[1])); | |
| 50 SkASSERT(approximately_equal(roots[1], -B) | |
| 51 || approximately_equal(roots[1], -C)); | |
| 52 } | |
| 53 } | |
| 54 } | |
| 55 } | |
| 56 } | |
| 57 | |
| 58 static void testOneCubic(bool limit, size_t aIndex, size_t bIndex, size_t cIndex
, size_t dIndex) { | |
| 59 const double A = mulA[aIndex]; | |
| 60 double B = rootB[bIndex]; | |
| 61 double C = rootC[cIndex]; | |
| 62 double D = rootD[dIndex]; | |
| 63 if (limit) { | |
| 64 B = (B - 6) / 12; | |
| 65 C = (C - 6) / 12; | |
| 66 D = (C - 2) / 6; | |
| 67 } | |
| 68 const double b = A * (B + C + D); | |
| 69 const double c = A * (B * C + C * D + B * D); | |
| 70 const double d = A * B * C * D; | |
| 71 double roots[3]; | |
| 72 const int rootCount = limit ? cubicRootsValidT(A, b, c, d, roots) | |
| 73 : cubicRootsReal(A, b, c, d, roots); | |
| 74 int expected; | |
| 75 if (limit) { | |
| 76 expected = B <= 0 && B >= -1; | |
| 77 expected += B != C && C <= 0 && C >= -1; | |
| 78 expected += B != D && C != D && D <= 0 && D >= -1; | |
| 79 } else { | |
| 80 expected = 1 + (B != C) + (B != D && C != D); | |
| 81 } | |
| 82 SkASSERT(rootCount == expected); | |
| 83 if (!rootCount) { | |
| 84 return; | |
| 85 } | |
| 86 SkASSERT(approximately_equal(roots[0], -B) | |
| 87 || approximately_equal(roots[0], -C) | |
| 88 || approximately_equal(roots[0], -D)); | |
| 89 if (expected <= 1) { | |
| 90 return; | |
| 91 } | |
| 92 SkASSERT(!approximately_equal(roots[0], roots[1])); | |
| 93 SkASSERT(approximately_equal(roots[1], -B) | |
| 94 || approximately_equal(roots[1], -C) | |
| 95 || approximately_equal(roots[1], -D)); | |
| 96 if (expected <= 2) { | |
| 97 return; | |
| 98 } | |
| 99 SkASSERT(!approximately_equal(roots[0], roots[2]) | |
| 100 && !approximately_equal(roots[1], roots[2])); | |
| 101 SkASSERT(approximately_equal(roots[2], -B) | |
| 102 || approximately_equal(roots[2], -C) | |
| 103 || approximately_equal(roots[2], -D)); | |
| 104 } | |
| 105 | |
| 106 static void cubicTest(bool limit) { | |
| 107 // (x - a)(x - b)(x - c) == x^3 - (a + b + c)x^2 + (ab + bc + ac)x - abc | |
| 108 for (size_t aIndex = 0; aIndex < mulACount; ++aIndex) { | |
| 109 for (size_t bIndex = 0; bIndex < rootBCount; ++bIndex) { | |
| 110 for (size_t cIndex = 0; cIndex < rootCCount; ++cIndex) { | |
| 111 for (size_t dIndex = 0; dIndex < rootDCount; ++dIndex) { | |
| 112 testOneCubic(limit, aIndex, bIndex, cIndex, dIndex); | |
| 113 } | |
| 114 } | |
| 115 } | |
| 116 } | |
| 117 } | |
| 118 | |
| 119 static void testOneQuartic(size_t aIndex, size_t bIndex, size_t cIndex, size_t d
Index, | |
| 120 size_t eIndex) { | |
| 121 const double A = mulA[aIndex]; | |
| 122 const double B = rootB[bIndex]; | |
| 123 const double C = rootC[cIndex]; | |
| 124 const double D = rootD[dIndex]; | |
| 125 const double E = rootE[eIndex]; | |
| 126 const double b = A * (B + C + D + E); | |
| 127 const double c = A * (B * C + C * D + B * D + B * E + C * E + D * E); | |
| 128 const double d = A * (B * C * D + B * C * E + B * D * E + C * D * E); | |
| 129 const double e = A * B * C * D * E; | |
| 130 double roots[4]; | |
| 131 bool oneHint = approximately_zero(A + b + c + d + e); | |
| 132 int rootCount = reducedQuarticRoots(A, b, c, d, e, oneHint, roots); | |
| 133 if (rootCount < 0) { | |
| 134 rootCount = quarticRootsReal(0, A, b, c, d, e, roots); | |
| 135 } | |
| 136 const int expected = 1 + (B != C) + (B != D && C != D) + (B != E && C != E &
& D != E); | |
| 137 SkASSERT(rootCount == expected); | |
| 138 SkASSERT(AlmostEqualUlps(roots[0], -B) | |
| 139 || AlmostEqualUlps(roots[0], -C) | |
| 140 || AlmostEqualUlps(roots[0], -D) | |
| 141 || AlmostEqualUlps(roots[0], -E)); | |
| 142 if (expected <= 1) { | |
| 143 return; | |
| 144 } | |
| 145 SkASSERT(!AlmostEqualUlps(roots[0], roots[1])); | |
| 146 SkASSERT(AlmostEqualUlps(roots[1], -B) | |
| 147 || AlmostEqualUlps(roots[1], -C) | |
| 148 || AlmostEqualUlps(roots[1], -D) | |
| 149 || AlmostEqualUlps(roots[1], -E)); | |
| 150 if (expected <= 2) { | |
| 151 return; | |
| 152 } | |
| 153 SkASSERT(!AlmostEqualUlps(roots[0], roots[2]) | |
| 154 && !AlmostEqualUlps(roots[1], roots[2])); | |
| 155 SkASSERT(AlmostEqualUlps(roots[2], -B) | |
| 156 || AlmostEqualUlps(roots[2], -C) | |
| 157 || AlmostEqualUlps(roots[2], -D) | |
| 158 || AlmostEqualUlps(roots[2], -E)); | |
| 159 if (expected <= 3) { | |
| 160 return; | |
| 161 } | |
| 162 SkASSERT(!AlmostEqualUlps(roots[0], roots[3]) | |
| 163 && !AlmostEqualUlps(roots[1], roots[3]) | |
| 164 && !AlmostEqualUlps(roots[2], roots[3])); | |
| 165 SkASSERT(AlmostEqualUlps(roots[3], -B) | |
| 166 || AlmostEqualUlps(roots[3], -C) | |
| 167 || AlmostEqualUlps(roots[3], -D) | |
| 168 || AlmostEqualUlps(roots[3], -E)); | |
| 169 } | |
| 170 | |
| 171 static void quarticTest() { | |
| 172 // (x - a)(x - b)(x - c)(x - d) == x^4 - (a + b + c + d)x^3 | |
| 173 // + (ab + bc + cd + ac + bd + cd)x^2 - (abc + bcd + abd + acd) * x + abcd | |
| 174 for (size_t aIndex = 0; aIndex < mulACount; ++aIndex) { | |
| 175 for (size_t bIndex = 0; bIndex < rootBCount; ++bIndex) { | |
| 176 for (size_t cIndex = 0; cIndex < rootCCount; ++cIndex) { | |
| 177 for (size_t dIndex = 0; dIndex < rootDCount; ++dIndex) { | |
| 178 for (size_t eIndex = 0; eIndex < rootECount; ++eIndex) { | |
| 179 testOneQuartic(aIndex, bIndex, cIndex, dIndex, eIndex); | |
| 180 } | |
| 181 } | |
| 182 } | |
| 183 } | |
| 184 } | |
| 185 } | |
| 186 | |
| 187 void QuarticRoot_Test() { | |
| 188 testOneCubic(false, 0, 5, 5, 4); | |
| 189 testOneQuartic(0, 0, 2, 4, 3); | |
| 190 quadraticTest(true); | |
| 191 quadraticTest(false); | |
| 192 cubicTest(true); | |
| 193 cubicTest(false); | |
| 194 quarticTest(); | |
| 195 } | |
| OLD | NEW |