Exercices

3. Exercices🔗

Note

Les exercices sont de difficulté variable, de ★ (simple) à ★★★ (complexe).

Quelques exercices sont également disponibles sous forme de Jupyter notebooks:

3.1. Introduction🔗

3.1.1. Intégration: méthode des rectangles ★🔗

La méthode des rectangles permet d'approximer numériquement l'intégrale d'une fonction f:

\[\int_a^b f(x)\,\mathrm{d}x \approx h \sum_{i=0}^{n-1} f(x_{i}) \quad\text{avec}\quad h = (b-a)/n \quad\text{et}\quad x_i = a + (i+1/2)h.\]

On définit la fonction sq renvoyant le carré d'un nombre par (cf. Fonctions):

def sq(x) :
    return x**2

Écrire un programme calculant l'intégrale de cette fonction entre a=0 et b=1, en utilisant une subdivision en n=100 pas dans un premier temps. En comparant au résultat analytique, quelle est la précision de la méthode, et comment dépend-elle du nombre de pas?

3.1.2. Jeu bête: Fizz Buzz ★🔗

Écrire un programme jouant au Fizz Buzz jusqu'à 99:

1 2 Fizz! 4 Buzz! Fizz! 7 8 Fizz! Buzz! 11 Fizz! 13 14 Fizz Buzz! 16...

Astuce

utiliser print(chaine, end=" ") pour éviter le retour à la ligne entre chaque affichage.

3.1.3. PGCD: algorithme d'Euclide ★★🔗

Algorithme d'Euclide.

Écrire un programme calculant le PGCD de deux nombres (p.ex. 756 et 306) par l'algorithme d'Euclide.

3.1.4. Tables de multiplication ★🔗

Écrire un programme affichant les tables de multiplication:

1 x 1 = 1
1 x 2 = 2
...
9 x 9 = 81

3.2. Manipulation de listes🔗

3.2.1. Triangle de Pascal ★🔗

Calculer le Triangle de Pascal d'ordre 5:

[
    [1],
    [1, 1],
    [1, 2, 1],
    [1, 3, 3, 1],
    [1, 4, 6, 4, 1],
    [1, 5, 10, 10, 5, 1]
]

3.2.2. Crible d'Ératosthène ★★🔗

Implémenter le Crible_d'Ératosthène pour afficher les nombres premiers compris entre 1 et un entier fixe, p.ex.:

Liste des entiers premiers <= 41
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41]

Attention

l'implémentation ne doit contenir aucune division ou opération modulo %!

3.2.3. Carré magique ★★★🔗

Un carré magique d’ordre n est un tableau carré n × n dans lequel on écrit une et une seule fois les nombres entiers de 1 à , de sorte que la somme des n nombres de chaque ligne, colonne ou diagonale principale soit constante.

Table 3.1 Carré magique d'ordre 5, où toutes les sommes sont égales à 65.🔗

11

18

25

2

9

10

12

19

21

3

4

6

13

20

22

23

5

7

14

16

17

24

1

8

15

Pour les carrés magiques d’ordre impair, on dispose de l’algorithme suivant – (i,j) désignant la case de la ligne i, colonne j du carré; on se place en outre dans une indexation « naturelle » commençant à 1:

  1. la case (n,(n+1)/2) contient 1 ;

  2. si la case (i,j) contient la valeur k, alors on place la valeur k+1 dans la case (i+1,j+1) si cette case est vide, ou dans la case (i-1,j) sinon. On respecte la règle selon laquelle un indice supérieur à n est ramené à 1.

Programmer cet algorithme pour pouvoir construire un carré magique d’ordre impair quelconque.

3.3. Programmation🔗

3.3.1. Suite de Syracuse (fonction) ★🔗

Écrire une fonction suite_syracuse(n) retournant la (partie non-triviale de la) suite de Syracuse pour un entier n. Écrire une fonction temps_syracuse(n, altitude=False) retournant le temps de vol (éventuellement en altitude) correspondant à l'entier n. Tester ces fonctions sur n=15:

