This behavior caught my eye the other day:
PROGRAM Sine
USE, INTRINSIC :: ISO_Fortran_env, dp=>REAL64 !modern DOUBLE PRECISION
IMPLICIT NONE
REAL(dp), PARAMETER :: PI = 3.1415926535897931 !PRINT *, ACOS(-1.0D+00)
PRINT *, SIN(PI / 4.0_dp)
PRINT *, COS(PI / 4.0_dp)
PRINT *, SQRT(2.0_dp) / 2.0_dp
END PROGRAM Sine
Unformatted output gives 0.70710679664085752, 0.70710676573223719, and 0.70710678118654757, respectively. I don’t understand why this doesn’t agree to double precision: SQRT(2.0D+00) / 2.0D+00 is exact, while SIN(PI / 4.0D+00) and COS(PI / 4.0D+00) should come from the Taylor series for those functions. I can get slightly better agreement out of a pocket calculator, so clearly I'm doing something wrong. Replacing the above intrinsic functions with DSIN, DCOS, and DSQRT produces the same result--though I wouldn't expect any different. Does anyone understand why this is happening?
