Hi!

While trying out some stuff in Julia, I came across the fact that cos(pi/2) does not give exact 0, but some type dependent output (6e-17 for Float64, -4e-8 for Float32, 5e-4 for Float16). Is this intended? And what is the best way to avoid this?

Not intended, but hard to avoid because of floating point error. See

You can use

```
julia> cospi(0.5)
0.0
```

If what you are wondering is why Julia does not know the cosine involving the `Irrational`

object `π`

exactly, this is because it is converted to a float when performing the operation `π/2`

. Base Julia does not implement `AffineFunctionOfIrrational`

or anything like that. Someone might have done so in a package.

Something similar to `UniformScaling`

for irrationals could maybe help fix that issue.

I see, thanks so much for your responses!

I’m stuck with one more question though, but that might be because I’m starting out with Julia:

I’ve had a look in the source code at `/opt/julia/share/julia/base/special/trig.jl`

```
function cos(x::T) where T<:Union{Float32, Float64}
absx = abs(x)
if absx < T(pi)/4
if absx < sqrt(eps(T)/T(2.0))
return T(1.0)
end
return cos_kernel(x)
elseif isnan(x)
return T(NaN)
elseif isinf(x)
cos_domain_error(x)
else
n, y = rem_pio2_kernel(x)
n = n&3
if n == 0
return cos_kernel(y)
elseif n == 1
return -sin_kernel(y)
elseif n == 2
return -cos_kernel(y)
else
return sin_kernel(y)
end
end
end
```

that is apparently using approximations for cos(x), -sin(x), -cos(x) or sin(x) depending on the interval x is in.

(I must admit I haven’t had a closer look at the `rem_pio2_kernel`

function, but I admit that it returns the div and rem value from the division by π/2)

I defined my own function based on the same idea of evaluating the trigonometric functions only in the interval [-π/4, π/4] around either their extremum (cosine) or root (sine) as suggested by `cos_kernel`

```
function mycos(x)
if x < 0
return mycos(-x)
else
j = div(x+pi/4,pi/2)%4
newx = x%2pi - j*pi/2
if j == 0 return cos(newx)
elseif j == 1 return -sin(newx)
elseif j == 2 return -cos(newx)
elseif j == 3 return sin(newx)
end
#if that should (in whatever case) not work, return normal value
return cos(x)
end
end
```

This does return the values at multiple of π/2 correctly (not talking about performance here). I don’t see too may differences to the function in `trig.jl`

and wonder why that one is not capable of doing that…

Probably because it’s not in fact correct (e.g. the cosine of the closest Float64 representation of `pi/2`

is not zero) and you get worse results everywhere. For example

```
julia> x = 1.571
1.571
julia> y1 = mycos(x)
-0.00020367320369523912
julia> y2 = cos(x)
-0.00020367320369517786
julia> y3 = cos(BigFloat(x))
-2.036732036951778709106180717800557926547406375863370495134343067126837115362263e-04
julia> abs(y1 - y3)
6.124442755166129148251053819710905768303105297194328731628846377367955634449768e-17
julia> abs(y2 - y3)
1.299519376970903891054505204264901252558765305671268371153622632044365550231594e-20
```

If you are specifically concerned about the cosine of multiples of pi, use the `cospi`

function. It’s there for exactly that purpose.

Alright! Thanks for your answer.