Simulateur de galaxie en python

Comment appliquer Newton efficacement en python

Pour ceux qui ne se sont pas rendu compte, j’adore les calculs un peu mathématiques (théorie du chaos, train lumineux, ordinateur quantique), avec comme langage de prédilection le python.

Aujourd’hui, je viens vous parler de mon dernier projet : un simulateur de gravité. Le but du programme est simple. On initialise des masses ponctuelles avec positions et vitesses. Ensuite on lance mon programme qui va calculer la force de Newton entre toutes les masses, et calculer les prochaines positions des astre.

Il peut donc simuler le système Terre-Lune, le système soleil ou même une collision de galaxies comme on peut voir dans la vidéo.

Je commence par vous donner l’adresse du programme sur mon github.

Le programme intéressant est le fichier processor_multi.pyx. Contrairement à mes autres algorithmes un peu rebutants celui-ci est très simple. La force de gravité de Newton qu’exerce la planètei sur la planètej s’écrit :

M1

Ensuite, pour calculer la nouvelle position de l’objet, on reste sur Newton avec sa troisième lois. Dans notre cas de simulation d’astres massifs, cela devient :

M2

Place au code

Pour i allant de 0 à n
    Fx = 0
    Fy = 0
    Pour j allant de 0 à n
        Si i différent de j
            Fx = Fx + G*m[i]*m[j]*rx/(rx²+ry²)^1.5
            Fy = Fy + G*m[i]*m[j]*ry/(rx²+ry²)^1.5
    x[t][i] = 2*x[t-1][i] - x[t-2][i] + Fx/m[i]
    y[t][i] = 2*y[t-1][i] - y[t-2][i] + Fy/m[i]

Et hope le programme et fini ! C’est tout

Ou alors il y a un petit problème. Ce programme fait intervenir deux boucles for, une complexité . Cela ne pose pas de problème pour la lune ou le système solaire, par contre pour simuler des galaxies, le programme va mettre plusieurs jours (22 jours dans mon cas). Il faut donc coder plus malin.

Premièrement, on remarque que la force est symétrique. Si on calcule la force M3 , on peut facilement calculer M4. J’ai mis en place cette technique par les matrices Fx et Fy qui représentent : M5 et M6. Ainsi je calcule que la moitié triangulaire haute des matrices, j’obtiens les matrices entières par M7. Il me reste plus qu’a sommer sur les lignes pour avoir la force externe pour chaque astre, M8.

Deuxième, il va falloir compiler. Pour ce faire j’ai utilisé la bibliothèque Cython qui permet de compiler du python. Elle apporte également un langage proche du python, mais fortement typé où toutes les variables doivent être déclarées au début. C’est ce langage que j’ai utilisé, et qui donne des fichiers *.pyx.

Au final après toutes ces optimisations, le programme de collision de galaxies est passé de 22 jours à 8 heures. C’est donc un gain très significatif 😉

A bientôt et amusez vous bien avec le programme.