Computing infinite series with iterated functions

I found this question (linked here) which asks to find what this infinite series converges to

$$\sum_{n=1}^{\infty} \int_0^{\pi} f_n(x) dx$$

where $$f_{n+1}(x) = \sin(f_n(x))$$ and $$f_1 = \sin(x)$$. I actually wrote a few scripts which I figure could be optimized but I'm looking for a better approach or if the type of problem is known. I've linked the script here and posted it

import scipy.integrate as integrate
import math

def repeated(f,  n):
if n < 1:
raise ValueError()
elif n == 1:
return f
else:
return lambda x: repeated(f, n-1)(f(x))

def get_sum(n):
### return nth sum

sum_ = 0
lower_bound = 0
upper_bound = math.pi

for n in range(1,n+1):
sum_ += integrate.quad(repeated(lambda x: math.sin(x), n), lower_bound, upper_bound)
return sum_

tolerance = 10**-4
initial_val = get_sum(1)
k = 2
while True:
new_val = get_sum(k)
absolute_error = abs(new_val-initial_val)
print("K : {k}\nKth sum : {s}\nAbsolute Error : {e}".format(k=k,s=new_val, e=absolute_error))
if absolute_error < tolerance:
break
initial_val = new_val
k+=1

I can run this until around k=500 and it starts to get slow. I modified it and increased the step size and I hit the max recursion depth for compose at around k=2500. I then attempted to interpolate the integral values with an exponential function instead. This seemed to be a way around the recursion issue but it doesn't seem accurate. Any suggestions would be great.

• Are you sure it converges? I appears that sequence formed by the function evaluated at $\pi/2$ (the maximum of the function) grows monotonically. Jun 21 '21 at 23:25
• Sorry, I meant "the sequence formed by the partial sums". I tried with 1000 000 terms. I think that sequence gives an estimate for a lower bound since you can come the area of the triangle that lies below your convex function as $f_{\max} \pi/2$. Jun 22 '21 at 0:03
• I didn't used a recursive function. I stored the last two values and accumulated them. Jun 22 '21 at 0:20
• I don't believe there is a finite limit here because for each term the integral $\int_0^{\pi} f_n(x) dx > \int_0^{\epsilon} f_n(x) dx$, for any positive $\epsilon < \pi$. So if we take $\epsilon \ll 1$ then $f_n(x) \to x$ and the sum becomes $\sum_1^n \epsilon^2/2 = n \epsilon^2/2$ which diverges for $n \to \infty$. Jun 22 '21 at 2:46
• On the function sequence involved, the iterated sine, or any other of the form $y_{k+1}=y_k-ay_k^3+...$, see math.stackexchange.com/questions/1449281/…, math.stackexchange.com/questions/1609995/…, math.stackexchange.com/questions/105452/… Jun 25 '21 at 7:18

There is a lot to unpack here, and probably this is a better question for math.SE. TL; DR version is: in exact arithmetic, this does not converge; see Maxim Umansky's answer. In FP arithmetic, it will converge. I don't know to what, see the long version for an attempt.

Let's assume that it is possible to compute $$\int_0^{\pi} f_n(x)\text{d}x$$ exactly. I don't want to deal with the extra arithmetic coming from there, and it doesn't change the final answer much anyways.

After a while $$\int_0^{\pi} f_n(x)\text{d}x$$ will become small enough that $$fl(\sum_{n=1}^N \int_0^{\pi} f_n(x) \text{d}x+ \int_0^{\pi} f_{N+1}(x)\text{d}x)= fl(\sum_{n=1}^N \int_0^{\pi} f_n(x)\text{d}x)$$ where $$fl$$ is the rounding function. This is a common phenomenon in computational mathematics that you should be careful about.

(A well-known example is $$\sum_{n=1}^{\infty} 1/n$$ which diverges in exact arithmetic, but converges in FP arithmetic. If you reverse the ordering of the sum $$\sum_{n=M}^{1} 1/n$$ for $$M\gg 1$$, then you get the expected divergence result.)

Also, $$f_n(x)$$ itself can be thought as a fixed point iteration, and the integral $$\int_0^{\pi} f_n(x)$$ can be interpreted as testing the convergence of the fixed point iterations for all the initial guesses in the interval $$[0,\pi]$$ at the $$n$$-th iteration. The fixed point iteration $$x=\sin(x)$$ has the solution $$x=0$$ for any initial guess $$x_0\in [0,\pi]$$ (not too hard to prove, so I will skip). Hence, the integral $$\int_0^{\pi} f_n(x)\text{d}x \to 0$$ as $$n\to\infty$$. The upper bound for this integral is $$\pi f_n(\pi/2)$$.

Now, we can naively bound $$\sin^n(\pi/2)$$ ($$\sin^n$$ means $$\sin$$ composited with itself $$n$$ times) from below by $$1/n$$. Here is a proof by induction:

$$n=1$$ case: $$\sin^1(\pi/2) \geq 1/1$$. Hence, we are good here.

Inductive hypothesis: Assume for $$n=k$$, $$\sin^n(\pi/2)\geq 1/n$$.

$$n=k+1$$ case: Consider \begin{align}\sin(\sin^k(\pi/2)) - 1/(k+1) &= \sin(\sin^k(\pi/2)) - 1/k + 1/(k^2+k) \\ &\text{by the inductive hypothesis and since sine is an increasing func. in } [0,1]\\ &\geq \sin(1/k) - 1/k + 1/(k^2+k) \\ &\text{by Taylor's remainder theorem}\\ &\geq 1/k -1/(6k^3) -1/k +1/(k^2+k)\\ &\text{for any } k\geq 1\\ &\geq 0. \end{align}

So, the lower bound for the sum $$\pi\sum_{n=1}^{\infty} f_n(\pi/2)$$ is $$\pi\sum_{n=1}^{\infty} 1/n$$, which will definitely not increase after $$n=4.5\times 10^{15}\approx 1/\varepsilon_{\text{mach}}$$ since $$fl(1+\varepsilon_{\text{mach}}) = fl(1)$$. I would expect $$\pi\sum_{n=1}^{\infty} f_n(\pi/2)$$ to "converge" earlier than that (the word "converge" is in quotation marks, because it is not a real convergence).

On the other hand, we can use $$1$$ as an upper bound to the upper bound $$f_n(\pi/2)$$ (since the fixed point iterations converge to zero). Then we see that $$4.5\pi\times 10^{15}\approx\pi\sum_{n=1}^{\infty} 1 \geq \pi\sum_{n=1}^{\infty} f_n(\pi/2) \geq \pi\sum_{n=1}^{\infty} 1/n \approx \pi(\ln(4.5\times 10^{15}) + \gamma),$$ where $$\gamma$$ is the Euler–Mascheroni constant. (Approximation signs in the equation above mean floating point approximation, and it is non-standard notation. I am abusing the notation here)

Hence, the sum $$\sum_{n=1}^{\infty} \int_0^{\pi} f_n(x)$$ definitely has a finite limit when calculated on a computer using floating point arithmetic and that value is no larger than $$4.5\pi\times 10^{15}$$. I would expect it to be slightly larger $$\pi\ln(4.5\times 10^{15})$$, but I cannot prove that rigorously. That is more of an educated guess.

Edit: If Maxim Umansky's answer is not clear to you, you can look for a lower bound for the integral $$\int_{0}^{\pi} f_n(x)\text{d}x$$ using quadrature rules. A simple to prove lower bound comes from the observation that $$f_n$$ is concave (see Theorems 1 and 2) for all finite $$n$$ so $$\int_{0}^{\pi} f_n(x)\text{d}x \geq f_n(\pi/2)\pi/2 \geq \pi/(2n)$$. The first inequality basically says that the composite trapezoidal rule applied to the integral using the subintervals $$[0,\pi/2]$$ and $$[\pi/2,\pi]$$ will be a lower bound to the value of the integral. Then using the crude lower bound to $$f_n(\pi/2)$$, we can say that $$\int_{0}^{\pi} f_n(x)\text{d}x\geq \frac{\pi}{2} \frac{1}{n}$$. Hence, $$\sum_{n=1}^{\infty}\int_{0}^{\pi} f_n(x)\text{d}x\geq \frac{\pi}{2} \sum_{n=1}^{\infty}\frac{1}{n}$$ which diverges.

Theorem 1: If $$h(x)$$ is concave and non-decreasing, and $$g(x)$$ is concave, then $$f(x) = h(g(x))$$ is concave.

Theorem 2: If $$h(x)$$ is concave and non-increasing, and $$g(x)$$ is concave, then $$f(x) = h(g(x))$$ is concave.

• Thank you. It wasn't that clear at first. This helps. Jun 22 '21 at 10:43
• No worries. You understand why your Python code, even it was lightning quick, wouldn't give you the correct answer, right? Jun 22 '21 at 15:01
• I do now. I probably jumped the gun trying to calculate it first but it's a lot more complex than I thought. Jun 22 '21 at 15:28

There is no finite limit for this series sum. Note that for each $$n$$ the function $$f_n$$ is positive definite, $$f_n(x) > 0$$ within the semi-open interval $$(0,\pi]$$, and we can construct the lower bound for the sum as follows. Consider some parameter $$\epsilon \in (0,\pi]$$. Then for each $$n$$ the integral $$\int_0^{\pi} f_n(x) dx \ge \int_0^{\epsilon} f_n(x) dx$$. But if we take $$\epsilon \ll 1$$ then for each $$n$$ the function $$f_n(x) \to x$$ and the sum $$\sum_1^n \int_0^{\epsilon} f_n(x) dx \to \sum_1^n \epsilon^2/2 = n \epsilon^2/2$$ which diverges for $$n \to \infty$$.