Initiation au FORTRAN 77

Pourquoi le FORTRAN

Cette page est ancienne, elle date de 2006 déjà. Pourtant, j'ai décidé de la conserver. Parce que le FORTRAN est toujours utilisé par les scientifiques, et qu'il est encore enseigné par quelques fondus...

Les raisons de cette longévité sont nombreuses :

FORTRAN est un langage structuré, de syntaxe simple et facile à apprendre, qui a su évoluer.
C'est un langage utilisé depuis 50 ans : les programmes et bibliothèques de routines disponibles se chiffrent en millions de lignes. On trouve toujours le bout de code que l'on cherche dans une bibliothèque existante.
C'est un langage universel : la plupart des scientifiques mathématiciens, physiciens, chimistes, géologues, biologistes, ingénieurs, etc. connaissent ou ont côtoyé FORTRAN.

Si vous êtes ici, c'est que vous avez entendu parler de FORTRAN ou même êtes-vous entrain de "sécher" sur votre premier programme. Je vous propose donc d'examiner les principales caractéristiques de ce langage et d'en étudier quelques applications simples au calcul scientifique.
Cette page n'a pas la prétention de faire du lecteur un programmeur confirmé mais de donner quelques indications pour aider les débutants en FORTRAN 77. Nous n'aborderons que les notions essentielles pour écrire des programmes de physique numérique. Elles seront introduites à mesure du besoin.

Un peu d'histoire

Le premier compilateur FORTRAN fut conçu et implémenté par John Backus (IBM Lab, San José Californie) en 1954. Il s'agissait du FORTRAN I. L'objectif de Backus était de concevoir un langage de programmation qui permettrait de réaliser des programmes en utilisant une syntaxe la plus proche possible de la syntaxe mathématicienne. Le terme FORTRAN traduit cette ambition, puisqu'il signifie FORmula TRANslation.
FORTRAN I était un peu primaire, mais quand même très supérieur aux assembleurs et langages machines utilisés au début des années 50. Il évolua rapidement en FORTRAN II et III.
FORTRAN prit son essor dans les années 60, surtout à partir de 1966, date à laquelle l'ANSI (organisme de standardisation américain) publia la norme ANSI66 qui normalisa le FORTRAN IV. A partir de cette date, tous les compilateurs utilisaient la même syntaxe. Un chercheur pouvait donc utiliser un code source écrit en FORTRAN IV sur n'importe quelle machine (bien qu'à l'époque, en calcul scientifique, en dehors des machines IBM, il n'y avait pas vraiment de salut...)
Par la suite, FORTRAN subit plusieurs évolutions, marquées par la publication par l'ANSI de normes :

1978, la publication de la norme ANSI 77, qui donna le FORTRAN 77, encore largement utilisé aujourd'hui (les exemples ci-dessous sont écrits en 77). Le FORTRAN 77 introduit la programmation structurée dans les habitudes des programmeurs FORTRAN.
1992, la publication de la norme ANSI 90, et pour la première fois la reconnaissance ISO (normalisation internationale). Le FORTRAN 90, outre quelques améliorations dans la gestion de la mémoire, amène essentiellement les outils nécessaires au calcul parallèle sur les machines multiprocesseurs. Il se décline en une nouvelle norme, qui n'est pas ANSI, et que l'on rencontre chez certains constructeurs, dont IBM, le HPF pour High Performance FORTRAN.
1997, la dernière norme en date, qui est une norme ISO, le FORTRAN 95.

Depuis la rédaction de cette page en 2006, le FORTRAN a encore évolué. Nous en sommes aujourd'hui (en décembre 2018) à la norme Fortran 2018 (ISO/CEI 1539-1:2018). POur ma part, je n'utilise plus FORTRAN 77 mais FORTRAN 90 ou 95. Néanmoins, les fondements du langage demeurent inchangés.

Le FORTRAN est donc un langage qui a plus d'un demi-siècle d'existence, mais qui reste un des plus utilisés en calcul scientifique.

Les bases du FORTRAN

Organisation d'un programme FORTRAN

Comme tous les programmes, un programme FORTRAN est un ensemble organisé d'instructions données par le programmeur à la machine. Je supposerais que le lecteur sait faire la différence entre un code source (écrit avec un éditeur, dans un langage de programmation, ici le FORTRAN) d'un code exécutable (construit par le compilateur et le linker, en langage machine). Nous parlerons ici essentiellement de code source (ou programme source, on dit aussi "source").
Un programme FORTRAN, du moins son code source, le plus simple est organisé comme suit :

 PROGRAM ProgName

              <liste d'instructions>

              END

Les mots PROGRAM et END sont les mots clé FORTRAN qui délimitent le programme principal. Le nom de ce programme sera ProgName.
La liste d'instructions bornée par PROGRAM et END forme un bloc fonctionnel. Tel quel (en enlevant <liste d'instructions>, tout de même!), ce bloc fonctionnel forme un programme FORTRAN tout à fait correct qui pourrait être compilé et exécuté. Evidement, il ne ferait pas grand chose..
On peut adjoindre à ce programme principal deux autres types de blocs fonctionnels:

 SUBROUTINE SubName (Arguments)

              <liste d'instructions>

              END

              et
              type FUNCTION FunctionName (Arguments)

              <liste d'instructions>

              END

Ces deux blocs sont utilisés pour alléger le programme principal de calculs ou de traitements répétitifs. Nous verrons plus loin comment les utiliser.
Un programme FORTRAN classique ressemble donc à quelque chose comme ça:

 PROGRAM ProgName

              < liste d'instructions >

              END

              SUBROUTINE SubName (Arguments)

              <liste d'instructions>

              END

              type FUNCTION FunctionName (Arguments)

              <liste d'instructions>

              END

Note: il est possible de disposer les blocs SUBROUTINE et FUNCTION avant le bloc PROGRAM. C'est une affaire d'habitude et de convenance personnelle. Une seule chose compte : la lisibilité de votre code (on y reviendra dans les exemples plus loin). Par exemple, évitez d'écrire un bloc PROGRAM de 5 000 lignes !
Il reste maintenant le plus gros du travail : identifier et comprendre ce que l'on entend par <liste d'instructions>
Comme tous les langages informatiques, FORTRAN utilise des données et des opérateurs. Nous commencerons par examiner les données, leur type et leur forme.