>>> suite_syracuse(15)
[15, 46, 23, 70, 35, 106, 53, 160, 80, 40, 20, 10, 5, 16, 8, 4, 2, 1]
>>> temps_syracuse(15)
17
>>> temps_syracuse(15, altitude=True)
11

3.3.2. Triangles de Pascal et Sierpińksi ★🔗

L'objectif est de calculer et afficher les Triangle de Pascal et Triangle de Sierpińksi, p.ex.:

Pascal's triangle of order 5:
            1
          1   1
        1   2   1
      1   3   3   1
    1   4   6   4   1
  1   5  10  10   5   1

Sierpinski's triangle of order 11:
     *
     **
    * *
    ****
   *   *
   **  **
  * * * *
  ********
 *       *
 **      **
* *     * *
****    ****
  1. Écrire une fonction next_row(row) retournant le rang suivant à partir du rang précédent (liste).

    >>> next_row([1, 1])
    [1, 2, 1]
    
  2. Écrire une fonction compute_pascal(n=5) retournant le triangle de Pascal d'ordre n (liste de listes)

    >>> compute_pascal(3)
    [[1], [1, 1], [1, 2, 1], [1, 3, 3, 1]]
    
  3. Écrire une fonction print_pascal(triangle) affichant un triangle de Pascal comme ci-dessus (utiliser str.center() pour centrer une chaîne de caractère).

  4. Écrire une fonction compute_sierpinski(n=5) calculant un triangle de Sierpiński d'ordre n à partir d'un triangle de Pascal de même ordre, en remplaçant les nombres impairs par 1, et les nombres pairs par 0.

  5. Écrire une fonction print_sierpinski(triangle) affichant un triangle de Sierpiński comme ci-dessus.

  6. BONUS: Écrire une fonction main() lisant un entier positif sur la ligne de commande (utiliser sys.argv), et affichant les triangles de Pascal et de Sierpiński de cet ordre. Appeler cette fonction au sein d'un bloc principal dans le script pascal.py afin de pouvoir l'exécuter sur la ligne de commande:

    $ python pascal.py 5
    Pascal's triangle of order 5:
            1
          1  1
         1  2  1
       1  3  3  1
      1  4  6  4  1
    1  5  10 10 5  1
    Sierpinski's triangle of order 5:
      *
      **
     * *
     ****
    *   *
    **  **
    

3.3.3. Flocon de Koch (programmation récursive) ★★★🔗

En utilisant les commandes left, right et forward de la bibliothèque graphique standard turtle dans une fonction récursive, générer à l'écran un flocon de Koch d'ordre arbitraire.

Flocon de Koch d'ordre 3.

Fig. 3.8 Flocon de Koch d'ordre 3.🔗

3.3.4. Jeu du plus ou moins (exceptions) ★🔗

Écrire un jeu de « plus ou moins »:

Vous devez deviner un nombre entre 1 et 100.
Votre proposition: 27
C'est plus.
[...]
Vous avez trouvé en 6 coups!

La solution sera générée aléatoirement par la fonction random.randint(). Le programme devra être robuste aux entrées invalides (« toto », 120, etc.), et aux lâches abandons par interruption (KeyboardInterrupt ou EOFError).

Dans le même genre: coder un mastermind (utiliser p.ex. la bibliothèque externe rich).

3.4. Programmation Orientée Objet🔗

3.4.1. Animaux (POO/TDD) ★★🔗

Téléchargez le fichier animaux.py et exécutez les tests prédéfinis (dans la seconde partie du fichier) via la ligne de commande (à éxecuter dans un terminal système, après installation de la librairie externe pytest):

$ pytest -v animaux.py

Dans un premier temps, les tests échouent, puisque le proto-code (dans la première partie du fichier) n'est pas encore correct. L'exercice consiste donc à modifier progressivement les classes Animal et Chien pour qu'elles passent avec succès tous les tests:

