{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Déterminant, cofacteurs, etc.\n",
    "\n",
    "Marc Lorenzi\n",
    "\n",
    "2 avril 2023"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import random\n",
    "import time\n",
    "compteur = 0"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il est bien connu que l'algorithme de calcul d'un déterminant par un développement selon une ligne ou une colonne est hautement inefficace. Cela n'empêche, c'est un excellent exercice de programmation. Et puis sait-on jamais ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. La formule de Lagrange"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.1 Introduction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Soit $A\\in\\mathcal M_n(\\mathbb K)$. Pour tous $i,j\\in[|0,n-1|]$, notons $A[i|j]$ la matrice obtenue en supprimant la ligne $i$ et la colonne $j$ de $A$. Remarquons que $A[i|j]\\in\\mathcal M_{n-1}(\\mathbb K)$.\n",
    "\n",
    "- Le mineur de $A$ ligne $i$, colonne $j$ est \n",
    "\n",
    "$${\\rm Min}(A,i,j)=\\det A[i|j]$$\n",
    "\n",
    "- Le cofacteur de $A$ ligne $i$, colonne $j$ est \n",
    "\n",
    "$${\\rm Cof}(A,i,j)=(-1)^{i+j}{\\rm Min}(A,i,j)$$\n",
    "\n",
    "La *formule de Lagrange* nous dit que pour tout $i\\in[|0,n-1|]$,\n",
    "\n",
    "$$\\det A=\\sum_{j=0}^{n-1}A_{ij}{\\rm Cof}(A,i,j)$$\n",
    "\n",
    "On a aussi pour tout $j\\in[|0,n-1|]$,\n",
    "\n",
    "$$\\det A=\\sum_{i=0}^{n-1}A_{ij}{\\rm Cof}(A,i,j)$$\n",
    "\n",
    "Ces égalités s'appellent aussi les formules de *développement* de $\\det A$ par rapport à la ligne $i$ ou la colonne $j$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.2 Afficher une matrice"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La fonction `print_mat` affiche de façon agréable une matrice."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def print_mat(A):\n",
    "    p = len(A)\n",
    "    q = len(A[0])\n",
    "    s = q * '+-----' + '+'\n",
    "    for i in range(p):\n",
    "        print(s)\n",
    "        for j in range(q):\n",
    "            print('|%4d ' % A[i][j], end='')\n",
    "        print('|')\n",
    "    print(s)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print_mat([[1, 2, 3],[4, 5, 6],[7, 8, 9]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.3 Mineurs, cofacteurs, comatrice"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La fonction `mineur` prend en paramètres une matrice carrée $A$ et deux entiers $i$ et $j$. Elle renvoie ${\\rm Min}(A,i,j)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def mineur(A, i, j):\n",
    "    A1 = supprimer_ligne_colonne(A, i, j)\n",
    "    return determinant(A1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le cofacteur de $A$ ligne $i$, colonne $j$ est égal à $(-1)^{i+j}$ fois le mineur correspondant."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def cofacteur(A, i, j):\n",
    "    m = mineur(A, i, j)\n",
    "    if (i + j) % 2 == 0: return m\n",
    "    else: return -m"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Tant que nous y sommes, la comatrice de $A$ est la matrice de ses cofacteurs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def comatrice(A):\n",
    "    n = len(A)\n",
    "    return [[cofacteur(A, i, j) for i in range(n)] for j in range(n)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Enfin, la transcomatrice de $A$, est la transposée de sa comatrice."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def transposee(A):\n",
    "    m = len(A)\n",
    "    n = len(A[0])\n",
    "    return [[A[j][i] for j in range(m)] for i in range(n)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print_mat(transposee([[1] , [2], [3]]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print_mat(transposee([[1, 2, 3] , [4, 5, 6], [7, 8, 9]]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def tilde(A):\n",
    "    return transposee(comatrice(A))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mais ça veut dire quoi, tilde ? Un peu de culture !\n",
    "\n",
    "__Étymol. et Hist. 1834 orth. esp. (Boiste). Mot esp. att. dep. 1433 (E. de Villena ds Cor.-Pasc.), issu du lat. titulus (titre*) qui signifiait propr. « ce qui désigne, signale ». Cf. a. fr. title, fr. mod. ti(l)tre « signe abréviatif » (v. titre, étymol. 1 b; FEW t. 13, 1, p. 360a).__"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ah, maintenant on comprend mieux."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.4 Déterminant"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le déterminant de $A$ peut être calculé par un développement par rapport à la première ligne, le cas d'un déterminant vide étant évident. Remarquons que nous allons écrire une fonction `determinant` qui se prépare à appeler `cofacteur` qui appelle `mineur` qui appelle ... `determinant`. Nous avons là un exemple de fonctions mutuellement récursives. Mais comme ces fonctions s'appellent sur des matrices de tailles strictement décroissantes et que l'on sait traiter le cas de la matrice « vide », nous sommes assurés de la terminaison de l'algorithme."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La variable `compteur` qui apparaît dans le code ci-dessous est une variable globale, qui compte le nombre d'additions et de multiplications effectuées par l'algorithme. C'est très mal de faire cela mais c'est pédagogique. Nous aurons ainsi une idée de la complexité du calcul d'un déterminant par cette mauvaise (?) méthode."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def determinant(A):\n",
    "    global compteur\n",
    "    n = len(A)\n",
    "    if n == 0: return 1\n",
    "    else:\n",
    "        s = 0\n",
    "        for j in range(n):\n",
    "            s = s + A[0][j] * cofacteur(A, 0, j)\n",
    "            compteur = compteur + 2 # 1 addition, 1 multiplication\n",
    "    return s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il reste à écrire la fonction supprimant une ligne et une colonne de la matrice $A$. Une petite fonction auxiliaire `phi` nous évite d'avoir à considérer 4 cas dans la fonction principale. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def phi(p, i):\n",
    "    if p < i: return p\n",
    "    else: return p + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def supprimer_ligne_colonne(A, i, j):\n",
    "    n = len(A)\n",
    "    rg = range(n - 1)\n",
    "    return [[A[phi(p, i)][phi(q, j)] for p in rg] for q in rg]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Testons."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = [[1, 2, 3], [4, 5, 6], [7, 8, 9]]\n",
    "print_mat(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "B = supprimer_ligne_colonne(A, 1, 0)\n",
    "print_mat(B)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(determinant(A))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print_mat(tilde(A))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.4 Matrices « aléatoires »"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour des tests un peu plus conséquents, définissons une fonction prenant en paramètre un entier $n$ et renvoyant une matrice $n\\times n$ dont les coefficients sont des entiers aléatoires entre $-9$ et $9$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def random_matrix(n):\n",
    "    rg = range(n)\n",
    "    return [[random.randint(-9, 9) for i in rg] for j in rg]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour une matrice $3\\times 3$, le calcul du déterminant demande 30 opérations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = random_matrix(3)\n",
    "print_mat(A)\n",
    "compteur = 0\n",
    "print('Déterminant : ', determinant(A))\n",
    "print(\"Nombre d'opérations : \", compteur)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Et pour une matrice $8 \\times 8$ ?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = random_matrix(8)\n",
    "print_mat(A)\n",
    "compteur = 0\n",
    "print('Déterminant : ', determinant(A))\n",
    "print(\"Nombre d'opérations : \", compteur)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il nous faut 219200 opérations. C'est beaucoup."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Factorielle"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.1 Complexité du calcul du déterminant"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Faisons un petit calcul. Calculer un déterminant $0 \\times 0$ est immédiat et demande 0 opération. Pour calculer un déterminant de taille $n$, il faut calculer $n$ déterminants de taille $n-1$, les multiplier par des nombres, puis additionner. Appelons $C_n$ le nombre d'opérations nécessaires au calcul d'un détermant de taille $n$. Nous avons $C_0= 0$ et, pour tout entier $n>0$, $C_n = n C_{n-1} + 2n$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def C(n):\n",
    "    if n == 0: return 0\n",
    "    else: return n * C(n - 1) + 2 * n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for k in range(10):\n",
    "    print(k, C(k))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "On est tout bons, $C_8=219200$. Et $C_n$, il vaut combien ? Le lecteur consciencieux montrera par récurrence sur $n$ que \n",
    "\n",
    "$$C_n = 2n!\\sum_{k=0}^{n-1}\\frac{1}{k!}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Juste pour le plaisir d'écrire du Python ...."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def fact(n):\n",
    "    p = 1\n",
    "    for k in range(1, n + 1): p *= k\n",
    "    return p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "[fact(n) for n in range(11)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def complexite(n):\n",
    "    p = 2 * fact(n)\n",
    "    s = 0\n",
    "    for k in range(n ):\n",
    "        s = s + 1.0 / fact(k)\n",
    "    return p * s"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "complexite(8)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Il est bien connu que \n",
    "\n",
    "$$\\sum_{k=0}^{n-1}\\frac 1 {k!}\\underset{n\\to\\infty}{\\longrightarrow}e$$\n",
    "\n",
    "Ainsi, un équivalent de la complexité du calcul d'un déterminant de taille $n$ est \n",
    "\n",
    "$$C_n\\sim 2en!$$\n",
    "\n",
    "Cela signifie que le temps de calcul d'un déterminant est $T_n\\sim 2Ken!$ où la constante $K$ dépend évidemment de l'ordinateur utilisé. Que vaut $K$ pour ma machine ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2 Tests"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = random_matrix(8)\n",
    "print_mat(A)\n",
    "compteur = 0\n",
    "t1 = time.time()\n",
    "d = determinant(A)\n",
    "t2 = time.time()\n",
    "t_mat8 = t2 - t1\n",
    "print(d)\n",
    "print(compteur)\n",
    "print('%5.3fs' % t_mat8)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Python me dit que pour calculer un déterminant de taille 8 il lui faut environ 0.5 seconde. Ainsi, $2ke8!\\simeq 0.5$, d'où la valeur de $K$ pour ma machine :"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "K = t_mat8 / (2 * fact(8) * 2.71828)\n",
    "print(K)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "C = 2 * K * 2.718281828 "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pour calculer un déterminant de taille $n$, ma machine met environ $1.24\\times 10^{-5}n!$ secondes. Voici donc le temps prévisionnel pour un déterminant de taille 10."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "C * fact(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "On essaye ?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = random_matrix(10)\n",
    "t1 = time.time()\n",
    "d = determinant(A)\n",
    "t2 = time.time()\n",
    "print(d)\n",
    "print('%5.3fs' % (t2 - t1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Et pour un déterminant de taille 17 ? Rappelons qu'il y a 86400 secondes dans une journée et 365 jours par an."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "1.08e-5 * fact(17) / 86400 / 365"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "121 ans. Eh bien non, on ne teste pas :-). Alors tout cela ne sert à rien ? Pas si sûr."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Inutile ?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.1 Calcul formel"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Soient $a,b,c\\in\\mathbb C$. La matrice \n",
    "\n",
    "$$A=\\begin{pmatrix}0&a&b&c\\\\a&0&c&b\\\\b&c&0&a\\\\c&b&a&0\\end{pmatrix}$$\n",
    "\n",
    "est-elle inversible ?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import *\n",
    "a, b, c = symbols('a b c')\n",
    "init_printing()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = [\n",
    "    [0, a, b, c],\n",
    "    [a, 0, c, b],\n",
    "    [b, c, 0, a],\n",
    "    [c, b, a, 0]\n",
    "]\n",
    "Matrix(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "e = determinant(A)\n",
    "e"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "factor(e)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "La matrice est inversible si et seulement si $a\\pm b\\pm c\\ne 0$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Tentons l'exo 412 du TD 820. Et ayons une pensée pour ceux qui n'ont toujours pas installé Jupyter."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "alpha, beta, gamma = symbols('alpha beta gamma')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = [\n",
    "    [alpha - beta - gamma, 2 * alpha, 2 * alpha],\n",
    "    [2 * beta, beta - alpha - gamma, 2 * beta],\n",
    "    [2 * gamma, 2 * gamma, gamma - alpha - beta]\n",
    "]\n",
    "Matrix(A)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "D = determinant(A)\n",
    "D"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "factor(D)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.2 Conclusion"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Tout n'est pas sombre en Terre du Milieu. Pour des matrices de taille raisonnable dont certains coefficients contiennent des paramètres, nous pouvons déterminer si oui ou non ces matrices sont inversibles, grâce à notre algorithme  de calcul de déterminant."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}