Visu in situ avec l’outil de visualisation VisIt

Simulation de la convection de Rayleigh-Bénard par Visit

Simulation de la convection Rayleigh-Bénard à l'aide de Visit

 

Equipe

Billon Florian
Cassadour Martin
Darricau Baptiste
Mignot Benjamin

 

Encadrant

Neau Hervé

Pedrono Annaïg

 

 


L'objectif de cette page est création d'un tutoriel pour la visualisation in situ à l'aide du logiciel Visit. On aborde l'écriture du code fortran (les subroutines à utiliser, forme du programme principale) pour envoyer les données que l'on veut observer à Visit, sa compilation en créant un Makefile associé et de l'utilisation de l'interface Visit pour tracer les champs voulus.


               

Introduction

Qu'est  ce que VisIt ? Qu'est in Situ ?

VisIt est un logiciel de visualisation et d'animation interactif de données. Outre le fait d'être un simple logiciel de visualisation, il est particulièrement adapté pour faire de la visualisation in situ.

Le principe du in Situ est de permettre à un logiciel de simulation de récupérer les données,et ce  pendant que le code tourne, sans toutefois enregistrer afin de les visualiser en temps réel. De plus il est possible de commander le code à travers le logiciel, ici VisIt, et de faire avancer ou stopper la routine de calcul du code afin d'avoir un bon contrôle de la simulation.

Il a notamment l'avantage de ne nécessiter aucune création de fichiers de données pour fonctionner, permettant ainsi de pouvoir faire tourner de gros calculs et d'afficher nos résultats sans occuper de mémoire, ni de nous faire perdre du temps avec l'enregistrement de ces fichiers.

De la même manière, il peut également modifier les paramètres de notre programme en temps réel, ce qui évite d'avoir à relancer le programme à chaque fois que l'on veut tester quelque chose de nouveau.

Ainsi le in Situ est très intéressant pour les calculs  particulièrement longs car on peut vérifier le bon déroulement du code en direct et sans stocker de données et modifier les paramètres du code en fonction de la visualisation.

Visualisation d'un modèle d'avion de chasse sur VisIt

Le tutoriel

L'utilisation du in Situ avec VisIt sur des codes fortran reste relativement méconnue. En effet il n'existe quasiment aucune documentation sur ce domaine. L'objectif de ce tutoriel est de vous guider pour adapter votre code Fortran afin de pouvoir l'utiliser en In Situ et le visualiser grâce à VisIt.

Attention : Il existe de nombreuses versions de VisIt. Ce tutoriel est basé sur VisIt 2.9.0. Rien n'est sûr quand à un fonctionnement du code fourni sur d'autres versions de VisIt.

Dans un premier temps nous allons vous présenter l'installation de VisIt et une prise en main rapide sur les fichiers exemples.

Ensuite nous vous expliquerons les modifications à faire sur votre code pour le faire fonctionner en in situ et le visualier grâce à VisIt. Pour cela nous utiliserons deux exemples: un code simple que nous avons créer qui résout une équation d'advection de la température dans une conduite puis un code en parallèle créer par Jean Favre, du CSCS, qui résout une équation de Laplace.

Enfin, nous vous indiquerons la démarche à suivre pour l'utilisation in situ de VisIt : la partie visualisation et la partie contrôle des paramètres.

 

Installation et prise en main de VisIt

Installation de VisIt

Afin d'installer le logiciel, nous détaillons ci-dessous la procédure à suivre :

- Télécharger le fichier compressé sur le site https://wci.llnl.gov/simulation/computer-codes/visit/ dans l'onglet Downloads puis Executables

- Utiliser l'outil gunzip (commande gzip option -d) dans le terminal et non pas dans l'éditeur de fichier. Nous avons ainsi modifié notre archive initiale qui est maintenant décompressible à l'aide de l'outil tar.

- Pour lancer maintenant VisIt, dans le dossier bin du dossier décompressé on va chercher l'éxecutable qui s'appelle tout simplement visit.

- Si vous voulez pouvoir l'ouvrir depuis n'importe où dans votre terminal, il va vous falloir créer un PATH. Pour cela, créer un fichier qui s'appelle .bashrc dans votre répertoire personnel dans lequel on va mettre la ligne suivante:

export PATH=$PATH:/.../visit2_9_0.linux-x86_64/bin

Avec le /.../ correspondant au chemin où se trouve le dossier d'installation de VisIt.

Prise en main de VisIt

Une fois l'installation de Visit terminée, il ne reste plus qu'à le prendre en main:

Dans le dossier d'installation de VisIt, il y a un dossier data dans lequel se trouvent plusieurs fichiers exemples. Pour vous expliquer comment utiliser VisIt, nous allons utiliser ces fichiers.

Ouvrir un fichier (par exemple curv3d.silo). Ensuite dans la zone Tracer, cliquez sur Ajouter et la vous pouvez voir tout une liste de tracés que VisIt peut afficher. Cette liste ne dépend pas du fichier sélectionné. A l'intérieur de ces tracés, se trouvent les données qui se trouvent dans le fichier qui sont affichables sous la forme choisie. Par exemple avec le fichier curv3d.silo, les données curvmesh3d ne peuvent s'afficher que sous forme de Maillage ou Subset mais pas sous forme de Pseudo couleur, contrairement à la pression p.

 

Maintenant pour afficher ces données cliquez sur Tracer ( il se peut qu'il soit cacher dans la barre Tracés, dans ce cas cliquez sur >>). On obtient par exemple, le tracé suivant:

 

Modification du code

Code simple

Pour pouvoir utiliser le code en VisIt in Situ, il est nécessaire de modifier le code et de rajouter un certain nombre de subroutines indispensables au fonctionnement de VisIt. La plupart de ces subroutines ont un code générique utilisable pour tous les codes en fortran.  De plus certaines subroutines sont en option, nécessaire suivant le type de maillage et suivant les options voulues pour commander le code depuis visIt.

Ajout de données au module

Il est nécessaire de rajouter un certain nombre de paramètres dont se sert VisIt et que nous utiliserons dans le programme principal, dans la boucle de temps  et dans certaines subroutines de Visit:

LOGICAL :: visUpdateFlag, visRunFlag
  INTEGER simUpdate, visitstate, blocking, visResult
  REAL*8 simtime

Ces paramètres peuvent être écrits pour plus de lisibilité dans un module.

 

Attention!!!

Les variables que vous voulez lire avec Visit doivent être imprétativement en real 8. Il faut le signaler lorsque vous les mettez dans le module. Toutefois vous pouvez aussi mettre une variable real 4 en real 8 grâce à la fonction dble() qui permet de doubler sa précision. Par exemple, vous vouler que la variable T soit égale à 2 en real 8 vous écriverez: T=dble(2)

Cependant, le maillage lu avec Visit doit être, quant à lui, en real 4, c'est-à-dire en variable real par défaut.

Si c'est deux recommendations ne sont pas suivies les valeurs affichées par Visit seront complétement fausses.

 

 

Subroutines obligatoires

 visitbroadcastintfunction et visitbroadcaststringfunction

Le code Fortran n'ayant pas d'outil de pointage comme en code C, il est nécessaire de rajouter ces deux subroutines qui permettent à VisIt d'avoir un pointeur et de pouvoir fonctionner en in situ.

Ces deux subroutines génériques sont obligatoires quelque soit le code.

 

c---------------------------------------------------------------------------
c visitbroadcastintfunction
c---------------------------------------------------------------------------
      integer function visitbroadcastintfunction(value, sender)
      implicit none
      integer value, sender
      visitbroadcastintfunction = 0
      end

c---------------------------------------------------------------------------
c visitbroadcaststringfunction
c---------------------------------------------------------------------------
      integer function visitbroadcaststringfunction(str, lstr, sender)
      implicit none
      character*8 str
      integer     lstr, sender
      visitbroadcaststringfunction = 0
      end

 

visitactivatetimestep et processvisitcommand

Ces deux subroutines permettent de faire avancer le code à partir de VisIt, ce qui est l'un des atout du In Situ.

 

c---------------------------------------------------------------------------
c visitactivatetimestep
c---------------------------------------------------------------------------
      integer function visitactivatetimestep()
      implicit none
      include "visitfortransimV2interface.inc"
      visitactivatetimestep = VISIT_OKAY
      end

c-----------------------------------------------------------------
c processvisitcommand
c-----------------------------------------------------------------

      integer function processvisitcommand()
      use m_type
      implicit none

      include "visitfortransimV2interface.inc"
      integer command, e, doloop, success, ret
      integer VISIT_COMMAND_PROCESS
      integer VISIT_COMMAND_SUCCESS
      integer VISIT_COMMAND_FAILURE
      parameter (VISIT_COMMAND_PROCESS = 0)
      parameter (VISIT_COMMAND_SUCCESS = 1)
      parameter (VISIT_COMMAND_FAILURE = 2)
      integer h, dl

      success = visitprocessenginecommand()

        if(success.gt.0) then
              command = VISIT_COMMAND_SUCCESS
              ret = 1
          else
              command = VISIT_COMMAND_FAILURE
              ret = 0
          endif

          
       processvisitcommand = ret
      end

 

Les subroutines obligatoires mais utiles seulement dans le cas d'un code en parallèle :

subroutine visitslaveprocesscallback , visitcommandcallback

Ces deux  subroutines ne sont utilisées que si le code est écrit en parallèle mais les autres subroutines font appel à elles quelque soit le type de code. visitslaveprocesscallback est donc une subroutine vide. Quant à visitcommandcallback il est toutefois nécessaire de la modifier si jamais on veut commander le pas de temps depuis VisIt.

c---------------------------------------------------------------------------
c visitslaveprocesscallback
c---------------------------------------------------------------------------
      subroutine visitslaveprocesscallback ()
      implicit none
      end

visitcommandcallback sans commande depuis visIt :

c---------------------------------------------------------------------------
c visitcommandcallback
c---------------------------------------------------------------------------
      subroutine visitcommandcallback (cmd, lcmd, args, largs)
      implicit none
      character*8 cmd, args
      integer     lcmd, largs
      end

visitcommandcallback avec commande depuis visIt :

c visitcommandcallback
c---------------------------------------------------------------------------
      subroutine visitcommandcallback(cmd, lcmd, args, largs)
      use m_type
      implicit none
      character*8 cmd, args
      integer     lcmd, largs, err
      include "visitfortransimV2interface.inc"

c     Handle the commands that we define in visitgetmetadata.
      if(visitstrcmp(cmd, lcmd, "halt", 4).eq.0) then
           visRunFlag = .false.
      elseif(visitstrcmp(cmd, lcmd, "step", 4).eq.0) then
           call simulate_one_timestep()
           err = visitupdateplots()
      elseif(visitstrcmp(cmd, lcmd, "run", 3).eq.0) then
          visRunFlag = .true.
      elseif(visitstrcmp(cmd, lcmd, "update", 6).eq.0) then
         err = visitupdateplots()
      endif
      end

 

 

 

 

Subroutine visitgetmetadata

La subroutine visitgetmetadata est la subroutine coeur pour VisIt. Il s'agit d'une subroutine boite où on va donner à VisIt tous les types d'information que l'on va lui fournir, sans toutefois les lui fournir dans cette subroutine. On va ainsi lui donner le type de maillage du code, le nombre de maillages, le nombre de variables et leur type. On peut aussi rajouter le type de matériau et les commandes de simulation  si on veut contrôler le code depuis VisIt.

La base de la subroutine :

Quelques soient les options que l'on veut fournir à VisIt,  les arguments (inexistants) et le début de la subroutine sont les mêmes.

c---------------------------------------------------------------------------
c visitgetmetadata
c---------------------------------------------------------------------------
      integer function visitgetmetadata()
      use m_type
      implicit none
      include "visitfortransimV2interface.inc"
c     Local variables
      integer md, mmd, vmd, cmd, err

      simtime = iter*0.1

      if(visitmdsimalloc(md).eq.VISIT_OKAY) then
          err = visitmdsimsetcycletime(md, iter, simtime)

      if( visRunFlag ) then
              err = visitmdsimsetmode(md, VISIT_SIMMODE_RUNNING)
          else
              err = visitmdsimsetmode(md, VISIT_SIMMODE_STOPPED)
          endif

Vient ensuite les différentes parties suivant les options que l'on veut :

Comme dit précédemment ces parties sont génériques, elles ne contiennent rien des éléments du code du programme principal

maillage 2D rectangulaire :

On peut toutefois modifier les unités en remplaçant ce qui est entre crochet sachant que le nombre à coté est le nombre de lettres.

c     Add a 2D rect mesh
          if(visitmdmeshalloc(mmd).eq.VISIT_OKAY) then
              err = visitmdmeshsetname(mmd, "mesh", 4)
              err = visitmdmeshsetmeshtype(mmd,
     .            VISIT_MESHTYPE_RECTILINEAR)
              err = visitmdmeshsettopologicaldim(mmd, 2)
              err = visitmdmeshsetspatialdim(mmd, 2)
              err = visitmdmeshsetxunits(mmd, "m", 1)
              err = visitmdmeshsetyunits(mmd, "m", 1)
              err = visitmdmeshsetxlabel(mmd, "Width", 5)
              err = visitmdmeshsetylabel(mmd, "Height", 6)
              err = visitmdsimaddmesh(md, mmd)
          endif

maillage 3D  :

c'est le même que le maillage 2D, il suffit juste de rajouter l'axe z dans les unités et les variables.

 if(visitmdmeshalloc(m2).eq.VISIT_OKAY) then
              err = visitmdmeshsetname(m2, "mesh3d", 6)
              err = visitmdmeshsetmeshtype(m2,
     .            VISIT_MESHTYPE_CURVILINEAR)
              err = visitmdmeshsettopologicaldim(m2, 3)
              err = visitmdmeshsetspatialdim(m2, 3)
              err = visitmdmeshsetxunits(m2, "cm", 2)
              err = visitmdmeshsetyunits(m2, "cm", 2)
              err = visitmdmeshsetzunits(m2, "cm", 2)
              err = visitmdmeshsetxlabel(m2, "Width", 5)
              err = visitmdmeshsetylabel(m2, "Height", 6)
              err = visitmdmeshsetzlabel(m2, "Depth", 5)
              err = visitmdmeshsetcellorigin(m2, 1)
              err = visitmdmeshsetnodeorigin(m2, 1)

              err = visitmdsimaddmesh(md, m2)
          endif

 

Variable sur les mailles en 2D :

c Add a zonal scalar variable on mesh2d.
          if(visitmdvaralloc(vmd).eq.VISIT_OKAY) then
              err = visitmdvarsetname(vmd, "zonal", 5)
              err = visitmdvarsetmeshname(vmd, "mesh2d", 6)
              err = visitmdvarsettype(vmd, VISIT_VARTYPE_SCALAR)
              err = visitmdvarsetcentering(vmd, VISIT_VARCENTERING_ZONE)

              err = visitmdsimaddvariable(md, vmd)
          endif

 

Variable sur les noeuds en 3D :

c Add a nodal scalar variable on mesh3d.
          if(visitmdvaralloc(vmd).eq.VISIT_OKAY) then
              err = visitmdvarsetname(vmd, "nodal", 5)
              err = visitmdvarsetmeshname(vmd, "mesh3d", 6)
              err = visitmdvarsettype(vmd, VISIT_VARTYPE_SCALAR)
              err = visitmdvarsetcentering(vmd, VISIT_VARCENTERING_NODE)

              err = visitmdsimaddvariable(md, vmd)
          endif

 

Affichage de courbes :

Add a curve variable
          if(visitmdcurvealloc(cmd).eq.VISIT_OKAY) then
              err = visitmdcurvesetname(cmd, "sine", 4)
              err = visitmdcurvesetxlabel(cmd, "angle", 5)
              err = visitmdcurvesetxunits(cmd, "radians", 7)
              err = visitmdcurvesetylabel(cmd, "amplitude", 9)

              err = visitmdsimaddcurve(md, cmd)
          endif

 

Ajout de matériau :

  Add a material
          if(visitmdmatalloc(mat).eq.VISIT_OKAY) then
              err = visitmdmatsetname(mat, "mat", 3)
              err = visitmdmatsetmeshname(mat, "mesh2d", 6)
              err = visitmdmataddmaterialname(mat, "Iron", 4)
              err = visitmdmataddmaterialname(mat, "Copper", 6)
              err = visitmdmataddmaterialname(mat, "Nickel", 6)

              err = visitmdsimaddmaterial(md, mat)
          endif

 

Ajout de commandes de simulation:

Si on veut commander le code depuis VisIt, et ainsi faire fonctionner le code pas à pas. Cette partie est nécessaire et très pratique, notamment sur les codes devants tourner pendant une longue période .

c     Add simulation commands
          if(visitmdcmdalloc(cmd).eq.VISIT_OKAY) then
              err = visitmdcmdsetname(cmd, "halt", 4)
              err = visitmdsimaddgenericcommand(md, cmd)
          endif
          if(visitmdcmdalloc(cmd).eq.VISIT_OKAY) then
              err = visitmdcmdsetname(cmd, "step", 4)
              err = visitmdsimaddgenericcommand(md, cmd)
          endif
          if(visitmdcmdalloc(cmd).eq.VISIT_OKAY) then
              err = visitmdcmdsetname(cmd, "run", 3)
              err = visitmdsimaddgenericcommand(md, cmd)
          endif
      endif
      visitgetmetadata = md
      end

 

 

Subroutines modifiables

Toutes les subroutines suivantes sont obligatoires, même si vous ne vous en servez pas, sinon le programme ne compilera pas. Pour les subroutines que vous n'utilisez pas, il suffit de mettre une commande error, qui s'affichera si vous cherchez à y accéder sur VisIt.

Les autres subroutines doivent être adaptées à votre code afin d'associer les bons tableaux aux bonnes variables.

 visitgetmesh

Cette subroutine permet de donner à VisIt le maillage, elle est donc particulièrement importante.

Pour un maillage 2D rectangulaire il suffit de donner le  nombre et  les vecteurs positions des noeuds suivant les différentes directions. Dans notre cas ces variables sont respectivement nx1 et ny1 pour les nombre et pnx et pny  pour les vecteurs position. Nous les avons placé dans le module m_type pour éviter de les appeler dans la subroutine.

integer function visitgetmesh(domain, name, lname)
     
      use m_type
      implicit none
      character*8 name
      integer     domain, lname
      include "visitfortransimV2interface.inc"

      integer h, x, y, err

      h = VISIT_INVALID_HANDLE
      if(visitstrcmp(name, lname, "mesh", 4).eq.0) then          
c Create a rect mesh here
        if(visitrectmeshalloc(h).eq.VISIT_OKAY) then
              err = visitvardataalloc(x)
              err = visitvardataalloc(y)
              err = visitvardatasetf(x,VISIT_OWNER_COPY, 1,
nx1, pnx)
              err = visitvardatasetf(y,VISIT_OWNER_COPY, 1,
ny1, pny)
              err = visitrectmeshsetcoordsxy(h, x, y)
          end if
      end if
      
      visitgetmesh = h

      
      end

 

visitgetvariable

Ici on associe la variable au tableau correspondant et à son nom. Pour nous le tableau est Temp et le nom de la variable est "Température". Comme pour la subroutine précédente les tableaux sont dans le module.

 

 integer function visitgetvariable(domain, name, lname)
      use m_type
      implicit none
      character*8 name
      integer     domain, lname, err
      include "visitfortransimV2interface.inc"
     
      integer h

      if(visitvardataalloc(h).eq.VISIT_OKAY) then
        if(visitstrcmp(name, lname, "
temperature", 11).eq.0) then
        err = visitvardatasetd(h, VISIT_OWNER_SIM, 1,
     .                         
nx1*ny1, Temp)
        endif
      endif

      visitgetvariable = h
      end

Le reste des subroutines n'a pas été utilisé dans notre programme, elles sont néanmoins obligatoires.  Les codes suivants permettent de faire tourner le programme. Toutefois si vous voulez utiliser une option supplémentaire, tel que l'ajout d'un matériau par exemple, vous devez modifier  la subroutine associée.

 

c---------------------------------------------------------------------------
c visitgetmixedvariable
c---------------------------------------------------------------------------
      integer function visitgetmixedvariable(domain, name, lname)
      implicit none
      character*8 name
      integer     domain, lname
      include "visitfortransimV2interface.inc"
      visitgetmixedvariable = VISIT_INVALID_HANDLE
      end

c---------------------------------------------------------------------------
c visitgetcurve
c---------------------------------------------------------------------------
      integer function visitgetcurve(name, lname)
      implicit none
      character*8 name
      integer     lname
      include "visitfortransimV2interface.inc"
      visitgetcurve = VISIT_INVALID_HANDLE
      end

c---------------------------------------------------------------------------
c visitgetdomainlist
c---------------------------------------------------------------------------
      integer function visitgetdomainlist(name, lname)
      implicit none
      character*8 name
      integer     lname
      include "visitfortransimV2interface.inc"
      visitgetdomainlist = VISIT_INVALID_HANDLE
      end

c---------------------------------------------------------------------------
c visitgetdomainbounds
c---------------------------------------------------------------------------
      integer function visitgetdomainbounds(name, lname)
      implicit none
      character*8 name
      integer     lname
      include "visitfortransimV2interface.inc"
      visitgetdomainbounds = VISIT_INVALID_HANDLE
      end

c---------------------------------------------------------------------------
c visitgetdomainnesting
c---------------------------------------------------------------------------
      integer function visitgetdomainnesting(name, lname)
      implicit none
      character*8 name
      integer     lname
      include "visitfortransimV2interface.inc"
      visitgetdomainnesting = VISIT_INVALID_HANDLE
      end

c---------------------------------------------------------------------------
c visitgetmaterial
c---------------------------------------------------------------------------
      integer function visitgetmaterial(domain, name, lname)
      implicit none
      character*8 name
      integer     domain, lname
      include "visitfortransimV2interface.inc"
      visitgetmaterial = VISIT_INVALID_HANDLE
      end

 

 

 

Modification du programme principal

Il faut adapter le code du programme principal afin qu'il fasse appel aux différentes subroutines utilisées par VisIt.

Pour cela, il faut tout d'abord faire appel à la bibliothèque VisIt, ainsi qu'à la commande. Il faut donc mettre le code suivant au début du code (après implicite none).

include "visitfortransimV2interface.inc"
  external processvisitcommand
  integer processvisitcommand

Vient ensuite la définition des variables. Une seule variable est à rajouter au premier code : la variable erreur : integer err.

Avant la l'attribution et l'allocation des différents paramètres, il est nécessaire de faire appel à la subroutine simulationarguments. Cela se fait sous la forme suivante :

call simulationarguments()
      err = visitsetupenv()
      err = visitinitializesim("prog", 4, &
     "Interface to prog", 17, &
     "/no/useful/path", 15, &
     VISIT_F77NULLSTRING, VISIT_F77NULLSTRINGLEN, &
     VISIT_F77NULLSTRING, VISIT_F77NULLSTRINGLEN, &
     VISIT_F77NULLSTRING, VISIT_F77NULLSTRINGLEN)

Cette subroutine est fournie avec VisIt dans le fichier simulationexemplev2.f qu'il ne faudra pas oublier de rajouter dans le makefile.

Toutefois si vous voulez l'intégrer à votre code, la voici :

c-----------------------------------------------------------------
c simulationarguments : The routine handles command line arguments
c-----------------------------------------------------------------
      subroutine simulationarguments()
      implicit none
      character (len=80) str
      integer err, i, N, len
      integer visitsetoptions, visitsetdirectory, visitopentracefile
      N = iargc()
      i = 1
      len = 80
5     if (i.le.N) then
          call getarg(i, str)
          if(str.eq."-dir") then
              call getarg(i+1, str)
              err = visitsetdirectory(str, len)
              i = i + 1
          elseif(str.eq."-options") then
              call getarg(i+1, str)
              err = visitsetoptions(str, len)
              i = i + 1
          elseif(str.eq."-trace") then
              call getarg(i+1, str)
              err = visitopentracefile(str, len)
              i = i + 1
          endif
          i = i + 1
          goto 5
      endif
      end

    

 

 

La suite du programme n'est pas modifié jusqu'à la fin de l'initialisation des différents paramètres du programme.

Viens ensuite la boucle principale de Visit qui va servir à connecter Visit et le code. Nous avons décider de faire une seule subroutine traitant toute cette boucle. Toutefois il s'agit de la solution que nous avons choisi, ce n'est pas l'unique solution.

Il suffit donc d'appeler cette subroutine : call mainloop() et le programme est terminé.

Subroutine mainloop :

Pour tous les programmes et subroutines qui vont faire appel à des subroutines de VisIt, il est nécessaire de faire appel à la bibiothèque VisIt. Le mainloop ne fait pas exception. Comme poour le programme principal il est nécessaire d'y faire appel avant le début du programme :

include "visitfortransimV2interface.inc"
  external processvisitcommand
  integer processvisitcommand

Initialisons tout d'abord le paramètre visRunFlag

visRunFlag = .FALSE.

Dans la boucle de temps (while), il est nécessaire d'informer VisIt de la possibilité de faire un pas de temps supplémentaire. De plus la dernière ligne est l'attente d'un signal de VisIt pour simuler une étape.

if(visRunFlag) then
    blocking = 0
  else
    blocking = 1
  endif

 visitstate = visitdetectinput(blocking, -1)

En cas d'input de la part de VisIt, on simule une étape :

  if (visitstate < 0) then
    print *, 'visitState < 0 in time stepping loop. Ignoring and continuing'
 
  elseif (visitstate == 0) then
call simulate_one_timestep()

elseif (visitstate == 1) then
    visRunFlag = .TRUE.
    visResult = visitattemptconnection()
    if (visResult == 1) then
       write (6,*) 'VisIt connected!'
    else
       write (6,*) 'VisIt did not connect!'
    endif

  elseif (visitstate == 2) then
    visRunFlag =.FALSE.
     if (processvisitcommand() == 0) then
      visResult = visitdisconnect()
      visRunFlag =.TRUE.
    endif
  endif

Cette boucle permet donc de connecter Visit avant de faire avancer le pas de temps, de faire attendre Visit une fois qu'il est connecté et de le déconnecter à le fin du calcul.

La subroutine simulate_one_timestep() appelé dans la boucle est composée du calcul et d'une avancé du pas de temps, qu'il faut signaler à VisIt. Notons que l'appel à la bibliothèque et la définition de err est ici aussi obligatoire. Une fois le pas de temps effectué, les dernières lignes de la subroutine simulate_one_timestep() sont :

iter=iter+1
err = visittimestepchanged()
err = visitupdateplots()

Ces lignes permettent de faire varier le tracé visualisé sous Visit.

 

Création du Makefile

Dans la création du Makefile, il faut ajouter deux items.

Le premier est include qui s'écrit sous la forme suivante:

INCLUDE = -I/.../visit/v2.9/visit2_9_0.linux-x86_64/2.9.0/linux-x86_64/libsim/V2/include

Le second est la librairie qui sert à trouver toutes les fonction et subroutines de visit. Il faut en plus lui ajouter une option de compilation. Cet item s'écrit sous la forme suivante:

LIB = -L/.../visit/v2.9/visit2_9_0.linux-x86_64/2.9.0/linux-x86_64/libsim/V2/lib -lsimV2 -lsimV2f \
-lstdc++ -Wl,--no-as-needed -ldl

Pour le reste, le Makefile reste inchangé mais il ne faut pas oublié de rajouter les includes et la librairie dans la création des objets comme dans l'exemple suivant:

mainloop.o:        mainloop.f90
$$    $(FC) $(OPT) $(INCLUDE) $(LIB) mainloop.f90 -c

 

Code en parallèle

Comme exemple de code en parallèle, nous prenons l'exemple d'un code de Jean Favre du CSCS (Swiss Center for Scientific Computing) qui permet la résolution d'une équation de Laplace.

Dans cette partie nous n'allons pas détailler les fonctions utilisées mais plutôt comment compiler ce code pour pouvoir l'utiliser ensuite avec VisIt.

Ces codes peuvent se trouver ici: http://calcul.math.cnrs.fr/spip.php?article173 tout en bas sur le lien TP.

Tout d'abord, comme indiqué dans le readme.txt, modifier le fichier CMakeLists.txt pour y mettre les chemins correspondants à votre installation de VisIt.

Ensuite, afin de rendre ce code compatible avec la nouvelle version de VisIt in situ, modifier le fichier PJacobi_InSitu.f qui se trouve dans le dossier F90 et y intégrer les lignes suivantes:

 

c---------------------------------------------------------------------------

c visitgetmixedvariable

c---------------------------------------------------------------------------

integer function visitgetmixedvariable(domain, name, lname)

implicit none

character*8 name

integer domain, lname

include "visitfortransimV2interface.inc"

visitgetmixedvariable = VISIT_INVALID_HANDLE

end

 

c---------------------------------------------------------------------------

 

Ensuite ouvrez un terminal et aller dans le dossier F90, d'ici lancer la commande:

ccmake -DCMAKE_Fortran_COMPILER=mpif90 -DCMAKE_EXE_LINKER_FLAGS="-lstdc++"

Appuyer sur [c] pour configurer chaque élément puis sur [g] pour générer les fichiers Makefile. Maintenant un simple make devrait suffire à compiler le code et à générer deux exécutables: pjacobi et pjacobi_visit. le premier correspond au programme seul qui résout l'équation de Laplace et le deuxième lui correspond au programme qui peut se connecter à VisIt.

Maintenant, lancer le programme qui peut se connecter à VisIt en utilisant la commande ./pjacobi_visit ou encore mpiexec -n 4 ./pjacobi_visit si vous voulez utiliser les quatre processeurs et passer à l'étape suivante.

 

 

Exemple de code

Afin de donner un aperçu d'un code complet utilisant VisIt en in situ, nous mettons ci-dessous le lien pour télécharger notre exemple de code.

/travaux/projnum/sites/default/files/users/fbillon/code_VisIt_InSitu.zip

Il fonctionne pour la version 2.9.0 de VisIt.

Ce code traite l'écoulement d'un fluide dans une conduite. On peut mettre au choix de l'utilisateur :

    - le type de transfert : advection, diffusion ou les deux en même temps.
    - le nombre de maille selon les x et les y
    - le choix du maillage: uniforme ou non uniforme
 

Utilisation de VisIt in situ

Visualisation

Maintenant que vous avez lancé votre programme, celui-ci aura créer un fichier à la racine de votre disque dans un dossier caché .visit.

Ouvrez VisIt et cliquer sur ouvrir et allez chercher ce dossier .visit dans lequel se trouve le dossier simulations. A l'intérieur de celui-ci, ouvrez le dernier fichier en .sim2.

Si tout s'est bien passé, sur le terminal ou vous avez lancer le programme devrait s'afficher la ligne:

VisIt connected !.

Vous pouvez donc aller sur VisIt et faire un tracer de ce que vous voulez comme indiquer dans la première partie de ce tutoriel. Normalement, votre tracé devrait s'afficher. Toutefois si cous demandez un tracé alors que vous n'avez pas codé l'option dans le code, un message d'erreur s'affichera.

Vous avez sûrement remarqué que votre programme ne fait rien, il est en attente. En réalité il est connecté à VisIt et attend ces commandes. Pour ce faire, sur VisIt, allez dans Fichier -> Simulations. A partir de là, vous avez un onglet Controls dans lequel se trouvent plusieurs commandes. Pour pouvez donc faire un step pour faire marcher votre calcul pas à pas, run pour le faire calculer en continu et halt pour l'arrêter.

 

 

 

Conclusion

A la base, ce projet était très ambitieux, la prise en main, la compréhension et l'application dans des cas concrets de VisIt in situ s'annonçait comme une énorme tache.

Il y a eu en effet de nombreux obstacles à notre avancement, notamment des obstacles sur lesquels nous n'avions pas les compétences pour les franchir. Grâce à l'aide de nos professeurs encadrant, nous avons pu arrivé à bout de ce projet, et le résultat est très gratifiant.

Nous avons néanmoins appris beaucoup avec ce projet. Tout d'abord la prise en main du logiciel VisIt, qui ne nous était pas familier. Ensuite l'apprentissage d'une nouvelle bibliothèque Fortran et surtout son implémentation dans un de nos code.

Enfin, ce tutoriel n'est qu'une initiation à VisIt in situ. Nous n'avons en effet utiliser qu'une petite partie des subroutines dont dispose la bibliothèque VisIt in situ, celles nécessaires à l'affichage d'un résultat. Il en reste de nombreuses autres qui permettent d'autres type d'affichage de résultats et il reste notamment la partie commande des paramètres avec VisIt in situ. C'est une partie très importante dont nous n'avons pas eu le temps de découvrir. Il reste donc bien du travail pour ceux qui souhaiteraient continuer ce projet.