$ pytest -v animaux.py
======================== test session starts =========================
[...]
collected 8 items
animauxSol.py::test_empty_init PASSED                          [ 12%]
animauxSol.py::test_wrong_init PASSED                          [ 25%]
animauxSol.py::test_init PASSED                                [ 37%]
animauxSol.py::test_str PASSED                                 [ 50%]
animauxSol.py::test_mort PASSED                                [ 62%]
animauxSol.py::test_lt PASSED                                  [ 75%]
animauxSol.py::test_mange PASSED                               [ 87%]
animauxSol.py::test_init_chien PASSED                          [100%]
========================= 8 passed in 0.03s ==========================

C'est le principe du Test Driven Development (voir Développement piloté par les tests).

3.4.2. Jeu de la vie (POO) ★★★🔗

On se propose de programmer l'automate cellulaire le plus célèbre, le Jeu de la vie.

Pour cela, vous créerez une classe Life qui contiendra la grille du jeu ainsi que les méthodes qui permettront son évolution. Vous initialiserez la grille aléatoirement à l'aide de la fonction random.choice(), et vous afficherez l'évolution de l'automate dans la sortie standard du terminal, p.ex.:

...#..#.....##.......
.....###.............
#........#...........
.....#...#...........
................##...
.....#.#......##..#..
..............##.##..
..............##.##..
................#....

Astuce

Pour que l'affichage soit agréable à l'oeil, vous marquerez des pauses entre l'affichage de chaque itération grâce à la fonction time.sleep().

Note

une version plus détaillée de cet exercice est disponible dans le TP5: jeu de la vie (POO).

Dans le même genre: Abelian sandpile model.

3.5. Manipulation de tableaux Numpy🔗

3.5.1. Manipulation de tableaux ★🔗

Écrire une fonction checkboard(n=8) générant un tableau $n×n$ en damier:

[[0, 1, 0, 1, ...],
 [1, 0, 1, 0, ...],
 [0, 1, 0, 1, ...],
 ...
]

Généraliser pour choisir la taille $m$ des cases, p.ex. checkboard(n=8, m=2):

[[1, 1, 0, 0, 1, 1, 0, 0],
 [1, 1, 0, 0, 1, 1, 0, 0],
 [0, 0, 1, 1, 0, 0, 1, 1],
 [0, 0, 1, 1, 0, 0, 1, 1],
 ...
]

Généraliser à des tableaux de rang arbitraire: checkboard(n=8, m=1, rank=2).

3.5.2. Inversion de matrice ★🔗

Créer un tableau carré réel \(\mathsf{r}\) aléatoire (numpy.random.randn()), calculer la matrice hermitienne \(\mathsf{m} = \mathsf{r} \cdot \mathsf{r}^T\) (numpy.dot() ou numpy.linalg.matmul()), l'inverser (numpy.linalg.inv()), et vérifier que \(\mathsf{m} \cdot \mathsf{m}^{-1} = \mathsf{m}^{-1} \cdot \mathsf{m} = \mathsf{1}\) (numpy.eye()) à la précision numérique près (numpy.allclose()).

3.5.3. Median Absolute Deviation🔗

En statistique, le Median Absolute Deviation (MAD) est un estimateur robuste de la dispersion d'un échantillon 1D: MAD = median(| x - median(x) |).

À l'aide des fonctions numpy.median() et numpy.abs(), écrire une fonction mad(x, axis=None) calculant le MAD d'un tableau, éventuellement le long d'un ou plusieurs de ses axes.

>>> a = np.arange(4 * 5).reshape(4, 5).astype(float) ** 2
>>> mad(a)
80.0
>>> mad(a, axis=0)
array([50., 60., 70., 80., 90.])
>>> mad(a, axis=1)
array([ 4. 15. 25. 35.])
>>> mad(a, axis=(0, 1))
80.0

3.5.4. Distribution du pull ★★★🔗

