+ (multiple-value-bind (q r)
+ (round (rational x) rational-pi-half)
+ (values q (coerce r 'double-float))))
+ (expected-val (op x)
+ "Calculate (OP X) precisely by shifting the argument by
+ an integral multiple of (/ PI 2) into the range from
+ (- (/ PI 4)) to (/ PI 4) and applying the phase-shift
+ formulas for the trigonometric functions. PI here is
+ precise enough that the result is exact to double-float
+ precision."
+ (labels ((precise-val (op q r)
+ (ecase op
+ (sin (let ((x (if (zerop (mod q 2))
+ (sin r)
+ (cos r))))
+ (if (<= (mod q 4) 1)
+ x
+ (- x))))
+ (cos (precise-val 'sin (1+ q) r))
+ (tan (if (zerop (mod q 2))
+ (tan r)
+ (/ (- (tan r))))))))
+ (multiple-value-bind (q r)
+ (round-pi-half x)
+ (precise-val op q r))))