https://gcc.gnu.org/bugzilla/show_bug.cgi?id=93871
--- Comment #33 from Thomas Henlich <thenlich at gcc dot gnu.org> --- Going back one step, I wonder if it would be good enough to perform a correctly rounded conversion from degrees to radians, and just use sin() as is. ... f = sgn * sin(deg2rad(arg)) end function fcn ! Compute correctly rounded value of x * pi / 180 for x = 0 ... 90 function deg2rad(x) result(y) real(sp) y real(sp), intent(in) :: x real(sp), parameter :: a = atan(1._4) / 45 real(sp), parameter :: b = atan(1._16) / 45 - a y = a * x + b * x end function deg2rad (The exact values of a and b still need some work)