Introduction to Numerical Problem Solving TX00BY09-3007

Assignment: 02 Approximations and round-off errors

Description: Solve problems 3, 4, and 5 from exercises 02.

Author: Joonas Forsberg

This document has been prepared for the second assignment for the course "Introduction to Numerical Problem Solving". It contains code examples and explanations for the written code. In addition, I've included reasoning behind the behavior where acceptable.

Problem 3

Modify your code to calculate the smallest number xmin, that differs from zero.

Problem 3 uses code from the previous exercise as a reference and it's included as a function exercise2 in the code example below.


In [65]:
def exercise2():
    
    eps = 1.0
    
    while eps + 1.0 > 1.0:
        eps = eps/2.0
    
    eps = 2.0 * eps
    print("Final value for eps is {}".format(eps))
    
def exercise3(start):
    
    x = start
    
    while x != 0.0:
        
        if x / 2 == 0.0:
            break
    
        x = x / 2
        
    print("Smallest possible number to differ from 0 is {}. Starting value set to {}".format(x, start))

exercise3(-0.5)
exercise3(0.1)


Smallest possible number to differ from 0 is -5e-324. Starting value set to -0.5
Smallest possible number to differ from 0 is 5e-324. Starting value set to 0.1

In this particular case, the smallest number to differ from 0.0 would be either 5e-324 or -5e-324. I've decided to break the while-loop in case x/2.0 is 0.0. If the value becomes 0.0, I would not be able to multiply it to receive the previous value the same way I did in the previous exericse.

There are two function calls with negative and positive starting values to demonstrate the difference being the same, despite of the starting value.

Problem 4

How does the (epsilon) change, when you increase or decrease the study point (e.g. 1.0 or 0.0 in the previous exercise)? What is the relationship of the study point and epsilon? Use values that differ several decades


In [140]:
%matplotlib notebook
import matplotlib.pyplot as plt
import numpy as np
import pylab
from matplotlib.pyplot import *
from numpy import *

# Create an array for the results to be saved
data = []

def exercise4(increment, count):
    y = 1.0
    
    i = 1
    while i <= count:
        eps = y
        while eps + 1.0 > 1.0:
            eps = eps/2.0
        
        eps = 2.0 * eps
        #print("Round {}, y = {}: Value for epsilon is {}".format(i, y, eps))
        i += 1
        y += increment
        
        data.append(eps)

exercise4(0.1,75)

# Create chart
fig, ax = plt.subplots()

# Set labels
ax.set_title('Value of Epsilon')
ax.set_ylabel('Epsilon')
ax.set_xlabel('Study Point')

# Draw chart
plt.plot(data)
plt.show()


I have decided to go with graphical represantion for this problem to give an easier way to interpret the data provided by the function exercise4. In the function, I am incrementing the base value (1.0) by the amount given in the function call, by x times, which is also given in the function call.

Data is then saved into an array which can be used to represent the data in graphical form to point of any relations between the study points and values for eps.

From the graph we can come up with couple of conclusions:

  • Value for eps is roughly between 2.15e-16 and 1.11e-16.
  • The value for eps decreases gradually until it reaches it's lowest value.
  • After reaching it's lowest value, eps starts to decrease from it's highest point
  • As the amount of study point increases, the amount of study points required for each cycle is roughly doubled.

If we take a look at the Wikipedia arcticle "Machine Epsilon" (https://en.wikipedia.org/wiki/Machine_epsilon), we can find some relation to the numbers we've seen in our graph. At the table, "Values for standard hardware floating point arithmetics", we can see the machine epsilon value for double precision being either 1.11e-16 or 2.22e-16, depending on which Professor you ask.

So based on this, I would say the code example above gives us an approximation that epsilon is somewhere in between 1.11e-16 and 2.22e-16, but is depended on the actual value which is being used to calculate the value. The search algorithm used in the process in "Linear Search"

Problem 5

The derivative of $f(x) = 1/(1 - 3x^2)$ is given by: $$f'(x) = \frac{6x}{(1-3x^2)^2}$$

Do you expect to have difficulties evaluating this function at x = 0.577? How does the value of the derivative change when you increase the significant figures in your calculations. Try using from 3-digit to 16- digit arithmetics. Start with the following code:


In [135]:
from decimal import *

def exercise5(x, prec):
    getcontext().prec = prec
    
    x = Decimal(x)
    
    x = 6*x / ((1-3*x**2)**2)
    
    return x

i = 3
while i <= 16:
    print("Value with {} decimals is {}".format(i, exercise5(0.577, i)))
    i += 1


Value with 3 decimals is 3.46E+6
Value with 4 decimals is 2.049E+6
Value with 5 decimals is 2.3646E+6
Value with 6 decimals is 2.35291E+6
Value with 7 decimals is 2352911
Value with 8 decimals is 2352910.8
Value with 9 decimals is 2352910.79
Value with 10 decimals is 2352910.793
Value with 11 decimals is 2352910.7926
Value with 12 decimals is 2352910.79260
Value with 13 decimals is 2352910.792602
Value with 14 decimals is 2352910.7926020
Value with 15 decimals is 2352910.79260199
Value with 16 decimals is 2352910.792601992

To answer the question, "Do you expect to have difficulties evaluating this function at x = 0.577?", my first answer would be no. Why would I have problems with this? Once we take a look at the results, it's clear that I was wrong. By using 3 digit accuracy, the answer is almost 1,5 times larger than the actual result is at 16 digits.

The integer part keeps changing it's value until we are using the accuracy of 8 decimals. This of course depending on application might not be suitable and we would still have to increase the accuracy to get the desired accuracy we are looking for. For example at 16 digit accuracy, the first 5 decimals will remain unchanged as the accuracy increases.