From: ph@pixar.UUCP (Paul Heckbert) Date: Mon Jul 20 12:42:47 PDT 1987 The "best" computational formula I've seen for ray-sphere intersection testing is what you called the "brute force" method, and goes something like this: ------------------------------------------ #define EPSILON 1e-7 /* tolerance to roundoff errors */ /* * SphereIntersect: find intersection of ray X=P+tD (where |D|=1) * with sphere |C-X|^2 = r^2 by solving the quadratic equation: * |V-tD|^2 = t^2 - 2t(V.D) + (V.V-r^2) = 0 * where V=C-P. * The roots are: * t = V.D +- sqrt((V.D)^2 - V.V + r^2) */ SphereIntersect(P, D, C, r2, tp) double P[3], D[3]; /* ray origin point and direction */ double C[3], r2; /* sphere center and radius_squared */ double t[2]; /* t at intersection points (if any) */ { double V[3], b, disc; s m a 0 0 3 V = vsub(C, P); /* vector subtract */ 0 3 2 b = vdot(V, D); /* dot product */ 0 4 4 disc = b*b - vdot(V, V) + r2; if (disc < 0.) return 0; /* ray misses sphere */ 1 0 0 disc = sqrt(disc); 0 0 1 t[0] = b-disc; /* first root */ 0 0 1 t[1] = b+disc; /* second root */ return 2; /* return #roots */ } ------------------------------------------ The numbers under "s m a" are the #sqrts, #multiplies, and #add/subtracts, respectively. With this scheme, it takes 7 multiplies and 9 adds to determine if the ray hits the sphere, and 1 sqrt, 2 additions more to find t at both intersection points. I'd wager that this is "optimal" given this set of inputs, outputs, and operators. [ General tip to fellow optimization freaks: I've found that I can often knock off a few superfluous multiplies and divides if I use the quadratic formula: roots of x^2-2b*x+c=0 are x=b+-sqrt(b*b-c) instead of the tired old formula: roots of a*x^2+b*x+c=0 are x=(-b+-sqrt(b*b-4*a*c))/(2*a) ]