{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Points périodiques de la fonction $x\\mapsto 4x(1-x)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "from math import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le but de ce notebook est l'étude de points remarquables de la fonction $f:[0,1]\\to[0,1]$ définie par $f(x)=4x(1-x)$ : ses __points périodiques__. Nous allons montrer et illustrer le théorème suivant :\n",
    "\n",
    "__Théorème__ : Pour tout entier $n\\ge 1$, $f$ possède un point périodique de période $n$.\n",
    "\n",
    "Mieux, nous allons, pour tout entier $n$, estimer le nombre $\\tau_n$ de points périodiques de $f$ de période $n$.\n",
    "\n",
    "__Théorème__ : $\\tau_n\\sim 2^n$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mais allons-y petit à petit ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def f(x): return 4 * x * (1 - x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "xs = [k / 500 for k in range(500)]\n",
    "plt.plot(xs, [f(x) for x in xs])\n",
    "plt.axis([0, 1, 0, 1])\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Période d'un point"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dans tout ce qui suit, $f^k$ désigne la $k$ième __itérée__ de $f$ : $f^k=f\\circ \\ldots\\circ f$, $k$ fois. Bien entendu, $f^0=id$ et $f^1=f$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Définition__ : Soit $x\\in[0,1]$. On dit que $x$ est __périodique__ lorsqu'il existe $k>0$ tel que $f^k(x)=x$. Le plus petit entier $k$ vérifiant cette propriété est appelé la __période__ de $x$, et est noté $T(x)$, ou simplement $T$ s'il n'y a pas de confusion à craindre."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Par exemple, un point est de période 1 si et seulement si c'est un __point fixe__ de $f$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Proposition 1__ : Soit $x\\in[0,1]$ un point périodique. Soit $n\\in\\mathbb N$. Alors $f^n(x)=x$ si et seulement si $T(x)| n$.\n",
    "\n",
    "__Démonstration__ : Appelons $T$ la période de $x$. On a évidemment $f^0(x)=id(x)=x$. Soit $k\\in\\mathbb N$. Supposons que $f^{kT}(x)=x$. On a alors $f^{(k+1)T}(x)=f^{T+kT}(x)=f^T(f^{kT}(x))=f^T(x)=x$. On vient de montrer par récurrence que pour tout $n$ multiple de $T$, $f^n(x)=x$.\n",
    "\n",
    "Soit maintenant $n\\in\\mathbb N$. Supposons que $f^n(x)=x$. Effectuons la division euclidienne de $n$ par $T$ : $n=kT + r$ où $0\\le r<T$. On a $x=f^n(x)=f^r(f^{kT}(x))=f^r(x)$. Par minimalité de $T$, $r=0$ et donc $n=kT$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Notations__ : \n",
    "\n",
    "- Pour tout entier $n\\ge 1$, on note $P(n)$ l'ensemble des points périodiques de période $n$. Le but premier de ce notebook est de prouver que $P(n)\\ne\\emptyset$. \n",
    "- On note également $Q(n)$ l'ensemble des points périodiques dont la période divise $n$. On a évidemment $P(n)\\subset Q(n)$.\n",
    "- On note enfin $\\tau_n$ le cardinal de $P(n)$, c'est à dire le nombre de points périodiques de période $n$. Le but second de ce notebook est de trouver un équivalent de $\\tau_n$ lorsque $n$ tend vers l'infini.\n",
    "- Le but troisième de ce notebook est de nous amuser :-)\n",
    "\n",
    "D'après la proposition précédente, $Q(n)=\\{x\\in[0,1], f^n(x)=x\\}$.\n",
    "\n",
    "__Proposition 2__ : Pour tous $m,n\\ge 1$, $P(m)\\cap P(n)=\\emptyset$. Et $Q(m)\\cap Q(n)=Q(m\\land n)$.\n",
    "\n",
    "__Démonstration__ : laissée au lecteur. La première égalité est triviale, la seconde utilise la proposition 1."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Proposition 3__ : $P(n)= Q(n)\\setminus \\bigcup_{d|n, d\\ne n} P(d)$.\n",
    "\n",
    "__Démonstration__ : Facile aussi. Pour obtenir $P(n)$ on enlève à $Q(n)$ les éléments dont la \"vraie\" période n'est pas $n$, mais seulement l'un de ses diviseurs stricts."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Corollaire 4__ : $\\tau_n = |Q(n)| - \\sum_{d|n, d\\ne n} \\tau_d$.\n",
    "\n",
    "__Démonstration__ : En effet, la réunion qui apparaît ci-dessus est une union disjointe.\n",
    "\n",
    "Il nous reste évidemment à calculer $Q(n)$ et son cardinal. Nous disposerons alors d'une relation de récurrence permettant de calculer $P(n)$ et $\\tau_n$ pour tout entier $n\\ge 1$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Un calcul explicite des points périodiques de $f$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2. Passage aux sinus"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Proposition 5__ : Soit $x\\in[0,1]$. On pose $x=\\sin^2\\theta$ où $\\theta\\in[0,\\frac\\pi 2]$. On a alors\n",
    "\n",
    "$$\\forall n\\in\\mathbb N, f^n(x)=\\sin^2(2^n\\theta)$$\n",
    "\n",
    "__Démonstration__ : Récurrence sur $n$. C'est évident pour $n=0$. Supposons cette égalité vraie pour $n$. On a alors $f^{n+1}(x)=f(\\sin^2(2^n\\theta))=4\\sin^2(2^n\\theta)(1-\\sin^2(2^n\\theta))=(2\\sin(2^n\\theta)\\cos(2^n\\theta))^2$ d'où le résultat."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Avec les notations ci-dessus, soit $n\\ge 1$. On a $f^n(x)=x$ si et seulement si $\\sin(2^n\\theta)=\\pm\\sin\\theta$, c'est à dire $2^n\\theta\\equiv \\pm\\theta[\\pi]$, ou encore $(2^n\\pm 1)\\theta\\equiv 0[\\pi]$. On en tire $\\theta=\\frac{k\\pi}{2^n\\pm 1}$, où $k\\in\\mathbb Z$. N'oublions pas la condition $\\theta\\in[0,\\frac\\pi 2]$ ! Cela nous donne finalement :\n",
    "\n",
    "- $\\theta=\\frac{k\\pi}{2^n+ 1}$, avec $0\\le k\\le \\lfloor\\frac{2^n+1}{2}\\rfloor$, ou bien\n",
    "- $\\theta=\\frac{k\\pi}{2^n- 1}$, avec $0\\le k\\le \\lfloor\\frac{2^n-1}{2}\\rfloor$\n",
    "\n",
    "Seul $\\theta=0$ appartient aux deux listes ci-dessus. Ainsi,\n",
    "\n",
    "$$Q(n)=\\{\\sin^2(\\frac{k\\pi}{2^n+ 1}),0\\le k\\le \\lfloor\\frac{2^n+1}{2}\\rfloor\\}\\bigcup\\{\\\\sin^2(\\frac{k\\pi}{2^n- 1}), 0< k\\le \\lfloor\\frac{2^n-1}{2}\\rfloor\\}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Proposition 6__ : Pour tout $n\\ge 1$, le cardinal de $Q(n)$ est $2^n$.\n",
    "\n",
    "__Démonstration__ : Comptez ses éléments :-)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La fonction `Q` renvoie l'ensemble $Q(n)$ des points périodiques (enfin, des valeurs approchées) de période divisant l'entier $n$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2 L'ensemble $Q(n)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def Q(n):\n",
    "    a = 2 ** n + 1\n",
    "    s1 = [sin(k * pi / a) ** 2 for k in range(0, int(a / 2) + 1)]\n",
    "    a = 2 ** n - 1\n",
    "    s2 = [sin(k * pi / a) ** 2 for k in range(1, int(a / 2) + 1)]\n",
    "    s = s1 + s2\n",
    "    return set(s)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(Q(1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(Q(3))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(Q(4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "len(Q(4))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.3 L'ensemble P(n)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Un peu d'arithmétique"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour coder $\\tau_n$ nous aurons besoin de quelques fonctions arithmétiques.\n",
    "\n",
    "La fonction `smallest_divisor` renvoie le plus petit diviseur de $n$ supérieur ou égal à 2. Ce diviseur est évidemment premier."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def smallest_divisor(n):\n",
    "    p = 2\n",
    "    while p * p <= n and n % p != 0: p += 1\n",
    "    if p * p > n: return n\n",
    "    else: return p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "smallest_divisor(221)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La fonction `factor` renvoie la liste des couples $(p, c)$ où $p$ est un diviseur premier de $n$ et $c$ est le plus grand entier tel que $p^c|n$. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def factor(n):\n",
    "    s = []\n",
    "    while n != 1:\n",
    "        p = smallest_divisor(n)\n",
    "        c = 0\n",
    "        while n % p == 0:\n",
    "            n = n // p\n",
    "            c += 1\n",
    "        s.append((p, c))\n",
    "    return s"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "factor(2214)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "2 * 3**3 * 41"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Et enfin, la fonction `divisors` renvoie la liste des diviseurs de l'entier $n$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def divisors(n):\n",
    "    f = factor(n)\n",
    "    return divisors_aux(f)\n",
    "\n",
    "\n",
    "def divisors_aux(f):\n",
    "    if f == []: return [1]\n",
    "    else:\n",
    "        ds = divisors_aux(f[1:])\n",
    "        p, c = f[0]\n",
    "        return [d * p ** k for d in ds for k in range(c + 1)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "divisors(2214)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. L'ensemble P(n)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous avons maintenant tout ce qu'il faut pour écrire une fonction qui renvoie l'ensemble $P(n)$. Cette fonction est récursive, elle utilise l'égalité $P(n)= Q(n)\\setminus \\bigcup_{d|n, d\\ne n} P(d)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def P(n):\n",
    "    if n == 1: return set(Q(1))\n",
    "    else:\n",
    "        E = Q(n)\n",
    "        for d in divisors(n):\n",
    "            if d != n:\n",
    "                E = E.difference(P(d))\n",
    "        return E"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "P(3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Le nombre exact de points périodiques"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous avons montré que pour tout $n\\ge 1$, \n",
    "\n",
    "$$\\tau_n=2^n-\\sum_{d|n,d\\ne n}\\tau_d$$\n",
    "\n",
    "La fonction récursive $\\tau$ ne fait que lire la formule."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def tau(n):\n",
    "    s = 2 ** n\n",
    "    if n == 1: return 2\n",
    "    for d in divisors(n):\n",
    "        if d != n:\n",
    "            s = s - tau(d)\n",
    "    return s"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "[(n, tau(n)) for n in range(1, 21)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Proposition 7__ : Pour tout $n\\ge 1$, $\\tau_n \\le 2 ^ n$.\n",
    "\n",
    "__Démonstration__ : c'est évident, puisque $P(n)\\subset Q(n)$ et $|Q(n)|=2^n$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nous y sommes ! Voici le résultat annoncé au début du notebook."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Théorème 1__ : Pour tout $n\\ge 1$, $\\tau_n > 0$.\n",
    "\n",
    "__Démonstration__ : Tout d'abord, $\\tau_1=2>0$. Soit maintenant $n>1$. On a $\\tau_n=2^n-\\sum_{d|n,d\\ne n}\\tau_d$. Cela dit, un diviseur strict $d$ de $n$ vérifie $d\\le \\frac n 2$. Donc $\\sum_{d|n,d\\ne n}\\tau_d\\le \\sum_{d=1}^{\\lfloor n/2\\rfloor}\\tau_d$, la somme étant étendue à tous les entiers, et pas seulement aux diviseurs de $n$. Mais $\\sum_{d=1}^{\\lfloor n/2\\rfloor}\\tau_d\\le \\sum_{d=0}^{\\lfloor n/2\\rfloor}2^d = 2^{\\lfloor n/2\\rfloor+1}-1$. Ainsi, $\\tau_n\\ge 2 ^n -2^{\\lfloor n/2\\rfloor+1}+1$. Ce dernier minorant est clairement strictement positif.\n",
    "\n",
    "En fait, ce minorant est même très grand. Précisément,\n",
    "\n",
    "__Théorème 2__ : $\\tau_n \\sim 2^n$.\n",
    "\n",
    "__Démonstration__ : $2 ^n -2^{\\lfloor n/2\\rfloor+1}+1\\le \\tau_n\\le 2^n$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "__Théorème 3__ : $0\\le 2^n-\\tau_n\\le2 \\sqrt{2^n}$.\n",
    "\n",
    "__Démonstration__ : Reprendre la double inégalité de la preuve précédente."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Nos majorations et minorations ont été assez grossières. Il est fort probable que nous puissions faire mieux, comme nous allons nous en convaincre en faisant quelques simulations. Voici tout d'abord le nombre de points de période $n$ pour $1\\le n\\le 20$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 21\n",
    "plt.plot(list(range(1, N)), [tau(n) for n in range(1, N)])\n",
    "#plt.plot(list(range(1, N)), [2 ** n for n in range(1, N)])\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Et maintenant, $2^n-\\tau_n$ et son majorant $2\\sqrt{2^n}$ pour $90\\le n\\le 100$. On trace aussi, en bleu, $\\sqrt{2^n}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 101\n",
    "plt.plot(list(range(1, N)), [2 ** n - tau(n) for n in range(1, N)], 'k')\n",
    "plt.plot(list(range(1, N)), [2 * sqrt(2 ** n) for n in range(1, N)], 'r')\n",
    "plt.plot(list(range(1, N)), [sqrt(2 ** n) for n in range(1, N)], 'b')\n",
    "plt.axis([90, 100, 0, 1.2e15])\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Clairement, la courbe bleue a l'air \"mieux\" que la rouge.  On voit qu'il reste un petit travail à faire pour affiner notre majoration ... pour les courageux, voici l'affinage."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notons $d_1=1 < d_2 < \\ldots < d_r$ les diviseurs de $n$. On a $d_1\\ge 1,\\ldots,d_r\\ge r$, l'entier $r$ étant bien évidemment inférieur ou égal à $n$. Mais $\\frac n d_1,\\ldots,\\frac n d_r$ sont aussi les diviseurs de $n$, cette fois-ci rangés dans l'ordre strictement décroissant ! Et $\\frac n d_1 = n$, $\\frac n d_2 \\le \\frac n 2$,..., $\\frac n d_r \\le \\frac n r$. Ainsi, $\\sum_{d|n,d\\ne n}\\tau_d\\sum_{d|n,d\\ne n}2^d\\le\\sum_{k=2}^r 2^{n/k}\\le \\sum_{k=2}^n 2^{n/k}$. Cette dernière somme est inférieure à $2^{n/2} + (n-2)2^{n/3}=2^{n/2}+o(2^{n/2})$.\n",
    "\n",
    "Ainsi, $2^n\\ge \\tau_n\\ge 2^n-2^{n/2}+o(2^{n/2})$, ou encore \n",
    "\n",
    "$$0\\le 2^n-\\tau_n\\le 2^{n/2}+o(2^{n/2})$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La bonne courbe c'est la bleue !"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La fonction $f$ possède donc beaucoup de points périodiques. En fait, __vraiment__ beaucoup : je ne résiste pas à la tentation de donner un dernier résultat, en guise de conclusion.\n",
    "\n",
    "__Théorème $\\infty$__ : L'ensemble des points périodiques de $f$ est dense dans $[0,1]$.\n",
    "\n",
    "__Démonstration__ : Je vous laisse en exercice la première étape : l'ensemble $E=\\{\\frac {k\\pi}{2^n+1}, n\\ge 1, 0\\le k\\le \\lfloor \\frac{2^n+1}{2}\\rfloor\\}$ est dense dans $[0, \\frac\\pi 2]$. Inspirez-vous du cours sur les réels. On a fait des exercices similaires.\n",
    "\n",
    "Maintenant, soit $y\\in[0,1]$. Il existe $x\\in[0,\\frac\\pi 2]$ tel que $y=\\sin^2 x$. Par la densité de $E$, il existe une suite $(u_n)$ de points de $E$ telle que $u_n\\to x$ lorsque $n$ tend vers l'infini. Par continuité, $\\sin^2 u_n \\to \\sin^2 x=y$ lorsque $n$ tend vers l'infini. Mais les $u_n$ étant dans $E$, $\\sin^2 u_n$, comme on l'a vu plus haut, est un point périodique de $f$. Ainsi, $y$ est limite d'une suite de points périodiques."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Conclusion"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Peut-être regarderez vous dorénavant la stupide parabole du début du notebook d'un oeil moins désabusé :-) ?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}