# Méthode du cercle de séparation

En mathématiques, la méthode du cercle de séparation est un algorithme numérique de recherche des racines complexes d'un polynôme. Il fut présenté par Arnold Schönhage dans sa publication de 1982 le théorème fondamental de l'algèbre en termes de complexité de calcul (rapport technique, Mathematisches Institut der Universität Tübingen). Une application de l'algorithme réalisée par Xavier Gourdon est employée par le logiciel de calcul formel Magma.

## Description générale

L'idée fondamentale de la méthode du cercle de séparation est d'utiliser les méthodes de l'analyse complexe, plus précisément le théorème des résidus, pour construire des facteurs de polynômes. Avec ces méthodes, il est possible de construire un facteur d'un polynôme donné $p(x)=x^{n}+p_{n-1}x^{n-1}+\cdots +p_{0}$ pour toute région du plan complexe avec une frontière lisse par morceaux. La plupart de ces facteurs seront triviaux, c'est-à-dire des polynômes constants. Seules les régions qui contiennent des racines de p(x) donnent lieu à des facteurs non triviaux qui ont exactement ces racines de p(x) comme leurs propres racines, en préservant la multiplicité.
Dans la réalisation numérique de cette méthode, on utilise des disques D(c,r) (centre c, rayon r) dans le plan complexe comme régions. Le cercle limite d'un disque divise l'ensemble des racines de p(x) en deux parties, d'où le nom de la méthode. Pour un disque donné, on calcule des facteurs approchés en suivant la théorie analytique et on les affine en utilisant la méthode de Newton. Pour éviter l'instabilité numérique, il faut exiger que toutes les racines soient bien séparées du cercle limite du disque. Ainsi, pour obtenir un bon cercle de séparation, il doit être intégré dans un anneau sans racines A(c,r,R) (centre c, rayon intérieur r, rayon extérieur R) avec une grande largeur relative R/r.
En répétant ce processus pour les facteurs trouvés, on arrive finalement à une factorisation approximative du polynôme avec une précision requise. Les facteurs sont soit des polynômes linéaires représentant des zéros bien isolés, soit des polynômes d'ordre supérieur représentant des groupes de zéros.

## Détails de la construction analytique

Les identités de Newton sont une relation bijective entre les polynômes symétriques élémentaires d'un ensemble de nombres complexes et ses sommes de puissances. Par conséquent, il est possible de calculer les coefficients d'un polynôme $p(x)=x^{n}+p_{n-1}x^{n-1}+\cdots +p_{0}=(x-z_{1})\cdots (x-z_{n})$ (ou d'un de ses facteurs) à partir des sommes des puissances de ses zéros $t_{m}=z_{1}^{m}+\cdots +z_{n}^{m}\,,\quad m=0,1,\cdots ,n$ en résolvant le système triangulaire obtenu en comparant les puissances de u dans l'identité suivante des séries formelles de puissances ${\begin{aligned}&a_{n-1}+2\,a_{n-2}\,u+\cdots +(n-1)\,a_{1}\,u^{n-2}+n\,a_{0}\,u^{n-1}\\[1ex]&\quad =-(1+a_{n-1}\,u+\cdots +a_{1}\,u^{n-1}+a_{0}\,u^{n})\cdot (t_{1}+t_{2}\,u+t_{3}\,u^{2}+\dots +t_{n}\,u^{n-1}+\cdots ).\end{aligned}}$
Si $G\subset \mathbb {C}$ est un domaine dont la frontière C est lisse par morceaux et si les zéros de p(x) sont distincts deux à deux et ne sont pas sur la frontière C, alors d'après le théorème des résidus du calcul résiduel, on obtient ${\frac {1}{2\pi \,i}}\oint _{C}{\frac {p'(z)}{p(z)}}z^{m}\,\mathrm {d} z=\sum _{z\in G:\,p(z)=0}{\frac {p'(z)z^{m}}{p'(z)}}=\sum _{z\in G:\mathbb {N} ,p(z)=0}z^{m}.$
L'identité du côté gauche au côté droit de cette équation est également valable pour les zéros avec des multiplicités. En utilisant les identités de Newton, on peut calculer à partir de ces sommes de puissances le facteur $f(x):=\prod _{z\in G:\,p(z)=0}(x-z)$ de p(x) correspondant aux zéros de p(x) à l'intérieur de G. Par division polynomiale, on obtient également le second facteur g(x) dans p(x) = f(x)g(x).
Les régions couramment utilisées sont des cercles dans le plan complexe. Chaque cercle donne lieu à une division du polynôme p(x) en facteurs f(x) et g(x). En répétant cette procédure sur les facteurs en utilisant différents cercles, on obtient des factorisations de plus en plus fines. Cette récursion s'arrête après un nombre fini de divisions correctes, tous les facteurs étant des puissances non triviales de polynômes linéaires.
Le défi consiste maintenant à convertir cette procédure analytique en un algorithme numérique avec un bon temps d'exécution. L'intégration est approximée par une somme finie d'une méthode d'intégration numérique, utilisant la transformée de Fourier rapide pour l'évaluation des polynômes p(x) et p'(x). Le polynôme f(x) qui en résulte ne sera qu'un facteur approximatif. Pour s'assurer que ses zéros sont proches des zéros de p à l'intérieur de G et seulement de ceux-ci, il faut exiger que tous les zéros de p soient éloignés de la frontière C de la région G.

## Observation numérique basique

