S = 144 # Number to find sqrt of
x = 10 # Initial guess
for i in range(3):
x = 0.5 * (x + S/x)
print(x)12.00000011196177114 Numerical Methods
Babylonian Square Root Method
Input: number S
Choose initial guess \(x\)
Repeat until convergence:
\(x \gets \frac{1}{2}\left(x + \frac{S}{x}\right)\)
Output: \(x\)
In Section 2.3 we saw how a variable can be repeatedly overwritten. Now that we have loops we can write this process more elegantly.
Let’s use this algorithm to approximate the square root of 144, starting with an initial guess of \(x_0 = 10\). \[x_0 = 10\] \[x_1 = \frac{1}{2} \left(x_0 + \frac{12}{x_0}\right) = 12.2 \] \[x_2 = \frac{1}{2} \left(x_1 + \frac{12}{x_1}\right) = 12.002 \] \[x_3 = \frac{1}{2} \left(x_2 + \frac{12}{x_2}\right) = 12.00002 \]
14.0.1 Fixed Number of Runs
One approach is to run the algorithm a fixed number of times, for example three times:
14.0.2 Threshold
Another approach is to run the algorithm until a threshold is reached. Roughly speaking, if you iterate the algorithm and get the same guess back, you have converged. If you run the algorithm and get nearly the same guess back, then you have nearly converged. Our threshold value will be the value \(\epsilon\) where we call our approximation “close enough” if \[ |x_{i+1} - x_{i}| < \epsilon \]
S = 144 # Number to find sqrt of
x = 10 # Initial guess
while abs(0.5 * (x + S/x) - x) > 0.001:
x = 0.5 * (x + S/x)
print(x)12.000000111961771In both cases, there is some parameter that can be modified to run the algorithm fewer or more times. This particular algorith converges toward the correct answer quikcly (quadratically), even if our initial guess is bad.
import matplotlib.pyplot as plt
S = 144 # Number to find sqrt of
x = 100 # Initial bad guess
plt.figure()
plt.plot([0,10],[12,12]) # Known square root of 144
plt.plot(0,x, 'rx')
for iteration in range(10):
x = 0.5 * (x + S/x)
plt.plot(iteration+1, x, 'rx') # Plot each updated approximation
plt.show()