An Implicit Solution to a Nonlinear Difference Equation

By 苏剑林 | April 09, 2016

Problem Source

A thread on the Mathematics Research Forum, which I frequent, discussed the asymptotic solution for the following nonlinear difference equation:

The original post is here. I benefited a lot from that thread and learned many new techniques. The main idea is through cubing both sides and then setting $x_n=a_n^3$, transforming it into an equivalent recurrence problem:

Then, by using clever techniques, the asymptotic expansion can be obtained:

I will not mention the specific process here; readers can go to the aforementioned post to learn on their own. However, although this form of the solution is exquisite, there are some aspects that I am not entirely satisfied with:

1. The solution is an asymptotic series, which means the radius of convergence is actually 0;
2. It is a solution in the form of $n^{-k}$, which makes it difficult to calculate for small $n$, making high-precision calculation relatively difficult;
3. Of course, the original purpose of the problem was asymptotic calculation, but there seems to be no need for the asymptotic analysis to expand into so many terms;
4. It contains a limit value $a$ that is inherently difficult to calculate;
5. The solving process seems slightly unintuitive.

Of course, some of the above shortcomings are a bit like nitpicking. However, it is precisely these shortcomings that prompted me to look for a better form of the solution, which eventually led to this article.

Implicit Solution

The recurrence

The first-level asymptotic solution is obtained by omitting the subsequent $\frac{3}{x_n}+\frac{1}{x_n^2}$, yielding $x_{n+1}=x_n+3$, so $x_n=3n-2$. There are many subsequent ideas, such as iteration or undetermined coefficients, which can all find the asymptotic solution, but I tried another path: finding an implicit solution. First, we consider:

Expanding the last term at $x_n$ to the first order:

Taking $1/x_n$ as the order and omitting second-order terms, we get

If we let $f(x_n)=-\ln x_n$, then the square bracket becomes 0, so we obtain

The accuracy of this equality is $\mathcal{O}(x_n^{-2})$. To balance simplicity and accuracy, it is better to start from $x_2=8$, so

This is an implicit (approximate) solution with an accuracy of $\mathcal{O}(x_n^{-1})$. It maintains almost the same accuracy for nearly all $n$. The key here is to introduce a new term so that it becomes a linear recurrence (an arithmetic progression). This process can continue. Based on the previous work, by introducing a new function and continuing to expand it to the second order, we calculate:

It has $\mathcal{O}(x_n^{-3})$ accuracy, from which an implicit solution with $\mathcal{O}(x_n^{-2})$ accuracy can be solved.

Why Implicit Solution?

Why look for implicit solutions? There are roughly the following advantages. Under normal circumstances, traditional asymptotic solutions are obtained through Taylor series expansions. Taylor series expansions themselves have limitations (requiring differentiability and finite derivatives), so asymptotic solutions tend to appear; even if they are not asymptotic, the radius of convergence may be small. Implicit function solutions, generally speaking, are more stable and have better convergence. Although they may still be asymptotic, the divergence rate will be much lower. For example, $x=\sqrt{1+t}$, expanded as a power series, has a radius of convergence of only 1, but it can be expressed in the implicit function form $x^2-1=t$, maintaining both accuracy and simplicity.

In short, if an explicit expansion can do it, an implicit function can basically do it too; and for things an explicit expansion cannot do, an implicit function might still succeed. Therefore, implicit functions clearly perform better. In addition, for the recurrence in this article, the implicit function solution is more concise and does not involve the difficult-to-calculate $a$ everywhere. Readers may feel that solving a nonlinear equation every time $x_n$ is needed is complicated. In fact, if one wishes, the explicit solution can be obtained from the implicit solution directly through the inverse function, but this returns to asymptotic series and is not very necessary.

Recurrence Format Conducive to Programming

The two terms calculated earlier were just an introductory demonstration. Since the order was not high, the calculation was not difficult. However, to calculate further, especially for programming purposes, it is necessary to construct a recursive format that is easy to understand, much like the recursive calculation of the perturbation method. Suppose the introduced $f(x_n)$ is exact, then clearly, for an exact $f(x)$, it must satisfy:

At this time, the original recurrence becomes

Solving becomes much easier. After analysis, a parameter $q$ can be artificially introduced, and $f(x)$ can be treated as a binary function $f(x,q)$ of $x$ and $q$. After solving for the series solution of $f(x,q)$, we let $q=1$ to get $f(x)=f(x,1)$. The introduced format is as follows:

The idea behind this introduction is: take $x^{-1}$ as the infinitesimal order and let $f^{(n)}(x)$ also have the order of $x^{-n}$. Therefore, inside $f()$, $x^{-n} \to q^{n+1}x^{-n}$, and outside $f()$, there is $x^{-n} \to q^n x^{-n}$. At this time, software like Mathematica can be used to expand it quickly:

Series[f[x + 3*q + 3*q^2/x + q^3/x^2, q] - f[x, q] + 3*q/x + q^2/x^2, {q, 0, 5}]

The result is:

Setting the coefficients of each order of $q$ to 0, we solve successively:

Therefore,

That is,

resulting in

This has $\mathcal{O}(x_n^{-3})$ accuracy. Slightly changing the idea, we can write an automatic calculation program in Mathematica:

g[x_] = 0; 
ff[x_] = 3*1/x + 1/x^2 
Do[e = q^n; 
g[x_] = g[x] + 
q^n*Integrate[-D[ 
f[x + 3*q*e + ff[x/q]*e*q, q] - f[x, q] + 
g[x + 3*q + ff[x/q]*q] - g[x] + ff[x/q], {q, 
n + 1}] /. {q -> 0, f -> 0} // Simplify, x]/ 
3/(n + 1)!;, {n, 0, 5}] 
h[x_] = g[x] /. {q -> 1} // Expand

By modifying 5 in {n, 0, 5}, the accuracy can be increased or decreased. The above code gives:

Asymptotic Results and Limit Value

At this point, we solve

Let

we get

To obtain an explicit asymptotic expression, an iterative formula can be built based on the above equation:

Here $x^{(k)}_n$ refers to the $k$-th order approximation of $x_n$. Taking $x^{(0)}_n=3n+a-\ln 3$ as the initial value, iterating and expanding gives:

This is the same as the asymptotic expression given by mathe on the research forum. The Mathematica code for the iteration is (continued from before):

p[n_] = 3*n + a - Log[3]; 
Do[p[n_] = 
Normal[Series[-h[p[n]] + 3*n + a - Log[3], {n, Infinity, 5}]], {i, 
0, 5}] 
Series[p[n], {n, Infinity, 5}]

From this, it is found that $a$ here is exactly the limit value $a$:

This also provides an asymptotic expression for calculating $a$ at the same time:

Using $x_{10^9}$ calculated by wayne and substituting it into the above formula, we get