Changeset 3766
 Timestamp:
 Jan 12, 2015, 7:38:52 PM (7 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/tools/lolremez/solver.cpp
r3765 r3766 180 180 for (int i = 0; i < m_order + 1; i++) 181 181 { 182 struct { real value, error; } left, right, mid;183 184 left.value = m_control[i];185 left.error = EvalEstimate(left.value)  EvalFunc(left.value);186 right.value = m_control[i + 1];187 right.error = EvalEstimate(right.value)  EvalFunc(right.value);182 struct { real value, error; } a, b, c; 183 184 a.value = m_control[i]; 185 a.error = EvalEstimate(a.value)  EvalFunc(a.value); 186 b.value = m_control[i + 1]; 187 b.error = EvalEstimate(b.value)  EvalFunc(b.value); 188 188 189 189 static real limit = ldexp((real)1, 500); 190 190 static real zero = (real)0; 191 while (fabs( left.value  right.value) > limit)191 while (fabs(a.value  b.value) > limit) 192 192 { 193 mid.value = (left.value + right.value) / 2; 194 mid.error = EvalEstimate(mid.value)  EvalFunc(mid.value); 195 196 if ((left.error <= zero && mid.error <= zero) 197  (left.error >= zero && mid.error >= zero)) 198 left = mid; 193 /* Interpolate linearly instead of taking the midpoint, this 194 * leads to far better convergence (6:1 speedups). */ 195 real t = abs(b.error) / (abs(a.error) + abs(b.error)); 196 real newc = b.value + t * (a.value  b.value); 197 198 /* If the third point didn't change since last iteration, 199 * we may be at an inflection point. Use the midpoint to get 200 * out of this situation. */ 201 c.value = newc == c.value ? (a.value + b.value) / 2 : newc; 202 c.error = EvalEstimate(c.value)  EvalFunc(c.value); 203 204 if (c.error == zero) 205 break; 206 207 if ((a.error < zero && c.error < zero) 208  (a.error > zero && c.error > zero)) 209 a = c; 199 210 else 200 right = mid;211 b = c; 201 212 } 202 213 203 m_zeroes[i] = mid.value;214 m_zeroes[i] = c.value; 204 215 } 205 216
Note: See TracChangeset
for help on using the changeset viewer.