Algoritmo di rappresentazione dei polinomi simmetrici negli $e_i$

Lo scopo di questo notebook è implementare l'algoritmo impiegato nella dimostrazione del Teorema fondamentale dei polinomi simmetrici, che, dato un campo $F$, asserisce il seguente isomorfismo:

$$\operatorname{Sym}[X_n] \cong F[e_1, \ldots, e_n],$$

dove $X_n$ è l'insieme $\{x_1, \ldots, x_n\}$.

Innanzitutto, si sceglie un $n \in \mathbb{N}^+$, che rappresenta il numero di variabili di cui è composto il polinomio simmetrico, e si crea l'anello di polinomi $\mathbb{Q}[x_1, \ldots, x_n]$, nel quale vale il degree lexicographic order (deglex). Vengono infine create delle variabili simboliche a rappresentare i polinomi simmetrici elementari.

%display latex
n = 5

P = PolynomialRing(QQ, ",".join(f"x_{i}" for i in range(1, n+1)), order='deglex')
P.inject_variables()

var(",".join(f"e_{i}" for i in range(n+1)))

Defining x_1, x_2, x_3, x_4, x_5

$$\newcommand{\Bold}[1]{\mathbf{#1}}\left(e_{0}, e_{1}, e_{2}, e_{3}, e_{4}, e_{5}\right)$$

Si definisce il valore che deve assumere la funzione $e(d)$, che rappresenta il polinomio simmetrico $e_d$, nel seguente modo:

$$e(d) = \begin{cases}1 & \text{se }d=0, \\ \displaystyle \sum_{1\leq i_1 < i_2 < \cdots < i_d \leq n} \underbrace{x_{i_1} \cdots x_{i_d}}_{d\text{ volte}} & \text{altrimenti.} \end{cases}$$

from functools import reduce
from itertools import combinations
from operator import mul

def e(d: Integer):
    if d == 0:
        return 1
    
    return sum(reduce(mul, [P(f"x_{j}") for j in i], 1) for i in combinations(range(1, n+1), d))

Si definisce una funzione $\beta(\operatorname{exp\_lt})$ che prende in ingresso una tupla ordinata $\operatorname{exp\_lt} \in \mathbb{N}^n$ e restituisce una tupla dello stesso tipo definita nel seguente modo:

$$\beta(\operatorname{exp\_lt})_i = \begin{cases} \operatorname{exp\_lt}_i - \operatorname{exp\_lt}_{i+1} & \text{se }1 \leq i < n, \\ \operatorname{exp\_lt}_i & \text{altrimenti.} \end{cases}$$

def beta(exp_lt):
    return [e - exp_lt[i+1] if i != n-1 else e for i, e in enumerate(exp_lt)]

Si definiscono due funzioni analoghe, dette $\operatorname{e\_prod}$ e $\operatorname{e\_prod\_value}$. La seconda restituisce lo stesso polinomio di $\operatorname{e\_prod}$ sostituendo ai simboli $e_i$ i corrispettivi valori $e(i)$. Pertanto si definisce solo il valore della prima funzione.

Tale funzione $\operatorname{e\_prod}$ prende in ingresso una tupla $b \in \mathbb{N}^n$ e restituisce $e^b = {e_1}^{b_1} {e_2}^{b_2} \cdots {e_n}^{b_n}$.

def e_prod(b):
    return reduce(mul, [eval(f"e_{i+1}")^k for i, k in enumerate(b)], 1)


def e_prod_value(b):
    return reduce(mul, [e(i+1)^k for i, k in enumerate(b)], 1)

Si definisce infine la funzione $\operatorname{combination}(\operatorname{poly}(x))$, che prende in ingresso un polinomio simmetrico $\operatorname{poly}(x)$ e ne restituisce la rappresentazione in $F[e_1, \ldots, e_n]$, secondo il seguente algoritmo.

• se $\operatorname{poly}(x)$ è $0$ o ha grado nullo, la sua rappresentazione in $F[e_1, \ldots, e_n]$ è già $\operatorname{poly}(x)$, e quindi la funzione restituisce il polinomio in ingresso senza modificarlo,
• altrimenti, si considera, secondo il deglex, il leading term di $\operatorname{poly}(x)$, detto $\operatorname{lt}$:
  • detta $\alpha$ la tupla ordinata degli esponenti di $\operatorname{lt}$, si calcola $\beta(\alpha)$, detto $b$,
  • detto $c$ il coefficiente razionale di $\operatorname{lt}$, si restituisce la rappresentazione simbolica $c \cdot \operatorname{e\_prod(b)}$, a cui si aggiunge ricorsivamente la rappresentazione del polinomio $\operatorname{poly(x)} -\, c \cdot \operatorname{e\_prod\_value(b)}$, ottenuta reiterando l'algoritmo su di esso.

def combination(poly):
    if isinstance(poly, Integer) or isinstance(poly, Rational) or poly == 0 or poly.degree() == 0:
        return poly
    
    lt = sorted(list(poly), reverse=True)[0]
    
    try:
        b = beta(lt[1].exponents()[0])
        return lt[0] * e_prod(b) + combination(poly - lt[0] * e_prod_value(b))
    except (AttributeError, TypeError):
        raise TypeError("The given polynomial is not of symmetric kind")

Si sceglie inoltre un polinomio $\operatorname{poly}(x)$ che deve essere simmetrico -- qualora non fosse tale, l'algoritmo non terminerà con successo (se terminasse con successo, il polinomio sarebbe combinazione di polinomi simmetrici, e sarebbe dunque anch'esso simmetrico, ↯).

# Il polinomio deve essere simmetrico...
poly = x_1 + x_2 + x_3 + x_4 + x_5 + (x_1^2 + x_2^2 + x_3^2 + x_4^2 + x_5^2 - x_1*x_2*x_3*x_4*x_5)^2
poly

$$\newcommand{\Bold}[1]{\mathbf{#1}}x_{1}^{2} x_{2}^{2} x_{3}^{2} x_{4}^{2} x_{5}^{2} - 2 x_{1}^{3} x_{2} x_{3} x_{4} x_{5} - 2 x_{1} x_{2}^{3} x_{3} x_{4} x_{5} - 2 x_{1} x_{2} x_{3}^{3} x_{4} x_{5} - 2 x_{1} x_{2} x_{3} x_{4}^{3} x_{5} - 2 x_{1} x_{2} x_{3} x_{4} x_{5}^{3} + x_{1}^{4} + 2 x_{1}^{2} x_{2}^{2} + 2 x_{1}^{2} x_{3}^{2} + 2 x_{1}^{2} x_{4}^{2} + 2 x_{1}^{2} x_{5}^{2} + x_{2}^{4} + 2 x_{2}^{2} x_{3}^{2} + 2 x_{2}^{2} x_{4}^{2} + 2 x_{2}^{2} x_{5}^{2} + x_{3}^{4} + 2 x_{3}^{2} x_{4}^{2} + 2 x_{3}^{2} x_{5}^{2} + x_{4}^{4} + 2 x_{4}^{2} x_{5}^{2} + x_{5}^{4} + x_{1} + x_{2} + x_{3} + x_{4} + x_{5}$$

Si calcola adesso la rappresentazione di $\operatorname{poly}(x)$ in funzione dei polinomi simmetrici elementari secondo l'implementazione di $\operatorname{combination}$:

combination(poly)

$$\newcommand{\Bold}[1]{\mathbf{#1}}e_{1}^{4} - 4 \, e_{1}^{2} e_{2} - 2 \, e_{1}^{2} e_{5} + 4 \, e_{2}^{2} + 4 \, e_{2} e_{5} + e_{5}^{2} + e_{1}$$

---

Data una somma di potenze simmetrica, è possibile implementare un algoritmo di rappresentazione negli $e_i$ secondo le identità di Newton-Girard. Si definisce allora la funzione $\operatorname{newton\_girard}(k)$, che prende in ingresso un $k \in \mathbb{N}$ e restituisce la rappresentazione di $\sum_{i=1}^n {x_i}^k$ nel seguente modo:

$$\operatorname{newton\_girard}(k) = \begin{cases} n & \text{se } k = 0, \\ (-1)^{k+1} \cdot (k e_{k} + \sum_{i=1}^{k-1} (-1)^i \operatorname{newton\_girard}(i) \, e_{k-1}) & \text{altrimenti}, \end{cases}$$

dove si pone $e_i = 0$ per ogni $i > n$.

La funzione è implementata con l'ausilio di un sistema di caching, affinché non vengano ricalcolati i valori già calcolati.

from functools import lru_cache

@lru_cache(maxsize=None)
def newton_girard(k):    
    if k == 0:
        return n
    
    return ((-1)**(k+1) * ((k * eval(f"e_{k}") if k <= n else 0) +
                           sum((-1)**i * newton_girard(i) * eval(f"e_{k - i}") for i in range(max(1, k-n), k)))).expand()

Si definisce infine la funzione $\operatorname{sum\_of\_powers}(k)$, che, dato in ingresso un $k \in \mathbb{N}$, restituisce la rappresentazione di $\sum_{i=1}^n (x_i)^k$ secondo la funzione $\operatorname{combination}$ già definita.

def sum_of_powers(k):
    return combination(sum(eval(f"x_{i}")^k for i in range(1, n+1)))

Si sceglie infine un $k \in \mathbb{N}$ al quale elevare tutti le variabili della somma.

k = 8

Si computa la rappresentazione di $\sum_{i=1}^n (x_i)^k$ prima secondo $\operatorname{newton\_girard}$, e poi secondo $\operatorname{sum\_of\_powers}$, e si verifica che siano uguali.

a = newton_girard(k)
a

$$\newcommand{\Bold}[1]{\mathbf{#1}}e_{1}^{8} - 8 \, e_{1}^{6} e_{2} + 20 \, e_{1}^{4} e_{2}^{2} + 8 \, e_{1}^{5} e_{3} - 16 \, e_{1}^{2} e_{2}^{3} - 32 \, e_{1}^{3} e_{2} e_{3} - 8 \, e_{1}^{4} e_{4} + 2 \, e_{2}^{4} + 24 \, e_{1} e_{2}^{2} e_{3} + 12 \, e_{1}^{2} e_{3}^{2} + 24 \, e_{1}^{2} e_{2} e_{4} + 8 \, e_{1}^{3} e_{5} - 8 \, e_{2} e_{3}^{2} - 8 \, e_{2}^{2} e_{4} - 16 \, e_{1} e_{3} e_{4} - 16 \, e_{1} e_{2} e_{5} + 4 \, e_{4}^{2} + 8 \, e_{3} e_{5}$$

b = sum_of_powers(k)
b

$$\newcommand{\Bold}[1]{\mathbf{#1}}e_{1}^{8} - 8 \, e_{1}^{6} e_{2} + 20 \, e_{1}^{4} e_{2}^{2} + 8 \, e_{1}^{5} e_{3} - 16 \, e_{1}^{2} e_{2}^{3} - 32 \, e_{1}^{3} e_{2} e_{3} - 8 \, e_{1}^{4} e_{4} + 2 \, e_{2}^{4} + 24 \, e_{1} e_{2}^{2} e_{3} + 12 \, e_{1}^{2} e_{3}^{2} + 24 \, e_{1}^{2} e_{2} e_{4} + 8 \, e_{1}^{3} e_{5} - 8 \, e_{2} e_{3}^{2} - 8 \, e_{2}^{2} e_{4} - 16 \, e_{1} e_{3} e_{4} - 16 \, e_{1} e_{2} e_{5} + 4 \, e_{4}^{2} + 8 \, e_{3} e_{5}$$

a - b

$$\newcommand{\Bold}[1]{\mathbf{#1}}0$$

---

Sia ora:

$$f(x) = x^n + a_{n-1} x^{n-1} \ldots + a_0,$$

con radici $x_1$, $x_2$, ..., $x_n$. Allora si definisce discriminante di $f$ la seguente produttoria:

$$\Delta(f) = \displaystyle \prod_{1 \leq i < j \leq n} (x_i - x_j)^2.$$

In particolare, vale la seguente proprietà:

$$\exists\, i \neq j \mid x_i = x_j \iff \Delta(f) = 0.$$

Si verifica facilmente che $\Delta(f)$ è un polinomio simmetrico. Poiché ogni permutazione è una composizione di trasposizioni, è sufficiente verificare che invertendo due radici $x_i$ e $x_j$ il discriminante rimanga invariato:

• i fattori che contengono solo uno tra $x_i$ e $x_j$ mantengono la somma della base invariata o al più cambiano di segno, e dunque, elevando al quadrato, rimangono invariati,
• il fattore che contiene sia $x_i$ che $x_j$ cambia la propria base di segno, ma rimane invariato elevando al quadrato.

Poiché $\Delta(f)$ è un polinomio simmetrico nelle variabili $x_1$, ..., $x_n$, per il Teorema fondamentale dei polinomi simmetrici, si scrive in modo unico in funzione dei polinomi simmetrici elementari $e_i(x_1, \ldots, x_n)$. Tuttavia tali polinomi simmetrici elementari altro non sono che i coefficienti di $f(x)$, a meno del segno. In particolare vale che:

$$e_i(x_1, \ldots, x_n) = (-1)^{i} a_{n-i}, \quad \text{per } 0 \leq i \leq n.$$

Si definiscono quindi i simboli $a_0$, ..., $a_n$ e si costruisce una funzione $\Delta(n)$ che restituisce il discriminante di un generico polinomio di $n$-esimo grado dato in ingresso in funzione degli $a_i$ mediante $\operatorname{combination}$.

var(",".join(f"a_{i}" for i in range(n+1)))

def delta(n):
    
    f = reduce(mul, (eval(f"(x_{i}-x_{j})")^2 for i, j in combinations(range(1, n+1), 2)), 1)
    c = combination(f)

    for i in range(1, n+1):
        if i % 2:
            c = eval(f"c.substitute(e_{i}=-a_{n-i})")
        else:
            c = eval(f"c.substitute(e_{i}=a_{n-i})")

    return c

Dunque, per l'$n$ scelto in partenza, il discriminante del generico polinomio di $n$-esimo grado è il seguente:

delta(n)

$$\newcommand{\Bold}[1]{\mathbf{#1}}a_{1}^{2} a_{2}^{2} a_{3}^{2} a_{4}^{2} - 4 \, a_{0} a_{2}^{3} a_{3}^{2} a_{4}^{2} - 4 \, a_{1}^{3} a_{3}^{3} a_{4}^{2} + 18 \, a_{0} a_{1} a_{2} a_{3}^{3} a_{4}^{2} - 27 \, a_{0}^{2} a_{3}^{4} a_{4}^{2} - 4 \, a_{1}^{2} a_{2}^{3} a_{4}^{3} + 16 \, a_{0} a_{2}^{4} a_{4}^{3} + 18 \, a_{1}^{3} a_{2} a_{3} a_{4}^{3} - 80 \, a_{0} a_{1} a_{2}^{2} a_{3} a_{4}^{3} - 6 \, a_{0} a_{1}^{2} a_{3}^{2} a_{4}^{3} + 144 \, a_{0}^{2} a_{2} a_{3}^{2} a_{4}^{3} - 27 \, a_{1}^{4} a_{4}^{4} + 144 \, a_{0} a_{1}^{2} a_{2} a_{4}^{4} - 128 \, a_{0}^{2} a_{2}^{2} a_{4}^{4} - 192 \, a_{0}^{2} a_{1} a_{3} a_{4}^{4} + 256 \, a_{0}^{3} a_{4}^{5} - 4 \, a_{1}^{2} a_{2}^{2} a_{3}^{3} + 16 \, a_{0} a_{2}^{3} a_{3}^{3} + 16 \, a_{1}^{3} a_{3}^{4} - 72 \, a_{0} a_{1} a_{2} a_{3}^{4} + 108 \, a_{0}^{2} a_{3}^{5} + 18 \, a_{1}^{2} a_{2}^{3} a_{3} a_{4} - 72 \, a_{0} a_{2}^{4} a_{3} a_{4} - 80 \, a_{1}^{3} a_{2} a_{3}^{2} a_{4} + 356 \, a_{0} a_{1} a_{2}^{2} a_{3}^{2} a_{4} + 24 \, a_{0} a_{1}^{2} a_{3}^{3} a_{4} - 630 \, a_{0}^{2} a_{2} a_{3}^{3} a_{4} - 6 \, a_{1}^{3} a_{2}^{2} a_{4}^{2} + 24 \, a_{0} a_{1} a_{2}^{3} a_{4}^{2} + 144 \, a_{1}^{4} a_{3} a_{4}^{2} - 746 \, a_{0} a_{1}^{2} a_{2} a_{3} a_{4}^{2} + 560 \, a_{0}^{2} a_{2}^{2} a_{3} a_{4}^{2} + 1020 \, a_{0}^{2} a_{1} a_{3}^{2} a_{4}^{2} - 36 \, a_{0} a_{1}^{3} a_{4}^{3} + 160 \, a_{0}^{2} a_{1} a_{2} a_{4}^{3} - 1600 \, a_{0}^{3} a_{3} a_{4}^{3} - 27 \, a_{1}^{2} a_{2}^{4} + 108 \, a_{0} a_{2}^{5} + 144 \, a_{1}^{3} a_{2}^{2} a_{3} - 630 \, a_{0} a_{1} a_{2}^{3} a_{3} - 128 \, a_{1}^{4} a_{3}^{2} + 560 \, a_{0} a_{1}^{2} a_{2} a_{3}^{2} + 825 \, a_{0}^{2} a_{2}^{2} a_{3}^{2} - 900 \, a_{0}^{2} a_{1} a_{3}^{3} - 192 \, a_{1}^{4} a_{2} a_{4} + 1020 \, a_{0} a_{1}^{2} a_{2}^{2} a_{4} - 900 \, a_{0}^{2} a_{2}^{3} a_{4} + 160 \, a_{0} a_{1}^{3} a_{3} a_{4} - 2050 \, a_{0}^{2} a_{1} a_{2} a_{3} a_{4} + 2250 \, a_{0}^{3} a_{3}^{2} a_{4} - 50 \, a_{0}^{2} a_{1}^{2} a_{4}^{2} + 2000 \, a_{0}^{3} a_{2} a_{4}^{2} + 256 \, a_{1}^{5} - 1600 \, a_{0} a_{1}^{3} a_{2} + 2250 \, a_{0}^{2} a_{1} a_{2}^{2} + 2000 \, a_{0}^{2} a_{1}^{2} a_{3} - 3750 \, a_{0}^{3} a_{2} a_{3} - 2500 \, a_{0}^{3} a_{1} a_{4} + 3125 \, a_{0}^{4}$$

Infine si costruisce la funzione $\operatorname{evaluate\_delta}$ che, dato in ingresso un polinomio $f(x)$ di $n$-esimo grado, restituisce il valore di $\Delta(f)$, sostituendo agli $a_i$ generici di $\operatorname{delta}(n)$ i valori dei coefficienti di $f$.

def evaluate_delta(f):
    d = delta(n)
    
    for i, a in enumerate(f.list()):
        d = eval(f"d.substitute(a_{i}=a)")
        
    return d

Si verifica adesso che per un polinomio con radici multiple il valore di $\operatorname{evaluate\_delta}$ è precisamente zero:

f = (x-2)^2*(x-3)*(x-4)*(x-5)
f

$$\newcommand{\Bold}[1]{\mathbf{#1}}{\left(x - 2\right)}^{2} {\left(x - 3\right)} {\left(x - 4\right)} {\left(x - 5\right)}$$

f.expand()

$$\newcommand{\Bold}[1]{\mathbf{#1}}x^{5} - 16 \, x^{4} + 99 \, x^{3} - 296 \, x^{2} + 428 \, x - 240$$

evaluate_delta(f)

$$\newcommand{\Bold}[1]{\mathbf{#1}}0$$

Infine, si fa lo stesso con un polinomio con radici distinte, per verificare che $\operatorname{evaluate\_delta}$ è strettamente diverso da zero.

g = (x+2)*(x+I)*(x-I)*(x-1)*(x-2)
g

$$\newcommand{\Bold}[1]{\mathbf{#1}}{\left(x + 2\right)} {\left(x + i\right)} {\left(x - i\right)} {\left(x - 1\right)} {\left(x - 2\right)}$$

g.expand()

$$\newcommand{\Bold}[1]{\mathbf{#1}}x^{5} - x^{4} - 3 \, x^{3} + 3 \, x^{2} - 4 \, x + 4$$

evaluate_delta(g)

$$\newcommand{\Bold}[1]{\mathbf{#1}}-1440000$$

---

(c) 2022, ~videtta