(Schönhage 1982) Soit $p\in \mathbb {C} [X]$ un polynôme de degré n qui a k zéros à l'intérieur du cercle de rayon 1/2 et les n – k zéros restants à l'extérieur du cercle de rayon 2. Avec N= O(k) suffisamment grand, l'approximation des intégrales de contour en utilisant N points résulte en une approximation f$_0$ du facteur f avec l'erreur $\|f-f_{0}\|\leq 2^{2k-N}\,nk\,100/98,$ où la norme d'un polynôme est la somme des modules de ses coefficients.
Puisque les zéros d'un polynôme dépendent continûment de ses coefficients, on peut rendre les zéros de f$_0$ aussi proches que possible des zéros de f en choisissant N suffisamment grand. Cependant, on peut améliorer cette approximation plus rapidement en utilisant une méthode de Newton. La division de p avec le reste donne une approximation g$_0$ du facteur restant g. Maintenant $p-f_{0}g_{0}=(f-f_{0})g_{0}+(g-g_{0})f_{0}+(f-f_{0})(g-g_{0}),$ en écartant le dernier terme du second ordre, on doit résoudre $p-f_{0}g_{0}=f_{0}\Delta g+g_{0}\Delta f$ en utilisant n'importe quelle variante de l'algorithme d'Euclide étendu pour obtenir les approximations incrémentées $f_{1}=f_{0}+\Delta f$ et $g_{1}=g_{0}+\Delta g$. Cette opération est répétée jusqu'à ce que les incréments soient nuls par rapport à la précision choisie.

## Itération de Graeffe

L'étape cruciale de cette méthode est de trouver un anneau de largeur relative 4 dans le plan complexe qui ne contient aucun zéro de p et qui contient approximativement autant de zéros de p à l'intérieur qu'à l'extérieur de celui-ci. Tout anneau de cette caractéristique peut être transformé, par translation et mise à l'échelle du polynôme, en l'anneau compris entre les rayons 1/2 et 2 autour de l'origine. Mais tous les polynômes n'admettent pas un tel anneau de séparation.
Pour remédier à cette situation, l'tération de Graeffe est appliquée. Elle calcule une suite de polynômes $p_{0}=p,\qquad p_{j+1}(x)=(-1)^{\deg p}p_{j}({\sqrt {x}})\,p_{j}(-{\sqrt {x}}),$ où les racines de $p_{j}(x)$ sont les 2$^j$-èmes puissances dyadiques des racines du polynôme initial p. En divisant $p_{j}(x)=e(x^{2})+x\,o(x^{2})$ en parties paires et impaires, le polynôme suivant est obtenu par des opérations purement arithmétiques comme $p_{j+1}(x)=(-1)^{\deg p}(e(x)^{2}-x\,o(x)^{2})$. Les rapports des modules absolus des racines augmentent de la même puissance 2$^j$ et tendent donc vers l'infini. En choisissant j suffisamment grand, on trouve finalement un anneau de séparation de largeur relative 4 autour de l'origine.
La factorisation approximative de $p_{j}(x)\approx f_{j}(x)\,g_{j}(x)$ doit maintenant être ramenée au polynôme original. Pour ce faire, on utilise une alternance de pas de Newton et d'approximations de Padé. Il est facile de vérifier que ${\frac {p_{j-1}(x)}{g_{j}(x^{2})}}\approx {\frac {f_{j-1}(x)}{g_{j-1}(-x)}}$ est valable. Les polynômes du côté gauche sont connus à l'étape j, les polynômes du côté droit peuvent être obtenus en tant qu'approximant de Padé des degrés correspondants pour le développement en série entière de la fraction du côté gauche.

## Trouver un bon cercle

En utilisant l'itération de Graeffe et toute estimation connue de la valeur absolue de la plus grande racine, on peut trouver des estimations R de cette valeur absolue de n'importe quelle précision. On calcule maintenant les estimations des plus grandes et plus petites distances $R_{j}>r_{j}>0$ de toute racine de p(x) à l'un des cinq points centraux 0, 2R, –2R, 2iR, -2iR et on choisit celui qui a le plus grand rapport $R_{j}/r_{j}$ entre les deux. Cette construction permet de garantir que $R_{j}/r_{j}>\mathrm {e} ^{0,3}\approx 1.35$ pour au moins un centre. Pour un tel centre, il doit y avoir un anneau sans racine de largeur relative e$^{0, 3/n}$≈1+0, 3/n {\textstyle \mathrm {e} ^{0, 3/n}\approx 1+{\frac {0, 3}{n}}}. Après 3+log$_2$⁡(n){\textstyle 3+\log _{2}(n)} itérations de Graeffe, l'anneau correspondant du polynôme itéré a une largeur relative supérieure à 11 > 4, comme requis pour la division initiale décrite ci-dessus (Schönhage 1982). Après 4+log$_2$⁡(n)+log$_2$⁡(2+log$_2$⁡(n)){\textstyle 4+\log _{2}(n)+\log _{2}(2+\log _{2}(n))} itérations de Graeffe, l'anneau correspondant a une largeur relative supérieure à 2$^{13, 8}$⋅n$^{6, 9}$>(64⋅n$^3$)$^2${\textstyle 2^{13{, }8}\cdot n^{6, 9}>(64\cdot n^{3})^{2}}, ce qui permet un fractionnement initial très simplifié (Malajovich et Zubelli 1997)
Pour localiser le meilleur anneau sans racine, on utilise une conséquence du théorème de Rouché : Pour k = 1, ..., n – 1 l'équation polynomiale $\mathbb {N} ,0=\sum _{j\neq k}|p_{j}|u^{j}-|p_{k}|u^{k},u>0$ a, par la règle des signes de Descartes, zéro ou deux racines positives $u_{k}<v_{k}$. Dans ce dernier cas, il y a exactement k racines à l'intérieur du disque (fermé) $D(0,u_{k})$ et $A(0,u_{k},v_{k})$ est un anneau (ouvert) sans racine.