Numerical approximation using Halley's method. I'm working on an R exercise, but I'm having difficulty grasping the math behind the exercise in order to implement it properly.
The exercise requires me to approximate a function using the first order recurrence relation 
$$X_{n+1} = X_n - \frac{ f(X_n)} {f'(X_n) - \frac{f(X_n)f''(X_n)}{2f'(X_n)}}$$
aka Halley's/Bailey's method. Specifically I am asked to approximate 59^(1/7) to 9 decimal places. (I'm assuming 9 decimal places also means within a tolerance of 1e-9, is this correct?) What I need help with is understanding how to get each successive approximation by hand until I reach 9 accurate decimal places using this method.
 A: The first step is to realize, that computing $59^\frac{1}{7}$ can be done by computing the root of $f(x) = x^7 - 59\;$ in the inteval $(1,2)$ because $2^7 =128 > 59$.
Let's write the Halley/Bailey formula in the form
$$x_{n+1} = x_n - d(x_n)$$
$$d(x) = \frac{ f(x)} {f'(x) - \frac{f(x)f''(x)}{2f'(x)}}$$
From this you can easily get the actual changes for the iteration process and stop if $|d_k| = |d(x_k)| < 10^{-9}.$ 
Using the definition of $f(x)$ you can simplify $d(x)$ to get
$$d(x) = \frac{x(x^7-59)}{4x^7+177}\cdot$$
Now start the iteration with the initial guess $x_0=2.\;$ Here the resulting list of $d_n, x_{n+1}$
d0 = 0.2002902757619739      x1 = 1.799709724238026
d1 = 0.009190072516764635    x2 = 1.790519651721261
d2 = 0.9611796341215367E-6   x3 = 1.790518690541627

d3 = 0.1300618903541134E-15  x4 = 1.790518690541627

The last step is actually not necessary because the convergence order is 3 in this case (you can expect $d_3 \approx 10^{-18},\;$ it is a bit larger because I compute with 16 digits, and so $d_3$ is about $10^{-16}$).
