2022-01-14 Function errors
Contents
2022-01-14 Function errors¶
Summing series
Conditioning and well posedness
Relative and absolute errors
Activity
using Plots
default(linewidth=4)
What is floating point arithmetic?¶
fuzzy arithmetic
exact arithmetic, correctly rounded
the primary focus of numerical analysis
#(0.1 + 0.2) - 0.3
2^(-53)
1.1102230246251565e-16
plot(x -> (1 + x) - 1, xlims=(-1e-15, 1e-15))
plot!(x -> x)
Machine epsilon¶
We approximate real numbers with floating point arithmetic, which can only represent discrete values. In particular, there exists a largest number, which we call
The notation
eps = 1
while 1 - eps != 1
eps = eps / 2
end
eps
5.551115123125783e-17
eps = 1.f0
while 1 + eps != 1
eps = eps / 2
end
eps
5.9604645f-8
Beating exp
¶
Suppose we want to compute
f1(x) = exp(x) - 1
y1 = f1(1e-8)
9.99999993922529e-9
f2(x) = x + x^2/2 + x^3/6
y2 = f2(1e-8)
1.000000005e-8
Which answer is more accurate?
@show (y1 - y2) # Absolute difference
@show (y1 - y2) / y2; # Relative difference
y1 - y2 = -1.1077470910720506e-16
(y1 - y2) / y2 = -1.1077470855333152e-8
Conditioning¶
We say that a mathematical function
The function
Conditioning¶
A function may also represent something more complicated, implementable on a computer or by physical experiment.
Find the positive root of the polynomial
Find the eigenvectors of the matrix
Find much the bridge flexes when the truck of mass
drives over it.Find how long it takes to clean up when the toddler knocks the glass off the counter, as a function of the chair location
.Find the length of the rubber band when it finally snaps, as a function of temperature
during manufacturing.Find the time at which the slope avalanches as a function of the wind speed
during the storm.Find the probability that the slope avalanches in the next 48 hours as a function of the wind speed
during the storm.Find the probability that the hurricane makes landfall as a function of the observations
.
Specification¶
Some of these problems are fully-specified
Others involve sophisticated models and ongoing community research problems.
In some cases, the models that are computable may incur greater uncertainty than the underlying system. In such cases, an analog experiment might produce smoother variation of the output as the problem data
are varied.In others, the model might be better behaved than what it seeks to model.
Some of these problems may be ill-posed
Well-posedness¶
A problem is said to be well-posed if
a solution exists,
the solution is unique, and
the solution depends continuously on the problem specification.
Mathematically, continuous variation in part 3 can be arbitrarily fast, but there may be measurement error, manufacturing tolerances, or incomplete specification in real-world problems, such that we need to quantify part 3. This is the role of conditioning.
Computing ¶
function myexp(x)
sum = 1
for k in 1:100
sum += x^k / factorial(k)
end
return sum
end
myexp(1) - exp(1)
OverflowError: 21 is too large to look up in the table; consider using `factorial(big(21))` instead
Stacktrace:
[1] factorial_lookup
@ ./combinatorics.jl:19 [inlined]
[2] factorial
@ ./combinatorics.jl:27 [inlined]
[3] myexp(x::Int64)
@ Main ./In[9]:4
[4] top-level scope
@ In[9]:8
[5] eval
@ ./boot.jl:373 [inlined]
[6] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)
@ Base ./loading.jl:1196
function myexp(x)
sum = 0
term = 1
n = 1
while sum + term != sum
sum += term
term *= x / n
n += 1
end
sum
end
myexp(1) - exp(1)
4.440892098500626e-16
What’s happening?¶
We’re computing
for values of near zero.This function is well approximated by
.Values of
near 1 cannot represent every value.After rounding, the error in our computed output
is of order .
Absolute Error¶
Relative Error¶
Suppose I want to compute ¶
plot([x -> myexp(x) - 1 , x -> x],
xlims=(-1e-15, 1e-15))
What now?¶
We’re capable of representing outputs with 16 digits of accuracy
Yet our algorithm
myexp(x) - 1
can’t find themWe can’t recover without modifying our code
Modify the code¶
function myexpm1(x)
sum = 0
term = x
n = 2
while sum + term != sum
sum += term
term *= x / n
n += 1
end
sum
end
myexpm1 (generic function with 1 method)
plot(myexpm1, xlims=(-1e-15, 1e-15))
Plot relative error¶
function relerror(x, f, f_ref)
fx = f(x)
fx_ref = f_ref(x)
max(abs(fx - fx_ref) / abs(fx_ref), 1e-17)
end
badexpm1(t) = exp(t) - 1
plot(x -> relerror(x, badexpm1, expm1), yscale=:log10, xrange=(-2, 2))
Activity: 2022-01-12-taylor-series¶
Use Julia, Jupyter, Git
Look at how fast series converge when taking only finitely many terms
Explore instability, as is occuring for large negative
x
above, but not for the standard libraryexpm1