Le pull est une quantité statistique permettant d'évaluer la conformité des erreurs par rapport à une distribution de valeurs (typiquement les résidus d'un ajustement). Pour un échantillon \(\mathbf{x} = [x_i]\) et les erreurs associées \(\mathrm{d}\mathbf{x} = [\sigma_i]\), le pull est défini par:

  • moyenne optimale (pondérée par la variance): \(E = (\sum_{i} x_i/\sigma_i^2)/(\sum_i 1/\sigma_i^2)\);

  • erreur sur la moyenne pondérée: \(\sigma_E^2 = 1/\sum(1/\sigma_i^2)\);

  • définition du pull: \(p_i = (x_i - E_i)/(\sigma_{E_i}^2 + \sigma_i^2)^{1/2}\), où \(E_i\) et \(\sigma_{E_i}\) sont calculées sans le point i.

Si les erreurs \(\sigma_i\) sont correctes, la distribution du pull est centrée sur 0 avec une déviation standard de 1.

Écrire une fonction pull(x, dx) calculant le pull de tableaux 1D.

>>> pull(np.arange(3), np.arange(3) + 1)
array([-0.67356452,  0.36140316,  0.57498891])

3.6. Méthodes numériques🔗

3.6.1. Quadrature et zéro d'une fonction ★🔗

À l'aide des algorithmes disponibles dans scipy:

  • calculer numériquement l'intégrale \(\int_0^\infty \frac{x^3}{e^x-1}\mathrm{d}x = \pi^4/15\);

  • résoudre numériquement l'équation \(x\,e^x = 5(e^x - 1)\).

3.6.2. Schéma de Romberg ★★🔗

Écrire une fonction integ_romberg(f, a, b, epsilon=1e-6) permettant de calculer l'intégrale numérique de la fonction f entre les bornes a et b avec une précision epsilon selon la Méthode_de_Romberg.

Tester sur des solutions analytiques et en comparant à scipy.integrate.romberg().

3.6.3. Méthodes de Runge-Kutta ★★🔗

