Transcript Document
MA/CS 375
Fall 2002
Lecture 6
MA/CS375 Fall 2002
1
Exercise: Blending Two Images
• Find one partner to work with (i.e. your neighbor)
• Write a script which reads in two images.
• Use the min function to find the minimum number of rows for
both pictures (MINrows)
• Use the min function to find the minimum number of columns
for both pictures (MINcolumns)
• Crop the images using
> pic1 = pic1(1:MINrows,1:MINcols,1:3)
> pic2 = pic2(1:MINrows,1:MINcols,1:3)
• Create a third picture (pic3) which is the average of pic1 and
pic2
• Plot pic1,pic2,pic3 in the same figure using subplot
•MA/CS375
Put Fall
your
names on the title, print it out and hand it in.
2002
2
Finite Precision Effects
• Recall that we are computing in a finite
precision computing paradigm.
• i.e. all non integer numbers are
represented by a finite approximation.
• e.g. 1/2 is represented exactly, but 1/3
is not.
MA/CS375 Fall 2002
3
Non-Special IEEE Floating Point
Double Precision
• A double precision (Matlab default) number is
represented by a 64 bit binary number:
( 1) 1.m 2
s
e 1023
• Here s is the sign bit
• m is the mantissa (length to be determined)
• e is the exponent and 0 < e < 2047
MA/CS375 Fall 2002
4
A Convergent Binary
Representation of Any Number
Between 0 and 1
n
1
m mn n where mn 1 or 0
2
n 1
a similar representation in base 10:
n
1
d d n n where d n 0,1,2,3,4,5,6,7,8 or 9
10
n 1
Volunteer ?
5
MA/CS375 Fall 2002
Finite Binary Approximations of a
Real Number
n N
1
m mn n +TN where mn 1 or 0
2
n 1
1
We can easily show that TN is bounded as: TN N
2
(think of Zeno’s paradox)
MA/CS375 Fall 2002
6
Review
Definition of a Derivative
• Suppose we are given a function: f :
– where f is smooth and well defined at
every point on the real numbers.
– we can define the derivative of f at a point
x as:
df x
f x f x
lim
dx
0
– later on today we will try to use this
definition (and see what goes wrong )
MA/CS375 Fall 2002
7
Monster Functions
Due to some of these finite precision effects
there are some odd behaviors of Matlab
documented in:
“Numerical Monsters”, Christopher Essex,
Matt Davison and Christian Schulzky.
MA/CS375 Fall 2002
8
Team Tasks
We are going to investigate several examples
from:
“Numerical Monsters”, Christopher Essex,
Matt Davison and Christian Schulzky.
A.K.A “When Good Formulas Go Bad in Matlab”
MA/CS375 Fall 2002
9
Teams
• Form teams of four
• Choose one of the following monster
exercises
• Remember to put your names in the
titles
MA/CS375 Fall 2002
10
Monster #1
• Consider:
• What should its behavior be as:
f x
1 x 1
x 0
x
• Use subplots for the following 2 plots
• Plot this function at 1000 points in:
x 1,1
• Plot this function at 1000 points in:
4 ,4
• Label everything nicely, include your name in the title.
• In a text box explain what is going on, print it out and hand it
in
MA/CS375 Fall 2002
11
Monster #2
• Consider:
f x e log(1 e )
x
x
• What should its behavior be as:
x
• Plot this function at 1000 points in: x 0,50
• Explain what is going on in a text box, label
everything, print it out and hand it in.
MA/CS375 Fall 2002
12
Monster #3
• Consider:
f x 2 log 2 (1 2 )
x
x
• What should its behavior be as:
x
• Plot this function at 1000 points in:
x 0,60
• Explain what is going on in a text box and in
particular what happens at x=54, label everything,
print it out and hand it in.
MA/CS375 Fall 2002
13
Monster #4
sin x sin x
g x cos x
• What should its behavior be as:
0
• Consider:
• Plot four subplots of the function at 1000 points in:
x 0,2 * for 1e 4,1e 7,1e 12,1e 15
• Now fix x=0.5 and plot this as a function of
0.1,0.01,0.001,...,1e 20
for
• Explain what is going on, print out and hand in.
MA/CS375 Fall 2002
14
Explanations of Team Examples
MA/CS375 Fall 2002
15
Recall
Monster #1
• Consider:
f x
1 x 1
x
• What should its behavior be as: x 0
• Plot this function at 1000 points in:
• Plot this function at 1000 points in:
x 1,1
4 ,4
• Explain what is going on.
MA/CS375 Fall 2002
16
Monster #1
((large+small)-large)/small
MA/CS375 Fall 2002
17
Monster #1
((large+small)-large)/small
when we zoom in we see that the large+small operation is introducing
order eps errors which we then divide with eps to get O(1) errors !.
MA/CS375 Fall 2002
18
Monster #1
((large+small)-large)/small
Each stripe is a region
where 1+ x is a constant
(think about the gaps between
numbers in finite precision)
Then we divide by x and the
stripes look like hyperbola.
The formula looks like
(c-1)/x with a new c for each
stripe.
when we zoom in we see that the large+small operation is introducing
order eps errors which we then divide with eps to get O(1) errors !.
MA/CS375 Fall 2002
19
Recall
Monster #2
• Consider:
f x e log(1 e )
x
x
• What should its behavior be as:
x
• Plot this function at 1000 points in: x 0,50
• Explain what is going on in a text box, label
everything, print it out and hand it in.
MA/CS375 Fall 2002
20
Limit of
f x e log(1 e )
lim f x lim e log 1 e
x
x
x
x
x
x
log 1 e x
lim
x
x
e
e x
1 e x
lim
by l'Hopital's rule
x
x
e
1
lim
1
x
x 1 e
MA/CS375 Fall 2002
21
Monster #2
(finite precision effects from large*small)
f x 1 as x
As x increases past 30 we see
that f deviates from 1 !!
MA/CS375 Fall 2002
22
Monster #2 cont
(finite precision effects from large*small)
As x increases past ~=36 we
see that f drops to 0 !!
MA/CS375 Fall 2002
23
Recall
Monster #3
• Consider:
f x 2 log 2 (1 2 )
x
• What should its behavior be as:
x
x
• Plot this function at 1000 points in: x 0,60
• Explain what is going on. What happens at
x=54?
MA/CS375 Fall 2002
24
Monster 3
(finite precision large*small with binary stripes)
MA/CS375 Fall 2002
25
Monster 3
(finite precision large*small with binary stripes)
As we require more than 52 bits to represent
1+2^(-x) we see that the log term drops to 0.
MA/CS375 Fall 2002
26
Recall
Monster #4
sin x sin x
g x cos x
• What should its behavior be as:
0
• Consider:
• Plot four subplots of the function at 1000 points in:
x 0,2 * for 1e 4,1e 7,1e 12,1e 15
• Now fix x=0.5 and plot this as a function of
0.1,0.01,0.001,...,1e 20
for
• Explain what is going on, print out and hand in.
MA/CS375 Fall 2002
27
Monster 4 cont
Behavior as delta 0 :
sin x sin x
sin x cos cos( x )sin sin x
lim
lim
0
0
sin x cos 1 cos( x )sin
lim
0
cos x
or if you are feeling lazy use the definition of derivative, and remember:
d(sin(x))/dx = cos(x)
MA/CS375 Fall 2002
28
Monster 4 cont
(parameter differentiation, delta=1e-4)
OK
MA/CS375 Fall 2002
29
Monster 4 cont
(parameter differentiation, delta=1e-7)
OK
MA/CS375 Fall 2002
30
Monster 4 cont
(parameter differentiation, delta=1e-12)
Worse
MA/CS375 Fall 2002
31
Monster 4 cont
(parameter differentiation, delta=1e-15)
Bad
When we make the delta around about machine precision we see
O(1) errors !.
MA/CS375 Fall 2002
32
Monster 4 cont
(numerical instablitiy of parameter differentiation)
As delta gets smaller we see that the approximation improves, until delta ~= 1e-8 when
it gets worse and eventually the approximate derivate becomes zero.
MA/CS375 Fall 2002
33
Approximate Explanation of
Monster #4
f x f x f ' x O 2
1) Taylor’s thm:
yˆ xˆ ˆ
x O ( )
2) Round off errors
3) Round off in computation
of f and x+delta
fˆ yˆ f yˆ O ( )
f x O ( )
fˆ yˆ fˆ xˆ f y f x O
4) Put this together:
MA/CS375 Fall 2002
f ' x O 2 O
34
ˆf yˆ fˆ xˆ f ' x O 2 O
fˆ yˆ fˆ xˆ
O O
i.e. for
or equivalently
10
8
approximation error decreases as delta
decreasesize.
BUT for
MA/CS375 Fall 2002
10
8
round off dominates!.
35
Summary
• Ok – so in the far limits of the range of finite
precision really dodgy things happen.
• But recall, the formula for the derivative of sine
worked pretty well for a large range of numbers.
• Try to avoid working with two widely separated
numbers.
MA/CS375 Fall 2002
36
Summary of Lecture 6
• Today’s lecture was all about the way Matlab
(and in fact almost all languages) store
representations of real numbers
• We found out:
– because the set of possible double precision
numbers only covers a discrete subset of points
on the real line there are finite gaps between
consecutive numbers
– If A is “large” and B is very “small” the sum A+B
actually results in A+B ~= A
– and so on..
MA/CS375 Fall 2002
37