Les types de données

Les principaux types de donnée simples (nous verrons les tableaux plus loin) que nous utiliserons sont:

INTEGER : entier signé sur 4 octets (31 bits pour la valeur et un bit de signe). Nombre compris entre -2³¹ et 2³¹ -1
REAL : nombre "réel" codé en virgule flottante sur 4 octets. Nombre égal à +ou- 0.m*2^e , où la mantisse m est codée sur 23 bits et l'exposant e sur 8 bits signés. En gros, un REAL peut prendre en valeur absolue une valeur entre 1.4*10^-45 et 3.4*10³⁸
DOUBLE PRECISION : REAL dont la mantisse est codée sur 52 bits et l'exposant sur 11 bits. Un DOUBLE PRECISION peut prendre une valeur comprise entre 4.9*10^-324 et 1.8*10³⁰⁸
COMPLEX : 2 REAL dans un même nombre

Il en existe d'autres (CHARACTER,LOGICAL) que nous verrons plus loin.

Les constantes

Selon le type de la constante son écriture sera différente (attention, c'est un piège classique):

une constante INTEGER s'écrira simplement 1, -1, 3, 0
une constante REAL s'écrira obligatoirement avec un . ou le E de la notation en virgule flottante, par exemple: 0. , 1.24 , 31415E-4, 6.02E23
une constante DOUBLE PRECISION doit obligatoirement s'écrire en virgule flottante, avec un D en lieu et place du e : 31415D-4 , 0D0, 1.D0
une constante complex s'écrira sous forme d'un couple entre parenthèse en respectant les conventions des REAL : (0.0,1.0) , (1.1E-10, 3.5E-5)

On peut (c'est même une bonne habitude) nommer ses constantes en utilisant le mot clé PARAMETER comme dans :

 REAL pi

              PARAMETER (pi = 3.14)

Notez que j'ai d'abord précisé que pi est un REAL.

Les variables

Pour construire une variable dans un programme, on énonce tout simplement le type de la variable et son nom. Par exemple :

 INTEGER toto

              COMPLEX z

              REAL x,y,z

On définit les variables d'un programme, d'une SUBROUTINE ou d'une fonction entre le mot clé de début de bloc (PROGRAM,par exemple) et la première instruction exécutable.
Il y a plusieurs règles à connaître à propos des variables, dont certaines sont spécifiques à FORTRAN :

Le nom d'une variable peut contenir des chiffres, des lettres, un blanc ou un underscore (_). Il doit commencer par une lettre. Le FORTRAN est indifférent aux minuscules/majuscules. J'ai pris l'habitude d'écrire mes noms de variables en minuscules pour les distinguer des mots clé. Ce n'est pas très conforme au standard FORTRAN des puristes, mais c'est personnel...
Seuls les 6 premiers caractères d'une variable sont significatifs. Attention au piège : MinimumX et MinimumY représentent la même variable ! Bug assuré !!
Sauf indication contraire, en FORTRAN les variables commençant par i,j,k,l,m et n sont implicitement considérées comme des INTEGER. Les autres lettres désignent un REAL. Cela permet de ne pas déclarer explicitement les variables. Je vous déconseille formellement cette pratique. C'est la raison pour laquelle vous verrez en tête de tous mes programmes la déclaration IMPLICIT NONE, qui oblige la déclaration de chaque variable.
Chaque variable est locale par rapport à son bloc fonctionnel. Si vous déclarez un INTEGER toto dans le programme principal, il sera inconnu de toutes les subroutines et fonctions de ce programme. Pour partager des variables entre plusieurs blocs, on utilise l'instruction COMMON (voir son utilisation plus loin).

Les opérateurs et les fonctions

FORTRAN supporte les opérateurs arithmétiques suivants:

Addition +
Soustraction -
Multiplication *
Division /
Puissance **

Attention à la conversion de type dans une opération! Par exemple, si vous écrivez a**2, le résultat sera du type de a, car 2 est une constante entière. Si vous écrivez b**2 - 4*a*c et que a et c sont REAL alors que b est DOUBLE PRECISION, le résultat sera DOUBLE PRECISION. Mais le calcul 4*a*c sera fait en REAL. Ce genre de piège, qui se traduit par une perte de précision, est assez courant. Attention donc à la cohérence de vos typages de variables.
FORTRAN supporte un grand nombre de fonctions mathématiques: SIN, COS ( en REAL, DOUBLE PRECISION et même COMPLEX), TAN, ATAN, ASIN, ACOS, LOG10, LOG, EXP, SQRT, etc.. On les retrouvera dans nos programmes.
Signalons aussi les fonctions MIN et MAX, la fonction ABS, les fonctions INT (partie entière) et MOD (reste de la division entière), dont les noms varient selon le type des variables.

Les entrées/sorties

Voyons les instructions qui nous permettent d'afficher un résultat sur l'écran ou de saisir un paramètre depuis notre clavier. A ce stade nous en resterons à ce genre d'entrée/sortie. Nous parlerons des fichiers plus loin.

Lire

Pour lire les informations saisies depuis un clavier, il existe deux méthodes :

- la simple (et toujours applicable sauf aux chaînes de caractères), qui est la lecture non formatée. La donnée saisie prend le type de la variable dans laquelle on la stocke. Cela s'écrit:

 REAL x

              READ(*,*) x

- un peu plus compliquée, la lecture formatée. On indique au programme un certaine nombre d'informations en particulier le nombre de chiffres avant et après la virgule et bien d'autres choses. cela s'écrit, par exemple:

 REAL x

              READ(*,'(e12.5)') x

ce qui indique qu'on lit un réel en précisant que l'on veut 12 caractères, dont 5 chiffres pour la mantisse. Nous ne nous servirons pas de cette possibilité pour l'instant.
Le premier * dans le READ indique le canal de lecture. Si on n'en utilise qu'un seul (le clavier, par exemple), il est inutile de le préciser. Sinon, on indique un nombre entier compris entre 10 et 99, qui identifie le canal pour chaque opération d'entrée/sortie (voir plus loin la cas d'utilisation des fichiers).

Ecrire

Les principes sont les mêmes. On utilisera l'écriture non formatée sur l'écran. Pour ce faire, la bonne instruction est:

WRITE(*,*) x

Deux instructions à connaître pour saisir une information depuis le clavier de manière sympathique:

 WRITE(*,'(a,$)') 'Rayon de la sphère : '

              READ(*,*) rayon

Notez la syntaxe du WRITE. Comme je l'ai dit plus haut, la lecture et l'écriture non formattée d'une chaîne de caractères est interdite. Vous utiliserez donc le format indiqué.

Les particularité d'un code source FORTRAN

Que faut-il savoir de plus pour commencer à écrire un petit programme FORTRAN ? Quelques spécificités propres à FORTRAN....
Historiquement, le code source FORTRAN était saisi sur des cartes perforées. Ce qui induisait certaines contraintes qui ont survécues jusqu'en FORTRAN 77. Signalons les principales:

une ligne de commentaire doit comporter un C ou un * en colonne 1
les instructions doivent commencer en colonne 7 et leur longueur ne doit pas dépasser la colonne 72 (dans les codes exemples ci-dessous, je ne suis pas la règle...)
les étiquettes (qu'on n'utilise plus tellement vu qu'on les trouve dans les GO TO, honnis, et les écritures formattées) doivent être écrites dans les colonnes 2 à 5
pour couper une instruction en plusieurs lignes, on positionne un caractère quelconque (je mets un !) dans la colonne 6 de chacune des lignes suivant la première ligne.

J'ai dit tout à l'heure que FORTRAN n'était pas sensible à la casse. Mais vous avez sans doute remarqué aussi que j'écrivais tous les mots clé en majuscule. Ce n'est pas obligatoire, mais c'est une coutume en FORTRAN (elle provient des premiers standards de FORTRAN).
Ces règles ont évolué avec les dernières moutures de FORTRAN (90 et 95). Dans ces versions, elles sont moins contraignantes. Je les cite parce que j'utilise ici un compilateur 77, qui est encore très largement utilisé, mais aussi parce qu'elles font partie du folklore...
Tout ça à l'air un peu compliqué et ringard, mais on s'y fait très vite. D'autant que les éditeurs de textes, pour la plupart, connaissent et gèrent les règles de FORTRAN.

Notre premier programme

Le programme VolumeSphere

Pour notre premier programme, nous calculerons le volume d'une sphère en fonction de son rayon que nous indiquerons au programme. J'ai préféré un programme de calcul, même trivial, plutôt que le classique "Hello the word" de C. Nous sommes des scientifiques que diable ! Et puis, ce genre de programme nous permet de mettre en oeuvre les notions introduites ci-dessus.

Le code source de ce programme:

 C Programme de calcul du volume d'une sphere

              PROGRAM VolumeSphere

              IMPLICIT NONE

              C Definition des parametres

              REAL PI

              PARAMETER (PI=3.14159265)

              C Definition des variables

              REAL rayon, volume

              C Saisie du rayon de la spher

              WRITE(*,'(a,$)')'Rayon de la sphere : '

              READ(*,*) rayon

              C Calcul du volume

              volume = (4*PI/3)*rayon**3

              C Affichage

              WRITE(*,*) 'Volume de la sphere = ', volume

              END

Il commence par une ligne de commentaire indiquant l'objet du programme. Bonne habitude à prendre que de commenter son code.
Puis le mot clé PROGRAM suivi de son nom (Attention, là aussi seuls les 6 premiers caractères sont significatifs).
Notez ensuite l'instruction IMPLICIT NONE, qui m'oblige à déclarer toutes mes variables. Ainsi, je n'en n'oublierai pas.
Viennent les déclarations des constantes et autres paramètres: déclarer ici tous vos paramètres symboliques (PI, constante de Plank, densité de l'eau, distance terre-lune, etc...). Attention à la syntaxe et à déclarer le type des constantes.
Déclaration des variables: simple et de bon goût.
Et enfin le corps du programme:

La saisie de la variable rayon. J'aurai pu écrire WRITE(*,'("Rayon de la sphere: ",$)'), ce qui aurait donné la même chose. Notez le $, qui empêche le retour à la ligne après l'affichage de la chaîne 'Rayon de la sphere : '.
Le calcul du volume: rien à dire tellement cela ressemble à la formule mathématique
Puis l'affichage, comme on l'a définit ci-dessus. En mode formaté, j'aurai pu écrire WRITE(*,'("Volume de la sphere = ",e15.8)') volume

qui se termine par END.
Voilà un tout petit programme, mais qui contient tous les germes d'un grand code. Voyons comment le saisir, le compiler et l'exécuter.

Installer un environnement de développement

Pour créer et exécuter un programme FORTRAN il faut un minimum d'outils. A commencer par un ordinateur.
Là, survient mon premier dilemme: dois-je m'orienter vers le monde Windows ou Linux? Personnellement je préfère et je travaille sous Linux (Ubuntu). Mais je suppose que la plupart d'entre vous sont plutôt sous Windows. Il faut savoir que FORTRAN se moque du système d'exploitation (OS): la syntaxe est la même sous Linux ou Windows. Seuls les outils différent, et encore ! Les environnements de développement FORTRAN existent souvent pour les deux OS. Quant au compilateur, j'utilise G77, un compilateur GNU disponible à la fois sous Windowss, Linux et Mac OSX.
Bref, j'ai opté, pour les exemples de cette page, de travailler sous Windows. Les linuxiens me pardonneront (je suis des leurs!) et se consoleront en se disant que tout ce que j'écris ici est directement transposable (sauf les commandes système, évidement... mais les linuxiens sont rarement des informaticiens débutants!).
Donc, il faut rassembler quelques programmes, que nous allons trouver gratuitement sur le net :

un environnement de développement (IDE) qui associe un éditeur, un compilateur et un linker. J'ai choisi VFort. VFort travaille en FORTRAN 77. Il comporte une aide et un guide de programmation. Il est téléchargeable ici. Attention : il semble que cet IDE ne soit plus supporté. Aujourd'hui (2018), il est préférable d'utiliser Eclipse Photran, qui existe sous Windows, Linux et MacOSX et qui permet l'usage du compilateur G77 qui existe toujours lui... Vous pouvez aussi, c'est plus efficace, travailler en mode commande, avec n'importe quel éditeur de texte.
un logiciel de tracé de courbes, pour visualiser nos résultats. J'ai choisi GnuPlot, qui existe sous Windows et Linux. Il est téléchargeable sur le site www.gnuplot.info

Puis je vous laisse installer ces deux programmes.
Il serait pratique de créer un répertoire pour nos petits travaux de programmation. Le mien s'appelle numlab\labphysique.

Saisir le programme

Lancer le programme VFort. La vue principale de l'IDE s'affiche.

La fenêtre de gauche liste les différents fichiers de votre projet.
La fenêtre de droite est notre fenêtre de travail. C'est là que nous saisirons notre code.
La fenêtre du bas est la console: c'est là que s'afficheront les différents messages de compilation (des messages d'erreur, souvent, mais aussi le résultat de la compilation) et de link (aussi des messages d'erreur...).

Pour saisir notre programme, il faut:

Créer un nouveau fichier (dans le menu FILE/New). Il porte un nom par défaut (Fort1, par exemple)
Saisir le texte du programme (ou faire un copier/coller). Attention au respect des colonnes : FORTRAN 77 est intransigeant sur le sujet...
Sauvegarder le fichier dans votre répertoire (FILE/Save), en le nommant, VolumeSphere. On donne généralement le nom du programme au fichier source du programme principal. Le nom du fichier portera une extension .for
Créer un projet (Project/Create project) dans votre répertoire, en le nommant, VolumeSphere par exemple ou tout autre nom
Attacher le fichier que vous venez de créer (VolumeSphere dans notre exemple) au projet (Project/Add Files)
ATTENTION : sous Windows, pensez à ajouter à la variable d'environnement système PATH le chemin d'accès aux programmes de VFORT (le compilateur, le linker) sinon, vous auriez une erreur au build du projet. Dans mon cas, le chemin est C:\Program\VFort\Bin. Attention aussi en manipulant la variable PATH à ne pas effacer les chemins existants!

Le compiler - linker - exécuter

Et maintenant on peut compiler et linker le fichier source, en faisant Project/Build. Si vous n'avez fait aucune erreur de syntaxe, vous devez voir apparaître les messages suivants (ou approchant..) dans la fenêtre inférieure:

Compiling (command: g77 -c) ...
E:\NumLab\LabPhysique\VolumeSphere\VolumeSphere.for
Linking (command: g77) ...
Code generation: VolumeSphere.exe

Tout c'est bien passé. Sinon, corrigez vos erreurs indiquées par le compilateur et recommencez à l'étape Project/Build.
Vous noterez que le programme executable porte le nom donné dans l'instruction PROGRAM. Il se trouve dans le sous-répertoire Output de votre répertoire de travail (ici NumLab\LabPhysique\VolumeSphere\).
Pour le lancer, sous Windows, il faut ouvrir une fenêtre de commande (cmd.exe) et déclarer le répertoire output comme répertoire courant. Puis tapez VolumeSphere et Return. Le programme s'exécute.
Un conseil, gardez la fenêtre de commande ouverte, cela vous évitera des opérations fastidieuses.
Et voilà, votre premier programme fonctionne ! Passons à des choses un peu plus compliquées...

Un peu plus compliqué

Dans cette partie un peu longue, on va examiner:

les conditions
les boucles
les tableaux
les fonctions et subroutines
les fichiers

Les conditions

Pour aborder les conditions, il faut d'abord parler d'expression logique.
Il existe en FORTRAN un type de variable LOGICAL. Une variable de type LOGICAL peut prend deux valeurs: .TRUE. ou .FALSE. Il faut savoir que .TRUE. vaut 1 et .FALSE. vaut 0.
FORTRAN définit les expressions logiques suivantes, dans lesquelles x et y sont deux variables :

x .EQ. y signifie x = y
x .NE. y signifie x <> y
x .GT. y signifie x > y
x .LT. y signifie x < y
x .GE. y signifie x >= y
x .LE. y signifie x <= y

Il définit aussi les opérateurs logiques unaires:

.OR. OU
.AND. ET
.XOR. OU exclusif
.NOT. NON

Munis de ces expressions logiques, on peut définir les instructions conditionnelles suivantes :

 IF (condition) THEN

                <liste d'instructions>

              ENDIF

La signification de ce jeu d'instruction est simple : si l'expression condition est vraie, par exemple x .LT. 0, alors la liste d'instructions contenue entre le IF et le ENDIF est exécutée. Sinon, le programme continu et exécute l'instruction immédiatement suivante le ENDIF.
Il existe plusieurs variantes d'instructions conditionnelles. Par exemple, la condition alternative :

 IF (condition) THEN

                <liste d'instructions 1>

              ELSE

                <liste d'instructions 2>

              ENDIF

Si l'expression condition est vraie, alors la liste d'instructions 1 est exécutée (et pas la liste d'instructions 2). Si l'expression condition est fausse, c'est la liste d'instruction 2 qui est exécutée. Il s'agit en fait de la forme complète de l'instruction conditionnelle IF.
Le langage C possède une instruction très sympathique, le 'switch', qui permet de faire des traitements différents en fonction de la valeur d'une variable. Cette instruction n'existe pas en FORTRAN 77, mais on peut la simuler en utilisant le jeu d'instructions suivant :

 IF (condition1) THEN

                <liste d'instructions 1>

              ELSE IF (condition2) THEN

                <liste d'instructions 2>

                ELSE IF (condition3) THEN

              <liste d'instructions 3>

              ELSE

                <liste d'instructions n>

              ENDIF

Vous en devinez l'utilisation! Supposons que je doive exécuter un traitement différent selon qu'une variable x est égale à 1, 2, 3 ou différente de 1, 2, et 3. Je poserai :

 IF (x .EQ. 1) THEN

                <liste d'instructions 1>

              ELSE IF (x .EQ. 2) THEN

                <liste d'instructions 2>

              ELSE IF 

              (x .EQ. 3) THEN

                <liste d'instructions 3>

              ELSE

              C cas de x différent de 1,2 et 3

              <liste d'instructions 4>

              ENDIF

Si x = 1, seule la liste d'instructions 1 est exécutée. Si x=2, seule la liste d'instructions 2 est exécutée, etc...
C'est un jeu d'instruction très employé en traitement des données, mais qui peut générer des programmes un peu fouillis! Les commentaires sont indispensables.

Les boucles

Sans doute, ce sont les instructions les plus utilisées en calcul. FORTRAN 77 comporte deux types de boucles.

La boucle DO...ENDDO

Sa syntaxe

 DO var= debut, fin, pas

                <liste d'instructions>

              ENDDO

où var est la variable de boucle, qui évolue de debut à fin selon un incrément de pas. Notons que pas peut être omis, dans ce cas, il vaudra 1.

Son mode de fonctionnement

La liste d'instructions est exécutée un certain nombre de fois, déterminé par la valeur de debut et de fin. En principe, on ne doit pas sortir de la boucle autrement que lorsque var > fin, qui est la sortie naturelle par le ENDDO. C'est techniquement possible mais à proscrire rigoureusement.
var, debut, fin et pas sont des variables de type 2018 : Fortran 2018 (ISO/CEI 1539-1:2018.
pas peut être négatif. Dans ce cas debut doit être plus grand que fin.

Je vous laisse découvrir ce qui se passe si fin > debut, selon la valeur de pas (>0 ou <0)

Exemple

 somme = 0

              DO i = 1, 100

                somme = somme + i

              ENDDO

Cette boucle calcule la somme de 100 premiers entiers.
Curiosité: on peut écrire une boucle infinie en écrivant:

 DO

                <liste d'instructions>

              ENDDO

Dans ce cas, on sort de la boucle en utilisant l'instruction exit en un lieu approprié de la liste d'instructions. Ce genre de manip est à éviter dans un premier temps...

la boucle DO ... WHILE

Sa syntaxe

 DO WHILE(condition)

                <liste d'instructions>

              ENDDO

où condition est une expression logique.

Son mode de fonctionnement

Le programme entre dans la boucle si condition = .TRUE. Il reste dans la boucle tant que (i.e. la liste d'instruction est exécutée autant de fois que) condition vaut .TRUE.
Lorsque condition vaut .FALSE. l'exécution se poursuit à la première instruction qui suit le ENDDO. Cela implique qu'il existe dans la liste d'instructions une instruction susceptible de faire passer condition à .FALSE. Sinon le programme bouclera indéfiniment, ce qui rarement le but recherché (mais cela arrive parfois dans les traitements d'acquisition).

Exemple

 C déclaration des variables

              INTEGER n

              DOUBLE PRECISION somme, epsilon

              LOGICAL fin

              C Initialisation des variables

              epsilon = 1.d-10

              somme = 0d0

              fin = .FALSE.

              C Boucle de calcul. La condition est .NOT. fin de telle sorte que
              la condition soit vraie à l'entrée

              DO WHILE (.NOT. fin)

              somme = somme + 1d0/n

              n = n +1

              C On sort de la boucle si fin devient TRUE

              fin = (1d0/n .LT. epsilon*somme)

              ENDDO

Cette boucle calcule la suite des inverses jusqu'à ce que le terme 1/n soit inférieur à epsilon fois la somme partielle.
Attention, dans cet exemple à la gestion de la condition. C'est un cas très répandu mais un peu tordu. On pose fin = .FALSE. et la condition de test de la boucle est .NOT. fin. Cela signifie qu'on entre dans la boucle, car au premier tour (.NOT. fin) = .TRUE.
Lorsque la condition mathématique de convergence est remplie fin devient .TRUE. Le test de boucle nous donne donc .NOT. fin = .FALSE. et l'on sort de la boucle. Tout ça pour que le test est un sens mathématique usuel...

Les tableaux

Le dernier élément indispensable pour faire du calcul numérique! Je ne sais pas ce que l'on ferait sans tableau... J'ai l'habitude, comme beaucoup de mes collègues de les appeler des 'vecteurs' ou des 'matrices'. Vous allez vite comprendre !
On déclare un tableau en FORTRAN de la façon suivante:

type Nom(dim1, dim2, dim3, ...)>

type désigne le type des variables stockées dans le tableau nommé Nom. Le nombre de constantes dim indique la dimension du tableau. Un vecteur est donc un tableau de 1 dimension...
dim1, dim2, ... déterminent la taille du tableau. Ce sont des constantes entières. FORTRAN 77 ne gère pas les tableaux dynamiques, ni même le redimensionnement dynamique des tableaux. Il est de bon ton de fixer les valeurs dim1, dim2, ... sous forme de paramètres.

Par exemple:

 INTEGER n

              PARAMETER (n = 1000)

              REAL toto(n)

On vient de construire un vecteur de 1000 réels. Notons qu'un tel tableau occupe 4000 octets. Où encore une matrice carrée de 1000*1000 DOUBLE PRECISION :

DOUBLE PRECISION matcarre(n,n)

qui occupe elle 1000*1000*8 octets.....

Voilà pour la déclaration des tableaux. Voyons maintenant l'utilisation de ceux-ci :

Premier point important, en FORTRAN, un tableau s'adresse en partant de 1, contrairement au C où l'on part de 0.
Deuxième point important: l'indice d'un tableau ne doit jamais dépasser la dimension fixée à sa déclaration. Sinon, plantage assuré.

Pour stocker une variable dans un tableau, il faut préciser la position du tableau dans laquelle on veut écrire. Par exemple :

 i = 2

              toto(i) = cos(x)

 i = 1

              j = 1

              matcarre(i,j) = toto(i)

Cette dernière écriture signifie : je prends la valeur contenue dans la case i du tableau toto et je l'assigne à la case (i,j) de matcarre.
Voici quelques exemples superclassiques d'utilisation de tableaux.

La somme de deux matrices carrées :

 INTEGER n

              PARAMETER (n = 1000)

              DOUBLE PRECISION m1(n,n), m2(n,n), m3(n,n)

              DO i=1, n

              DO j=1, n

              m3(i,j) = m1(i,j) + m2(i,j)

              ENDDO

              ENDDO

Le produit scalaire de deux vecteurs :

 INTEGER n

              PARAMETER (n = 1000)

              DOUBLE PRECISION u(n), v(n)

              DOUBLE PRECISION Scalar

              Scalar = 0d0

              DO i=1, n

              Scalar = Scalar + u(i)*v(i)

              ENDDO

Les fonctions et subroutines

Les fonctions et les subroutines sont deux notions importantes de programmation. On emploie une fonction ou une SUBROUTINE lorsqu'on veut coder un traitement qui est souvent répété dans un programme. Il peut s'agir :

d'un calcul qui retourne une valeur simple à partir de variables. Par exemple, calculer le volume d'une sphère à partir de son rayon. Dans ce cas, on utilise une fonction. Sa principale caractéristique est de ne savoir retourner au programme qu'une seule valeur.
d'un traitement qui retourne plusieurs valeurs. Par exemple, on peut chercher les extrema et minima locaux dans un tableau. Dans ce cas on utilise une SUBROUTINE.

Les subroutines servent aussi fréquemment à matérialiser la décomposition fonctionnelle d'un problème. On découpe un problème compliqué en plusieurs blocs fonctionnels différents. Cela permet de mieux comprendre le code, mais aussi de pouvoir le faire réaliser par plusieurs programmeurs. Dans ce cas, on est vite confronté au besoin de partager des données entre subroutines.
ATTENTION: une fonction (resp. une SUBROUTINE) est un bloc fonctionnel différent du programme principal. Chaque bloc fonctionnel possède ses propres variables locales qui ne sont pas connues des autres. Ainsi, une variable peut être nommée toto dans le programme principal et une autre variable porter le même nom dans une fonction (resp. une SUBROUTINE) sans pour autant qu'il s'agisse de la même variable.
De même, le nom des paramètres dans l'en-tête d'une fonction (resp. d'une SUBROUTINE) n'a absolument pas de signification en dehors de cette fonction (resp.SUBROUTINE).

Les fonctions

La structure générale d'une fonction est la suivante :

 type FUNCTION NomFonction(par1, par2, ...)

              type par1

              type par2

              C déclaration des variables locales

              C Traitement effectué par la fonction

              <liste d'instructions dont une instruction sera de la forme
              NomFonction = ... >

              RETURN

              END

Depuis un programme, on appelle une fonction comme n'importe quelle fonction prédéfinie de FORTRAN. Par exemple :

 PROGRAM toto

              REAL R, V

              R = 1.

              C Notez que la variable dans toto s'appelle R, alors que la
              paramètre formel de VolumSphere s'appelle rayon

              V = VolumeSphere(R)

              END

              

              C définition de la fonction VolumeSphere

              REAL FUNCTION VolumeSphere(rayon)

              REAL rayon

              C Definition des parametres - Il s'agit de variables locales, PI
              est inconnu du programme toto

              REAL PI

              PARAMETER (PI=3.14159265)

              VolumeSphere = (4*PI/3)*rayon**3

              RETURN

              END

Attention, lorsque vous appelez une fonction depuis votre programme (on va voir comment), vous devez vous assurer que les paramètres utilisés dans l'appel sont bien strictement du même type que les paramètres déclarés dans l'en-tête de définition de la fonction.

Les subroutines

La structure générale d'une subroutine est la suivante :

 SUBROUTINE NomSub(par1, par2, ...)

              type par1

              type par2

              ....

              C déclaration des variables locales

              C Traitement effectué par la fonction

              <liste d'instructions>

              RETURN

              END

Depuis un bloc fonctionnel (programme, autre SUBROUTINE, fonction), on appelle une SUBROUTINE par un CALL, comme dans l'exemple ci-dessous. ATTENTION, FORTRAN 77, contrairement à FORTRAN 90, ne supporte pas les appels récursifs (une SUBROUTINE ne peut pas s'appeler elle-même).

 PROGRAM toto

              REAL R,V

              R = 1.

              CALL VolumeSphere(R,V)

              END

              C définition de la SUBROUTINE VolumeSphere

              SUBROUTINE VolumeSphere(rayon, volume)

              REAL rayon, volume

              C Definition des parametres - Il s'agit de variables locales, PI
              est inconnu du programme toto

              REAL PI

              PARAMETER (PI=3.14159265)

              volume = (4*PI/3)*rayon**3

              RETURN

              END

Attention, comme pour les fonctions, lorsque vous appelez une SUBROUTINE depuis votre programme, vous devez vous assurer que les paramètres utilisés dans l'appel sont bien strictement du même type que les paramètres déclarés dans l'en-tête de définition de la SUBROUTINE.
Nous verrons dans l'exemple du programme Euler l'utilisation de ces différents blocs fonctionnels.

Les fichiers

C'est un vaste sujet ! Nous nous limiterons sur cette page à examiner les instructions nécessaires à créer, lire et écrire les fichiers qui conviennent pour stocker nos résultats de calcul ou nos données d'entrée.
On ne considérera donc que les seuls fichiers séquentiels formattés (autrement appelés fichier texte). Ces fichiers ont la propriété intéressante de pouvoir être lus et modifiés avec n'importe quel éditeur de texte.

Créer et ouvrir un fichier séquentiel formaté

La syntaxe de l'instruction pour ouvrir un fichier est :

 OPEN(unité,

              FILE= chaine,

              FORM = chaine,

              STATUS = chaine

              ERR = numlabel)

où

unité : numéro d'identification de l'unité du fichier. C'est un nombre entier compris entre 10 et 99, choisi arbitrairement par le programmeur. Il existe des conventions historiques, dont on ne s'encombrera pas ici. Chaque fois qu'une instruction voudra adresser le fichier ouvert, elle le désignera par ce numéro d'unité.
FILE : chaîne de caractères contenant le nom et le chemin du fichier, par exemple 'd:\NumLab\LabPhysique\data.txt'
FORM: chaîne de caractères: 'formatted' pour les fichiers formatés (valeur par défaut) ou 'unformatted' pour les fichiers binaires (dont nous ne parlerons pas)
STATUS : chaîne de caractères :
- 'new' : le fichier n'existe pas, il est créé.S'il existe déjà le programme retourne une erreur
- 'old' : le fichier existe déjà. S'il n'existe pas, le programme retourne une erreur
- 'unknow' : ouvre le fichier, s'il n'existe pas il le créé, s'il existe il écrit à la fin de celui-ci. C'est le mode par défaut
ERR : numéro de label vers lequel le programme continuera son exécution s'il rencontre une erreur

Exemple : OPEN(10, FILE = 'd:\NumLab\LabPhysique\data.txt') ouvre le fichier data.txt, qu'il existe ou non (par défaut). En cas d'erreur, le programme se plante...

Lire un fichier séquentiel formaté

La syntaxe de l'instruction de lecture est :

 READ(unité,

              format,

              ERR = numlabel,

              END = numlabel) liste de données

où

unité : numéro d'identification de l'unité du fichier alloué à son ouverture
format : dans notre cas (fichier formaté), on indiquera *
ERR : numéro de label vers lequel le programme continuera son exécution s'il rencontre une erreur
END : numéro de label vers lequel le programme continuera son exécution lorsqu'il n'aura plus rien à lire dans le fichier (optionnel)

Exemple : READ (10,*,ERR=90, END =100)toto lit une donnée dans le fichier data.txt et la stocke dans la variable toto.

Ecrire une fichier séquentiel formaté

La syntaxe de l'instruction de lecture est:

 WRITE(unité,

              format,

              END = numlabel) liste de données

où

unité : numéro d'identification de l'unité du fichier alloué à son ouverture
format : dans notre cas (fichier formaté), on indiquera *
END : numéro de label vers lequel le programme continuera son exécution lorsqu'il n'aura plus rien à lire dans le fichier (optionnel)

Exemple : WRITE (10,*) toto écrit la variable toto dans le fichier data.txt.

Fermer un fichier

Facile... On utilise l'instruction CLOSE(unité).
En principe, tous les fichiers sont fermés automatiquement à la fin d'un programme. Mais c'est une bonne habitude à prendre que de les fermer par un CLOSE lorsqu'on a fini de s'en servir.

Exemple d'utilisation

Supposons que nous ayons stocké les résultats d'un calcul dans deux vecteurs x et y de dimension n. Nous voulons sauvegarder ces valeurs dans un fichier pour les plotter avec un logiciel comme Python, Matlab, Scilab ou GnuPlot.
Voilà le bout de code correspondant. On supposera que le fichier data.txt n'existe pas.

 INTEGER n, i

              PARAMETER (n = 100)

              DOUBLE PRECISION x(n), y(n)

              CHARACTER*80 nomfic

              nomfic = 'data.txt'

              OPEN(10, FILE=nomfic)

              DO i=1,n

              WRITE(10,*) REAL(x(i)), REAL(y(i))

              ENDDO

              CLOSE(10)

              END

Il y a plusieurs nouveautés dans ce code:

l'usage d'une variable CHARACTER et son assignation. On créé une variable contenant une chaîne de caractères, ici d'une longueur maximum de 80 caractères. On lui assigne simplement le nom du fichier que l'on écrit entre cotes.
la fonction REAL (ne pas confondre avec le type de variable!). La raison de l'usage de cette fonction, qui convertit un DOUBLE PRECISION en REAL est un peu curieux. Elle tient au logiciel de plot que vous utiliserez. certains, comme Matlab n'aime pas le 'D' utilisé par FORTRAN dans la notation double précision. Ils ne comprennent que le 'E' des REAL. Alors pour être sur de ne pas avoir d'ennuis, on convertit en REAL. Il risque d'y avoir quelques pertes de précision, mais pour faire de l'affichage graphique, cela a généralement peu d'importance.

Le premier programme de calcul

Simuler le mouvement d'un pendule simple

Nous voilà maintenant en possession de tous les éléments pour réaliser un programme de calcul en FORTRAN. Pour les mettre en oeuvre sur un exemple pratique, j'ai choisi de simuler le mouvement d'un pendule simple. Je ne reviens pas sur l'analyse physique du phénomène, que vous trouverez détaillée dans la page suivante de TangenteX.

L'équation différentielle non linéaire d'ordre 2 qui modélise ce mouvement est :

$ \dfrac{d^2\theta}{dt^2} = -\dfrac{g}{l}\sin(\theta) $

Pour résoudre cette équation différentielle, j'emploie la méthode d'Euler, décrite dans cette page de TangenteX, car c'est la plus simple à mettre en oeuvre, même si ce n'est pas la plus précise ni la plus efficace.

L'application du schéma d'Euler me donne les équations :

$\theta_{i+1} = \theta_i + h\omega_i $ (1)

$ \omega_{i+1} = \omega_i - h\dfrac{g}{l}\sin(\theta_i)$ (2)

Le programme

Nous retrouverons tout ce petit monde dans le programme suivant :

le tableau X1 contient les valeurs d'angle calculées
le tableau X2 contient les valeurs de vitesse angulaire calculées
le tableau t contient les pas temporels
les conditions initiales sont : angle initial 60° (attention à l'exprimer en rd) et vitesse initiale nulle.

 C Programme de demonstration de l'usage de la methode d'EULER
              

              C Calcul du mouvement d'un pendule simple 

              C Dominique Lefebvre mars 2006 

              

              PROGRAM PenduleEuler 

              IMPLICIT NONE 

              

              C Declaration des constantes 

              REAL g, PI, n, P 

              PARAMETER (g = 9.81, PI = 3.141592654) 

              PARAMETER (n = 2) ! Nb de cycles de calcul 

              PARAMETER (P = 1000) ! Nb de pas de calcul par cycle 

              

              C Declaration des variables 

              INTEGER i 

              REAL l, Theta, h 

              REAL X1(n*P), X2(n*P), t(n*P) ! resp. l'angle, la vitesse
              angulaire, le temps 

              

              C declaration des zone commune 

              COMMON /Systeme/l 

              

              C Initialisation des variables 

              l = g/(4*PI*PI) ! je choisis une longueur arbitraire qui m'arrange
              

              Theta = 2*PI*SQRT(l/g) ! la periode vaut 1 pour l'approximation
              des petits angles 

              h = Theta/P ! ! pas temporel pour le calcul 

              

              C Declaration des conditions initiales d'integration 

              X1(1) = 60*PI/180 ! angle initial en radian 

              X2(1) = 0 ! vitesse angulaire initiale 

              t(1) = 0.0 ! Origine des temps 

              

              C Calcul du mouvement selon l'algorithme d'Euler 

              DO i=1, n*P-1 

              t(i+1) = t(i)+h ! incrementation du temps 

              

              C Le systeme differentiel descriptif du pendule 

              CALL Derivee(X1(i), X2(i), f1, f2) 

              

              C L'algorithme d'EULER 

              CALL Euler(X1,X2,f1,f2,i,h) 

              ENDDO 

              

              C Sauvegarde des donnees dans le fichier Euler.dat 

              C On sauvegarde les données pour tracer la courbe theta = f(t)
              mais on pourrait 

              C tracer d'autre courbes, en particulier omega = f(theta). Il
              suffit de bien 

              C choisir les données à stocker 

              OPEN(10, FILE='Euler.dat') 

              DO i=1, n*P 

              WRITE(10,*) REAL(t(i)), REAL(X1(i)*180./PI) 

              ENDDO 

              CLOSE(10) 

              END 

              

              C Subroutine de calcul du systeme differentiel 

              C En cas d'utilisation d'un autre système différentiel, seule 

              C cette subroutine sera modifiee 

              SUBROUTINE Derivee(x,y,fx,fy) 

              REAL x,y,fx,fy 

              COMMON /Systeme/l 

              REAL g, l 

              PARAMETER (g = 9.81) 

              fx = y 

              fy = -(g/l)*SIN(x) 

              END 

              

              C Subroutine Algorithme d'Euler 

              SUBROUTINE Euler(x,y,f1,f2,i,dt) 

              REAL x(*), y(*), f1, f2, dt 

              INTEGER i 

              x(i+1) = x(i)+f1*dt 

              y(i+1) = y(i)+f2*dt 

              END

Le code source de ce programme est disponible dans la bibliothèque de codes de TangenteX.

Quelques remarques à propos de ce dernier programme :

L'usage de COMMON entre le programme principal et les subroutines. Cette instruction permet de partager des variables entre deux ou plusieurs blocs fonctionnels (je vous rappelle que les variables en FORTRAN sont locales). Attention, COMMON ne permet pas de partager les variables déclarées en PARAMETER, ce qui explique que j'ai recopié dans la subroutine Derivee le PARAMETER g. Le nom du COMMON est quelconque. Vous pouvez nommer un COMMON comme vous le souhaitez, à condition que ce nom soit unique dans le programme.
Les commentaires en fin de ligne. Vous avez noté le signe ! qui permet d'écrire des commentaires en fin de ligne. Pratique pour commenter ses instructions...

Utiliser GnuPlot pour tracer ses résultats

GnuPlot : présentation du traceur de courbes

J'utilise GnuPlot 3.7, que vous pouvez très facilement charger, que vous soyez sous OS X, Windows ou Linux (www.gnuplot.info). Je n'ai pas encore essayé la dernière version 4 ! Pour info, en 2018, nous en sommes à la version 5.2 et le logiciel s'est considérablement enrichi.

GnuPlot est un logiciel destiné au tracé et à l'impression des données en 2D et 3D. Il possède la majeure partie des fonctions utilisées quotidiennement pour faire du tracé de courbes.
Il n'est pas dans mon intention de faire de grands développements sur GnuPlot ici. Celui-ci possède une fonction d'aide (tapez help après le prompt gnuplot) qui vous donnera toutes les informations que vous voulez connaître. Nous allon voir ici les principales commandes pour tracer une courbe.

Les principales commandes de GnuPlot

Nous partirons du principe que nos données sont stockées dans un fichier texte, par exemple Euler.dat, produit par un programme de calcul.
Nous voulons tracer la courbe correspondante sur une feuille propre avec des légendes, un titre, etc. Voyons les principales commandes nécessaires. Je ne fais pas apparaître le prompt de GnuPlot..:

plot 'Euler.dat' : c'est la commande la plus simple. Je suppose que le répertoire de travail de gnuplot contient le fichier Euler.dat. Sinon, il faut le changer en utilisant la commande chdir du menu de Gnuplot. Vous obtenez un tracé avec les axes par défaut, sans légende, brut de fonderie.
set data style lines : pour modifier l'aspect de la courbe avant le tracé. On lui précise ici de relier les points par des segments, c'est plus joli! Un détail, cette option sera active tant que vous n'en changerait pas...
set xrange[0.0:2.0] ou set yrange[-8:8] pour fixer les bornes du tracé
set xlabel "Temps" pour indiquer la légende de l'axe des X. Evidement set ylabel "Angle" légende l'axe des Y.
set title "Variation de l'angle en fonction du temps" pour indiquer un titre à votre courbe...

Comment tracer un fichier de données produit par un code FORTRAN

Voyons comment tracer la courbe correspondante à l'exécution de notre programme Euler. Comme je suis un peu paresseux, je vais stocker toutes mes commandes gnuplot dans un fichier de commandes. Je créé donc un fichier avec un éditeur quelconque (Notepad par exemple), que j'appelle Euler.p. Voici son contenu :

 # Fichier de tracé de la courbe theta = f(t)selon la methode
              d'Euler

              set title "Variation de l'angle en fonction du temps"

              set xlabel "Temps"

              set ylabel "Angle"

              plot 'Euler.dat'

Comme vous l'avez deviné, le # désigne une ligne de commentaires
Je place ce fichier dans le répertoire qui contient le fichier Euler.dat, qui doit être le répertoire de travail de GnuPlot.
Maintenant, au prompt de GnuPlot, je tape load 'Euler.p' et j'obtiens une magnifique courbe !

Bibliographie

La programmation FORTRAN

Fortran 77 de Ficheux-Vapne
Numerical Recipes in Fortran: The Art of Scientific Computing de William H. Press, et al
Programmer en Fortran 90. Guide complet de Claude Delannoy
Fortran 90/95 for Scientists and Engineers de Stephen J. Chapman
Computer Applications in Physics With Fortran And Basic de Suresh Chandra

Gnuplot

le site www.gnuplot.info