Comment construire les polynômes de Taylor à n'importe quel ordre de fonctions usuelles sans calcul formel sur Haskell ?
Pour construire $T_{n,0}(x)= \frac{x^0}{0!}f^{(0)}(0) + \frac{x^1}{1!}f^{(1)}(0) + \frac{x^2}{2!}f^{(2)}(0) + \cdots $ on va créer deux listes infinies : $[\frac{x^0}{0!}, \frac{x^1}{1!} , \frac{x^2}{2!} ,...]$ et $[f^{(0)}(0), f^{(1)}(0), f^{(2)}(0),....]$.
La première liste est facile à créer avec scanl qui renvoie la liste des différents états de l'accumulateur d'un pliage:
Attention ! C'est une liste infinie...
*AnalyseNumerique> take 5 (vecTaylor 1) [1.0,1.0,0.5,0.16666666666666666,4.1666666666666664e-2]
Pour la deuxième liste, sans calcul formel, on va créer des listes pour chaque fonction élémentaire en cherchant une relation de récurrence.
Par exemple, pour l'exponentielle, c'est [1,1,1,1,....], pour $x \mapsto \frac{1}{1-x}$ c'est [1,1,2!,3!,4!,...], pour le sinus c'est [0,1,0,-1,0,1,0,-1,...].
On commence par donner des alias aux fonctions:
Ensuite, on crée un couple (fonction, liste de ses dérivées successives en 0) :
{- On renvoie le couple (fonctions, liste de ses dérivées successives en 0) On cherche le plus souvent une formule de récurrence et on utilise scanl Sinon, on répète des motifs tTruc :: (Double -> Double, [Double]) -} tExp = (exp, repeat 1.0) tInvMoins = (invMoins, scanl (\ acc k -> acc * k ) 1.0 [1 .. ]) tInvPlus = (invPlus, scanl (\ acc k -> acc * (-k)) 1.0 [1 .. ]) tRacine = (racinePlus, scanl (\ acc k -> acc * (-0.5) * (2*k - 1)) 1.0 [0 .. ]) tLogPlus = (logPlus, 0.0 : (snd tInvPlus)) tSin = (sin, cycle [0.0,1.0,0.0,-1.0]) tCos = (cos, cycle [1.0,0.0,-1.0,0.0])
Ce sont encore des listes infinies...
Par exemple pour $x \mapsto \frac{1}{1-x}$:
Il ne reste plus qu'à effectuer le produit scalaire des deux listes infinies avec sum et zipWith :
Pour l'afficher, on peut utiliser la simplissime bibliothèque EasyPlot qui permet d'interfacer GNUPlot sous Haskell.
On n'oubliera pas d'importer EasyPlot
import Graphics.EasyPlot
On peut alors utiliser plot X11 (Function2D [Title ..., Color ... ] [Range ... ..., Step ...] fonction) :
-- On trace sur un même graphique la fonction et son poly de Taylor sur [a,b] avec un pas de 1/2^8 compareDL :: Double -> Double -> ((Double -> Double) , [Double]) -> Int -> IO Bool compareDL a b tf n = plot X11 [ Function2D [Title ("n = " ++ (show n)) , Color Blue ] [Range a b, Step (2**(-8))] (polyTaylor (snd tf) n), Function2D [Title "Vraie fonction", Color Red ] [Range a b, Step (2**(-8))] (fst tf) ]
Ensuite, on demande l'affichage dans un terminal GNUplot avec l'option X11 (ou Windows ou Aqua selon les OS) où un export en PDF, PNG, JPG, etc.
*AnalyseNumerique> compareDL (-0.8) 1.2 tLogPlus 4
Renvoie, pour ${\rm ln}(1+x)$:
On peut en profiter pour créer une petite animation. Pour rester libre, on évitera Flash, GIF et cie....
Voici une petite animation PNG. On exporte les images en PNG et on utilise par exemple apngasm.
On commence par modifier compareDL pour pallier au problème de l'absence de contrôle sur l'échelle des y:
animIniDL :: Double -> Double -> ((Double -> Double) , [Double]) -> Int -> IO Bool animIniDL a b tf n = plot (PNG ("animIndi" ++ (show (n + 10)) ++ ".png")) ( [ Function2D [Title " " , Color White ] [Range a b, Step (2**(-8))] -- pour pallier l'abscence de contrôle sur l'échelle des y (polyTaylor (snd tf) k) | k <- [0..10] ] ++ [ Function2D [Title ("n = " ++ (show n)) , Color Blue ] [Range a b, Step (2**(-8))] (polyTaylor (snd tf) n), Function2D [Title "Vraie fonction", Color Red ] [Range a b, Step (2**(-8))] (fst tf) ] )
Ensuite on crée une liste d'objets graphiques :
Haskell doit communiquer avec le monde extérieur via ici la monade IO. On utilise l'intruction do qui permet d'exécuter des commandes IO séquentiellement : cela donne l'illusion de faire de la programmation impérative...
Il ne reste plus qu'à créer une fonction qui crée la suite de figures de manière transparente :
Par exemple pour les 11 premiers polynômes de Taylor de $x\mapsto {\rm ln} (1+x)$ sur $[-0.8,1.3]$:
anim (-0.8) 1.3 tLogPlus 10
Il reste, dans un shell, à regrouper les fichiers PNG
$ apngasm animLog.png animIndi*.png 1 1
pour créer cette animation:
