1.5e7 # a floating-point value 1.5 × 10⁷

x = 1/49 # division of two integers produces a floating-point value

x * 49

1 - x * 49

2.0^-52, eps(), eps(1.0), eps(Float64) # these are all the same thing

# 1.5 is exactly represented in binary floating point:
big(1.5), 1.5 == 3//2

# 0.1 is *not* exactly represented
big(0.1), 0.1 == 1//10

function my_cumsum(x)
    y = similar(x) # allocate an array of the same type and size as x
    y[1] = x[1]
    for i = 2:length(x)
        y[i] = y[i-1] + x[i]
    end
    return y
end

eps(Float32), eps(Float64)

x = rand(Float32, 10^7) # 10^7 single-precision values uniform in [0,1)
@time y = my_cumsum(x)
yexact = my_cumsum(float64(x)) # same thing in double precision
err = abs(y - yexact) ./ abs(yexact) # relative error in y

using PyPlot
n = 1:10:length(err) # downsample by 10 for plotting speed
loglog(n, err[n])
ylabel("relative error")
xlabel("# summands")
# plot a √n line for comparison
loglog([1,length(err)], sqrt([1,length(err)]) * 1e-7, "k--")
gca()[:annotate]("~ √n", xy=(1e3,1e-5), xytext=(1e3,1e-5))
title("naive cumsum implementation")

VERSION < v"0.3.6" && warn("cumsum is less accurate for Julia < 0.3.6")
@time y2 = cumsum(x)
err2 = abs(y2 - yexact) ./ abs(yexact)
loglog(n, err2[n])
ylabel("relative error")
xlabel("# summands")
title("built-in cumsum function")
loglog(n, sqrt(log(n)) * 1e-7, "k--")

1e300 # okay: 10³⁰⁰ is in the representable range

(1e300)^2 # overflows

realmax(Float64), realmax(Float32)

1e-300 # okay

(1e-300)^2 # underflows to +0

realmin(Float64), realmin(Float32)

-1e-300 * 1e-300 # underflows to -0

+0.0 == -0.0

1 / +0.0, 1 / -0.0

signbit(+0.0), signbit(-0.0)

1 / (1 / -Inf)

0 * Inf, Inf / Inf, 0 / 0

sqrt(-1.0)

sqrt(-1.0 + 0im)