Développer des fonctions ode_rk1(f, y0, t), ode_rk2(f, y0, t) et ode_rk4(f, y0, t) permettant d'intégrer numériquement l'équation différentielle du 1er ordre \(y' = f(t, y)\) avec \(y(t_0) = y_0\), en utilisant les Méthodes_de_Runge-Kutta d'ordre 1, 2 et 4.

Tester sur le cas analytique:

\[\begin{split}\begin{cases} f(t, y) = \sin^2(t) \times y(t) \\ y(t=0) = 1 \end{cases} \Leftrightarrow y(t) = e^{(t - \sin(t)\cos(t))/2}\end{split}\]

et en comparant à scipy.integrate.odeint().

3.7. Visualisation (matplotlib)🔗

3.7.1. Quartet d'Anscombe ★★🔗

Après chargement des données, calculer et afficher les propriétés statistiques des quatres jeux de données du Quartet d'Anscombe:

Table 3.2 Quartet d'Anscombe🔗

I

II

III

IV

x

y

x

y

x

y

x

y

10.0

8.04

10.0

9.14

10.0

7.46

8.0

6.58

8.0

6.95

8.0

8.14

8.0

6.77

8.0

5.76

13.0

7.58

13.0

8.74

13.0

12.74

8.0

7.71

9.0

8.81

9.0

8.77

9.0

7.11

8.0

8.84

11.0

8.33

11.0

9.26

11.0

7.81

8.0

8.47

14.0

9.96

14.0

8.10

14.0

8.84

8.0

7.04

6.0

7.24

6.0

6.13

6.0

6.08

8.0

5.25

4.0

4.26

4.0

3.10

4.0

5.39

19.0

12.50

12.0

10.84

12.0

9.13

12.0

8.15

8.0

5.56

7.0

4.82

7.0

7.26

7.0

6.42

8.0

7.91

5.0

5.68

5.0

4.74

5.0

5.73

8.0

6.89

Pour chacun des jeux de données, tracer y en fonction de x, ainsi que la droite de régression linéaire.

3.7.2. Diagramme de bifurcation: la suite logistique ★★🔗

Écrivez une fonction qui calcule la valeur d'équilibre de la Suite logistique pour un \(x_0\) (nécessairement compris entre 0 et 1) et un paramètre \(r\) (parfois noté \(\mu\)) donné.

Générez l'ensemble de ces points d'équilibre pour des valeurs de \(r\) comprises entre 0 et 4:

Diagramme de bifurcation

Fig. 3.9 Diagramme de bifurcation.🔗

N.B. Vous utiliserez la bibliothèque Matplotlib pour tracer vos résultats.

3.7.3. Ensemble de Julia ★★🔗

Représentez l'ensemble de Julia pour la constante complexe \(c = 0.284 + 0.0122j\):

Ensemble de Julia

Fig. 3.10 Ensemble de Julia pour \(c = 0.284 + 0.0122j\).🔗

On utilisera la fonction numpy.meshgrid() pour construire le plan complexe, et l'on affichera le résultat grâce à la fonction matplotlib.pyplot.imshow().

Voir également: Superposition d'ensembles de Julia

3.8. Packaging🔗

3.8.1. Partie 1 ★🔗

  1. Créer un package erathostene_pkg constitué du seul module erathostene_mod.py contenant la fonction erathostene (voir ci-dessus), avec l’arborescence suivante:

    Erathostene_repo/            # Repo
    ├── erathostene_pkg/         # Package
    │   ├── erathostene_mod.py   # Module
    │   └── __init__.py          # Initialisation
    ├── LICENSE                  # Licence
    ├── README                   # Texte libre
    └── pyproject.toml           # Configuration
    
  2. Ajouter un __main__.py à votre package pour pouvoir exécuter:

    $ python -m erathostene_pkg 50
    [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]
    

    Utiliser argparse (ou sys.argv pour ce cas simple) pour la gestion de l'argument sur la ligne de commande.

  3. Ajouter un point d'entrée (entry point) à votre fichier de configuration pour pouvoir exécuter:

    $ erathostene_pgm 50
    [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]
    

3.8.2. Intermède Git ★🔗

Héberger votre package sur un dépôt dédié sur le serveur https://gitlab.in2p3.fr. Faire tester l'installation et l'utilisation de votre package par un·e collègue, en lui transmettant l'adresse du dépôt.

3.8.3. Partie 2 ★🔗

Prendre exemple sur pyyc.

3.8.3.1. Documentation🔗

  1. Documenter la fonction erathostene, le module erathostene_mod et le package erathostene_pkg par des doc-strings.

  2. Ajouter un répertoire docs/ et y configurer et générer la documentation sphinx du package (sphinx-quickstart et sphinx-apidoc).

  3. Inclure dans une section dédiée de la documentation un notebook d'exemple (ajouté sous docs/notebooks/).

3.8.3.2. Tests🔗

  1. Inclure un doc-test dans la documentation de la fonction erathostene.

  2. Ajouter un répertoire tests/ et y ajouter un test unitaire, en utilisant pytest ou unittest.

  3. Évaluer la couverture des tests avec coverage.

3.8.3.3. Intégration continue🔗

  1. Ajouter un fichier de configuration d'intégration continue .gitlab-ci.yml pour une construction et une mise en ligne automatique de la documentation.

3.9. Manipulation de données (pandas)🔗

3.9.1. Visualisation ★🔗

À partir du jeu de données en ligne Data on CO2, générer une figure du genre:

Émission de CO2 par personne et cumulée.

Fig. 3.11 Évolution de l'émission de CO2 par personne et cumulée.🔗

Source: Our World In Data

3.9.2. Statistiques ★★🔗

  1. Chargez l'ensemble de données Fifa World Cup 2018.

  2. Combien de matchs et combien de pays ?

  3. Calculez le nombre moyen de buts par match. Même question pour les passes complétées et les distances parcourues.

  4. Trouvez l'équipe ayant couvert la plus grande distance dans un match.

  5. Calculez le score moyen.

  6. Trouvez les équipes avec la plus grande (resp. la plus petite) distance parcourue par match.

  7. Trouvez les matchs où les joueurs ont couru le plus (resp. le moins) au total.

Idées supplémentaires:

3.10. Exercices en vrac🔗

3.10.1. Équation différentielle: boulet de canon ★★🔗

À l'aide de la fonction scipy.integrate.odeint(), intégrer les équations du mouvement d'un boulet de canon soumis à des forces de frottement « turbulentes » (en \(v^2\)):

\[\ddot{\mathbf{r}} = \mathbf{g} - \frac{\alpha}{m}v\times\mathbf{v}.\]

Utiliser les valeurs numériques pour un boulet de canon de 36 livres:

g = 9.81        # Pesanteur [m/s2]
cx = 0.45       # Coefficient de frottement d'une sphère
rhoAir = 1.2    # Masse volumique de l'air [kg/m3]
rad = 0.1748/2  # Rayon du boulet [m]
rho = 6.23e3    # Masse volumique du boulet [kg/m3]
mass = 4./3.*np.pi*rad**3 * rho            # Masse du boulet [kg]
alpha = 0.5*cx*rhoAir*np.pi*rad**2 / mass  # Coeff. de frottement / masse
v0 = 450.       # Vitesse initiale [m/s]
alt = 45.       # Inclinaison du canon [deg]

Voir également: Équations_de_prédation_de_Lotka-Volterra.

3.10.2. Redshift de la galaxie hôte de SN 2005eu (astropy) ★★★🔗

L'objectif de l'exercice est 1. d'extraire le spectre de la galaxie hôte après soustraction du spectre du ciel, 2. de déterminer le redshift à partir des raies en émission.

Extraction du spectre

  1. Récupérez et ouvrez le cube \((x, y, \lambda)\) de la galaxie hôte de SN 2005eu (C06_326_085_006_17_R.fits) observée avec le spectrographe à champ intégral SNIFS.

    Il s'agit d'un fichier FITS avec NAXIS=3 (15×15 spaxels × 1402 longueurs d'onde), avec le signal en 1e extension et la variance en 2e:

    >>> hdu = astropy.io.fits.open("C06_326_085_006_17_R.fits")
    >>> signal = hdu[0].data * 1e16
    >>> variance = hdu[1].data * 1e32
    >>> y, x = np.mgrid[-7:+7.1, -7:+7.1]  # Coordonnées spatiales
    >>> hdr = hdu[0].header
    >>> wavelengths = hdr['CRVAL3'] + hdr['CDELT3'] * np.arange(hdr['NAXIS3'])
    
  2. Calculer et tracer le spectre intégré (et l'erreur associée) sur les deux dimensions spatiales.

  3. Construire et afficher l'image intégrée sur la dimension spectrale.

  4. Estimer le spectre du ciel nocturne à partir des zones externes du cube:

    ciel = mean(cube[r > r_out])
    

    Dans un premier temps, vous pouvez déterminer le centre et le rayon de l'ouverture de visu.

  5. Estimer le spectre de la galaxie à partir de la zone interne du cube auquel vous avez soustrait la contribution du ciel:

    galaxy = sum(cube[r < r_in])
    

Détermination du redshift

Note

Si vous n'avez pas réussi la partie précédente, vous pouvez utiliser le jeu de données SN2005eu_host.npz:

>>> archive = np.load("SN2005eu_host.npz")
>>> λ = archive['wavelength']
>>> signal = archive['flux']
>>> variance = archive['variance']

Estimer son redshift (et l'erreur associée) d'abord avec un simple ajustement gaussien de la raie Hα, puis à partir d'un ajustement conjoint des 5 raies en émission [Hα], [NIIa,b] et [SIIa,b]. Les longueurs d'onde de référence (dans l'air) sont disponibles en ligne.