Typage et précision en Fortran

Un problème récurrent, mais peu malheureusement trop ignoré des développeur⋅se⋅s Fortran, est le typage des flottants. Le⋅a développeur⋅se, un peu feignant⋅e, et faisant un peu trop confiance à son compilateur, aura tendance à faire des déclarations de ce type :

real(kind=8) :: x
x = 0.1

Iel pensera sûrement que tout se passera bien, et que x étant défini en double précision, c'est bien la précision qui sera obtenue après l'affectation précédente. Est-ce que c'est vraiment ce qui se passe ? Un moyen simple de le contrôler est d'écrire le programme suivant :

real(kind=8) :: x, y
x = 0.1; y = 0.1d0

print '("x = ", F20.18, ", y = ", F20.18)', x, y

Après compilation et exécution, voilà le résultat obtenu :

x = 0.10000000149011612, y = 0.100000000000000001

y a donc bien la précision attendu, par contre x est en fait en simple précision ! Des programmes entiers, par cette simple erreur, peuvent se retrouver corrompu et être, de fait, en simple précision et non en double.

Mais que se passe-t-il à l'affectation, pour qu'on obtienne ce résultat ? La norme Fortran stipule que les évaluations se font de droite à gauche, avec retypage éventuel si besoin. En clair, ça veut dire que le compilateur, arrivant à l'expression x = 0.1, va d'abord évaluer 0.1, qui est en simple précision. Il va ensuite l'affecter à x (x se retrouve ainsi égal à 0.1 en simple précision), puis, x étant déclaré comme étant en double précision, il fait un retypage de simple à double précision. Mais la promotion en double précision ne permet pas d'attendre effectivement la double précision, puisque ça n'était pas la précision originelle de 0.1. x se retrouve ainsi promu en double précision, mais en étant « complété » par des scories.

Il existe plusieurs possibilités pour contourner ce problème. Le premier est d'utiliser les options de compilations idoines pour s'assurer de la précision de ses flottants (-fdefault-real-8 avec gfortran, par exemple). Cette solution peut poser problème, les flottants en double précision étant alors promus en quadruple precision ! C'est possible de les forcer en double précision en utilisant d'autres options de compilation (-fdefault-double-8 avec gfortran), mais ça n'est pas toujours extrêmement simple à gérer.

Une solution est de relire attentivement son code, et se s'assurer que toutes les affectations sont en accord avec la norme... Ce qui peut être extrêmement fastidieux, et qui n'est pas non plus exempt d'erreurs.

La dernière solution est d'utiliser PCF, Precision Converter for Fortran, un script Python que j'ai écrit et qui permet d'harmoniser le typage des flottants dans l'ensemble d'un code Fortran. Il suffit de lui précisier les sources à convertir, et le prefixe désiré (d0, q0, _dp) pour qu'il se charge d'harmoniser l'ensemble des sources.