What is this function related with continued fractions?

Question

Playing with continued fractions, I came with the idea of iterating the limit of the simplest one:

$$1 + \cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1+\cdots}}}}}\ = \Phi$$

And then I thoght about iterating the result: $$\Phi + \cfrac{1}{\Phi+\cfrac{1}{\Phi+\cfrac{1}{\Phi+\cfrac{1}{\Phi+\cfrac{1}{\Phi+\cdots}}}}}\ = ?$$

And after that keep iterating again and again. Essentially, what I am doing is solving this:

$$x = a + \frac{1}{x}$$

for different values of $a$. Whose solution is:

$$ x = \frac{a \pm \sqrt{a^2 + 4}}{2}$$

And I'm just staying with the positive solution. At the end, I'm just exploring this sequence defined recursively:

$$ f(1) = \Phi \\ f(n+1) = \frac{f(n) + \sqrt{f(n)^2 + 4}}{2} $$

When you see the plot of the terms of the sequence, you get this picture:

I would like to know which is the function that generates that curve. It really seems to be a continuos function (notice that the picture is a set of points so close, that may seem to be continuous, but it's not).

So... What's the function underlying this curve and how can we find it?

Many thanks in advance!

Answer 1

This question is similar to MSE question 1072256 "$a_{n+1} = \log(1+a_n), a_1>0.$ ..." and to MSE question 3215 "Convergence of $\sqrt{n}x_n$ where $x_{n+1}=\sin(x_n)$". Because of the similarities, it is natural to assume that there is a power series $\, f(n) \approx f^*(n) := \sqrt{2n} \, (1 + c_1 x - c_2 x^2 + c_3 x^3 - c_4 x^4 + O(x^5)) \,$ where $\, x := 1/n \,$ and $\, c_k \,$ is a polynomial of degree $k$ in $\, y := \log(x). \,$ Using the defining recursion for $\, f(n) \,$ we can solve for the coefficients $\, c_k \,$ to get $$ c_1 = c + y/8, \quad c_2 = (1+4c)^2/32 + (1+4c)y/32 + y^2/128,$$ $$ c_3 \!=\! \frac{11 \!+\! 120 c \!+\! 384 c^2 \!+\! 384 c^3}{768} \!+\! \frac{5 \!+\! 32 c \!+\! 48 c^2}{256} y \!+\! \frac{1 \!+\! 3 c}{128} y^2 \!+\! \frac{y^3}{1024},$$ where $\, c\,$ is a constant depending on $\, f(0). \,$ If we use the initial values $\, f(0) = 0, f(1) = 1, \, f(2) = \phi, \dots, \,$ then $\, c \approx -0.291131527. \,$ Note that since $\, \sqrt{1/2}-1 \approx -0.292 \,$ the formula is already close for $\, n=1.$

To judge the accuracy, we get $\, f(1) = 1, \,$ $ f^*(1) \approx 1.00269, \,$ $f(2) \approx 1.61803, \,$ $ f^*(2) \approx 1.61826, \,$ $f(3) \approx 2.09529, \,$ and $ f^*(3) \approx 2.09535. \,$

Answer 2

I think you can set up a differential equation and solve it.

\begin{eqnarray} 2f'(n) &= -f(n)+\sqrt{f(n)^2+4}\\ \text{d}n&=\frac{2\, \text{d}\! f}{-f+\sqrt{f^2+4}} \\ n &= \frac{1}{4} (f (f + \sqrt{4 + f^2} + 4 \sinh^{-1}\left(\frac{f}{2}\right) \end{eqnarray} Solving for $f$ should result in the function you wanted. Wolfram plot

Answer 3

One way to solve this is to use the theory of continuous iterations, which is often best tackled using power series around fixed points. Let $g(z) = \frac{z + \sqrt{z^2 + 4}}2$. Then notice $$ g^{-1}(w) = \frac{w^2 - 1}w $$ This has no fixed points except at infinity, but we can actually make use of that if we conjugate by $w\to \frac1w$, since $h(w) = \frac1{g^{-1}(1/w)}$ will have a fixed point at $w=0$. We can see the series for $h(w)$ is given by $$ h(w) = \frac{w}{1-w^2} = w + w^3 + w^5 + ... $$ Formally composing this series with itself is fairly straightforward, for example: \begin{eqnarray} h(h(w)) &=& (w + w^3 + w^5 +...) + (w + w^3 + w^5 +...)^3 + (w + w^3 + w^5 +...)^5 + ... \\ &=& w + 2w^3 + 2w^5 + ... \end{eqnarray} In general, we will have $$ h^n(x) = w + p_3(n)w^3 + p_5(n)w^5 + ... $$ where $p_i(n)$ satisfy a recurrence given by \begin{equation} w + p_3(n+1) w^3 + p_5(n+1)w^5 + ... = h^{n+1}(x) = h(h^n(x)) \\ = (w + p_3(n)w^3 + p_5(n)w^5 + ...) + (w + p_3(n)w^3 + p_5(n)w^5 + ...)^3 + (w + p_3(n)w^3 + p_5(n)w^5 + ...)^5 + ... \\ = w + (p_3(n) + 1) w^3 + (p_5(n) + 3p_3(n) + 1) w^5 + ... \end{equation} We get recursions like $p_3(n+1) = p_3(n) + 1$. Clearly, this is solved by $p_3(n) = n$. The $w^5$ coefficient is $p_5(n+1) = p_5(n) + 3n + 1$, which is solved by $p_5(n) = n + \frac32 n(n+1)$. In general, you can actually find a unique polynomial for each $p_i(n)$ by using Faulhaber's formula since each one is defined by a recursion $p_i(n+1) = p_i(n) + q_i(n)$ for some appropriate polynomial $q_i(n)$. Notice this recursion will still be satisfied for noninteger $n$, so we have a function $H(t,w) = w + p_3(t)w^3 + p_5(t)w^5 +...$ satisfying: $$ H(t+1,w) = h(H(t,w)) $$ OK so what does this have to do with the original question? Well, we have \begin{eqnarray} H(t+1,w) &=& \frac1{g^{-1}(1/H(t,w))} \\ \frac1{H(t+1,w)} &=& g^{-1}\left(\frac1{H(t,w)}\right) \\ g\left(\frac1{H(t+1,w)}\right) &=& \frac1{H(t,w)} \end{eqnarray} substituting $-n = t+1$, and recalling the definition of $g$, we find: \begin{eqnarray} \frac1{H(-(n+1),w)} &=& g\left(\frac1{H(-n,w)}\right) \\ \frac1{H(-(n+1),w)} &=& \frac{\frac1{H(-n,w)} + \sqrt{\frac1{H(-n,w)^2} + 4}}2 \end{eqnarray} Then, as you can see, we will have $f(n) = \frac1{H(-n,w)}$, where $w$ is the starting value. I do not expect this function to have a nice closed form.

This formula for $f$ is only helpful on a small scale, but it does demonstrate we can draw a smooth (hence continuous) function through those points. If you want to know about the asymptotic behavior of $f(n)$ as $n\to\infty$, somos has already covered that in his answer.

Answer 4

Using what has been given in comments.

If we use $$2n=f^2+\log(f^2-1)\implies 2n-1=(f^2-1)+\log(f^2 -1)$$ let $t=(f^2-1)$ to make $$2n-1=t+\log(t)\implies t\, e^t=e^{2n-1}\implies t=W\left(e^{2 n-1}\right)$$ where appears Lambert function making the final solution to be
$$f=\sqrt{1+W\left(e^{2 n-1}\right)}\tag1$$

If we use $$n=\frac{1}{4} f \left(f+\sqrt{f^2+4}\right)+\sinh ^{-1}\left(\frac{f}{2}\right)$$ and expand the rhs as a series for infinitely large values of $f$, we should get $$n=\frac{f^2}{2}+\frac{1}{2}+\log \left({f}\right)+O\left(\frac{1}{f^2}\right)=\frac{f^2}{2}+\frac{1}{2}+\frac{1}{2}\log \left({f^2}\right)+O\left(\frac{1}{f^2}\right)$$ and, ignoring the higher order terms, the same process would lead to $$f=\sqrt{W\left(e^{2 n-1}\right)}\tag2$$

Now, for large values of $n$, you could use the asymptotic series given in the linked Wikipedia page $$W(x)=L_1-L_2+\frac{L_2}{L_1}+\frac{L_2(L_2-2)}{2L_1^2}+\frac{L_2(6-9L_2+2L_2^2)}{6L_1^3}+\cdots$$ where $L_1=\log(x)$ and $L_2=\log(L_1)$ and get, as a very first approximation $$f\approx \sqrt{(2n-1)-\log(2n-1)}$$

To allow you to compare with your numbers, I generated a table of numerical values $$\left( \begin{array}{cc} n & (1) \\ 500 & 31.5135 \\ 1000 & 44.6363 \\ 1500 & 54.6991 \\ 2000 & 63.1800 \\ 2500 & 70.6504 \\ 3000 & 77.4035 \\ 3500 & 83.6131 \\ 4000 & 89.3925 \\ 4500 & 94.8203 \\ 5000 & 99.9539 \end{array} \right)$$