Télécharger eqex.eso

Retour à la liste

eqex
C EQEX      SOURCE    CB215821  24/04/12    21:15:48     11897          
      SUBROUTINE EQEX
C***********************************************************************
C VERSION    : ????
C HISTORIQUE : 22/03/00: gounand
C Rajout des préconditionneurs ILUT (ILU with dual truncation) et d'une
C variante (ILUT2) qui remplit mieux la mémoire et des paramètres
C associés : ILUTLFIL (ILUT level of fill) et ILUTDTOL (ILUT drop
C tolerance)
C HISTORIQUE : 20/12/99: gounand
C Ajout des indices 'TYRENU' (type de renumérotation) et 'PCMLAG'
C (placement des multiplicateurs de Lagrange à la table d'indice
C 'METHINV'.
C HISTORIQUE : 08/04/04 : ajout ILUTP
C HISTORIQUE : 27/10/10: JCARDO: correction bug lié à IARG, qui pouvait
C                                parfois valoir zéro à tort (label 110)
C HISTORIQUE :
C***********************************************************************
      IMPLICIT REAL*8 (A-H,O-Z)
      IMPLICIT INTEGER (I-N)
      PARAMETER (NBM=21,NBL=4,NOPT=75)
      PARAMETER (NTB=2,NBH=7)
 
-INC PPARAM
-INC CCOPTIO
-INC CCNOYAU
-INC SMLREEL
-INC SMLENTI
      POINTEUR MLENT4.MLENTI
-INC SMLMOTS
      POINTEUR MINCO.MLMOTS
-INC SMELEME
      POINTEUR IGEOM.MELEME
-INC SMMODEL
-INC SMCHPOI
      POINTEUR MCHINI.MCHPOI
 
      LOGICAL XEQUA,TTRAN,TPROJ,LOG1
      INTEGER RESTRT
      INTEGER IPST
      CHARACTER*(LOCOMP) LMOTS(NBM),NOM,CHAI,MTYP,NOMI,NOML,NMACO
      CHARACTER*(LONOM) NOMZ
      CHARACTER*8 LSCHE(NBH)
      CHARACTER*20 NOMO,MEQUA,MNEFMD
      CHARACTER*8 TYPE,TYPC,TYPS
      CHARACTER*8 LOPTI(NOPT)
      CHARACTER*8 MOIMP(NBL)
      DIMENSION   KINCD(100)
      CHARACTER*8 TINCD(100)
 
      CHARACTER*8 LTAB(NTB)
      DIMENSION KTAB(NTB)
 
      DATA LMOTS /'ZONE    ','OPER    ','INCO    ','CLIM    ',
     &     'ITMA    ','ALFA    ','DTI     ','IIMP    ',
     &     'DUMP    ','OPTI    ','NOMVI   ','DOMINC  ',
     &     'TPSI    ','TFINAL  ','FIDT    ','NISTO   ',
     &     'NITER   ','OMEGA   ','EPS     ','IMPR    ',
     &     'EQUA    '/
      DATA MOIMP /'UIMP    ','VIMP    ','WIMP    ','TIMP    '/
 
      DATA LOPTI /
C indice KFORM 0 à 3
     &     'EFM1    ','EF      ','VF      ','EFMC    ','????????',
     &     'LINE    ','MACRO   ','QUAF    ','LINB    ','ISOQ    ',
C indice IDCEN
     &     'CENTREE ','SUPGDC  ','SUPG    ','TVISQUEU','CNG     ',
     &     'PSI     ','JOHNSON ','UPWIND  ','GODUNOV ','VANLEER ',
     &     'VLH     ','HUSVL   ','HUSVLH  ','AUSM    ','CG      ',
     &     'VSM     ','VSMCC   ','SUPGDCH ','SUPGH   ','????????',
C indice KPOIN
     &     'SOMMET  ','FACE    ','CENTRE  ','CENTREP0','CENTREP1',
     &     'MSOMMET ','????????','????????','????????','????????',
C indice KIMPL ---------------------->|indice ISCHT
     &     'IMPL    ','EXPL    ','SEMI    ','BDF2    ','BDF4    ',
     &     'DIV2    ','CMD     ','RIGIDITE','LIMITE  ','NODIV   ',
C indice IKOMP-------------->|    RNG        KMACO      ALE
     &     'CONS    ','NOCONS  ','CONS2   ','RNG     ','ALE     ',
C indice MTRMASS---------------------->|
     &     'MMPLEINE','MMDIAGO ','MMPG    ','MATCONS ','????????',
C indice IDEUL  ------------------------->|
     &     'EULER   ','EULERMS ','EULERMST','????????','????????',
C indice KPOIND->|
     &     'INCOD   ','INCOP   ','STABP   ','MUCONS  ','FTAU    ',
C
     &     'MUVARI  ','????????','????????','????????','????????'/
 
      DATA LSCHE /'EUL_EXPL','EUL_IMPL','TVISQ   ','SEMI    ',
     &            'CN      ','CNG     ','BDF2    '/
 
      DATA LTAB/'DOMAINE ','EQEX    '/
C***
C     WRITE(IOIMP,*) ' DEBUT EQEX '
      IPST=0
      NBIK=0
      IDP=0
      MTABD=0
      CALL INITI(KTAB,NTB,0)
 
C Définition des options par défaut
      CALL CRTABL(KOPT)
      CALL ECMM(KOPT,'SOUSTYPE','KOPT')
      CALL ECME(KOPT,'IDCEN',2)
      CALL ECME(KOPT,'RNG   ',1)
      CALL ECME(KOPT,'IKOMP',0)
      CALL ECME(KOPT,'KMACO',0)
      CALL ECMM(KOPT,'NMACO','xxxxxxxx')
      CALL ECME(KOPT,'KIMPL',0)
      CALL ECME(KOPT,'KFORM',0)
      CALL ECMF(KOPT,'AIMPL',1.D0)
      CALL ECME(KOPT,'ALE',0)
      CALL ECME(KOPT,'KMU',0)
      CALL ECME(KOPT,'KPOIND',99)
      CALL ECME(KOPT,'KPOIN',2)
      CALL ECME(KOPT,'MTRMASS ',1)
      CALL ECME(KOPT,'IDEUL   ',1)
      CALL ECME(KOPT,'ISCHT',0)
      CALL ECME(KOPT,'IDIV',0)
      CALL ECMF(KOPT,'CMD',0.2D0)
      CALL ECMF(KOPT,'STABP',1.D-2)
      CALL ECME(KOPT,'RIGIDITE',0)
      CALL ECME(KOPT,'LIMITE',0)
c     CALL ECMM(KOPT,'INEFMD','xxxxxxxx')
 
 
      CALL QUETYP(TYPE,0,IRET)
 
         MATABL=0
         MMODEL=0
      IF(TYPE.EQ.'MMODEL')THEN
         CALL LIROBJ('MMODEL',MMODEL,0,IRET)
         CALL LEKMOD(MMODEL,MTBLE,INEFMD)
         IF(MTBLE.EQ.0)RETURN
         KTAB(1)=MTBLE
         KTAB(2)=0
         CALL ECMM(KOPT,'INEFMD',LOPTI(5+INEFMD))
 
      ELSEIF(TYPE.EQ.'MOT')THEN
C Nouvelle directive EQUA
c        write(6,*)' Nouvelle directive EQUA'
         CALL LIRCHA(CHAI,1,LCHAR)
            KTAB(1)=-1
            KTAB(2)=0
           IF(CHAI(1:4).NE.'EQUA')THEN
            CALL ECRCHA(CHAI)
            GO TO 6
           ENDIF
 
         CALL CRTABL(MATABL)
         CALL ECMM(MATABL,'SOUSTYPE','EQEX')
         XEQUA=.TRUE.
         CALL ECML(MATABL,'XEQUA',XEQUA)
         CALL ECML(MATABL,'XRIG',.FALSE.)
C Lecture du nom de l'equation/inconnue
         NBIC=0
         JGN=4
         JGM=0
         NINCT=0
         SEGINI MLMOT2
         CALL ECMO(MATABL,'LISTINCO','LISTMOTS',MLMOT2)
         SEGDES MLMOT2
 3       CONTINUE
         CALL QUETYP(MTYP,0,IRET)
         IF(MTYP.EQ.'MMODEL  ')THEN
          GO TO 4
         ELSEIF(MTYP.EQ.'MOT     ')THEN
          CALL LIRCHA(CHAI,1,LCHAR)
          IF(CHAI.EQ.'RIGIDITE')THEN
          CALL ECML(MATABL,'XRIG',.TRUE.)
          GO TO 3
          ENDIF
 
          NBIC=NBIC+1
          SEGACT MLMOT2
          JGM=NBIC
          SEGADJ MLMOT2
          MLMOT2.MOTS(NBIC)=CHAI
          SEGDES MLMOT2
          GO TO 3
         ELSE
C On ne trouve pas d'objet de type %m1:8
          MOTERR( 1: 8) = '        '
          CALL ERREUR(38)
          RETURN
         ENDIF
 
 
 4       CONTINUE
         IF(NBIC.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
          MOTERR( 1: 8) = 'MOT     '
          CALL ERREUR(38)
          RETURN
         ENDIF
 
         CALL LIROBJ('MMODEL',MMODEL,0,IRET)
         CALL ECMO(MATABL,'MODELE','MMODEL  ',MMODEL)
          SEGACT MMODEL
*       Détermination de MACRO et INEFMD
          IMODEL = KMODEL(1)
          SEGACT IMODEL
        IF(NEFMOD.GE.129.AND.NEFMOD.LE.135)THEN
          INEFMD=1
        ELSEIF(NEFMOD.GE.136.AND.NEFMOD.LE.141)THEN
          INEFMD=2
        ELSEIF(NEFMOD.GE.143.AND.NEFMOD.LE.149)THEN
          INEFMD=3
        ELSEIF(NEFMOD.GE.158.AND.NEFMOD.LE.164)THEN
          INEFMD=4
        ELSE
C% Le type d'élément fini %m1:8 ne convient pas.
          WRITE(NOM,FMT='(I8)')NEFMOD
          MOTERR( 1: 8) = NOM
          CALL ERREUR(926)
          RETURN
        ENDIF
         SEGDES MMODEL,IMODEL
         CALL ECMM(KOPT,'INEFMD',LOPTI(5+INEFMD))
         CALL ECMM(MATABL,'INEFMD',LOPTI(5+INEFMD))
 
         JGN=8
         JGM=NBIC
         NINCT=0
         SEGINI MLMOT3
         CALL ECMO(MATABL,'TYPEINCO','LISTMOTS',MLMOT3)
         DO 5 I=1,NBIC
          CALL LIRCHA(CHAI,1,LCHAR)
           IF(LCHAR.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT     '
            CALL ERREUR(38)
            RETURN
           ENDIF
 
           CALL ECMM(MATABL,'TYPPRESS','MSOMMET ')
          IF(CHAI.EQ.'PRESSION')THEN
           CALL LIRCHA(NOM,1,LCHAR)
           CALL OPTLI(IPRE,LOPTI(33),NOM,4)
            IF(IPRE.EQ.0)THEN
C% On a lu : %m1:8 : , alors qu'on attend un des mots clés suivant :
C%  %m9:16  %m17:24  %m25:32  %m33:40
             MOTERR( 1: 8) = NOM
             MOTERR( 9:16) = LOPTI(33)
             MOTERR(17:24) = LOPTI(35)
             MOTERR(25:32) = LOPTI(36)
             CALL ERREUR(930)
             RETURN
            ENDIF
           CALL ECMM(MATABL,'TYPPRESS',NOM)
          ENDIF
 
          IF(CHAI.NE.'TEMPERAT'.AND.CHAI.NE.'VITESSE'.AND.
     &       CHAI.NE.'PRESSION')THEN
C Option %m1:8 incompatible avec les données
             MOTERR( 1: 8) =CHAI
             CALL ERREUR(803)
             RETURN
           ENDIF
          MLMOT3.MOTS(I)=CHAI
 5    CONTINUE
         SEGDES MLMOT3
 
C Lecture du Schema en temps
 
          CALL LIRCHA(CHAI,1,LCHAR)
           IF(LCHAR.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT     '
            CALL ERREUR(38)
            RETURN
           ENDIF
 
        IF(CHAI.NE.'PERM'.AND.CHAI.NE.'TRAN'.AND.CHAI.NE.'PROJ')THEN
C% On a lu : %m1:8 : , alors qu'on attend un des mots clés suivant :
C%  %m9:16  %m17:24  %m25:32  %m33:40
          MOTERR( 1: 8) = CHAI
          MOTERR( 9:16) = 'PERM'
          MOTERR(17:24) = 'TRAN'
          MOTERR(25:32) = 'PROJ'
          CALL ERREUR(930)
          RETURN
        ENDIF
 
        TPROJ=.FALSE.
        IF(CHAI.EQ.'PERM')TTRAN=.FALSE.
        IF(CHAI.EQ.'TRAN')TTRAN=.TRUE.
        IF(CHAI.EQ.'PROJ')THEN
        TTRAN=.TRUE.
        TPROJ=.TRUE.
        ENDIF
        CALL ECML(MATABL,'TRAN',TTRAN)
        CALL ECML(MATABL,'PROJ',TPROJ)
 
      IP=0
 
      IF(TTRAN)THEN
C Lecture des parametres du transitoire
      CALL QUETYP(MTYP,0,IRET)
           IF(IRET.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT ou R'
            CALL ERREUR(38)
            RETURN
           ENDIF
           IF(MTYP.EQ.'MOT')THEN
           CALL LIRCHA(CHAI,1,LCHAR)
           CALL ECMM(MATABL,'DELTAT',CHAI)
           ELSEIF(MTYP.EQ.'FLOTTANT')THEN
           CALL LIRREE(XVAL,1,IRET)
           CALL ECMF(MATABL,'DELTAT',XVAL)
           ELSE
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT ou R'
            CALL ERREUR(38)
            RETURN
           ENDIF
 
        CALL LIRCHA(NOM,1,IRET)
           IF(IRET.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT     '
            CALL ERREUR(38)
            RETURN
           ENDIF
 
      CALL OPTLI(IPST,LSCHE,NOM,NBH)
 
      IF(IPST.EQ.0)THEN
         WRITE(IOIMP,*)'Directive : ',NOM
         WRITE(IOIMP,*)'non trouvée dans la liste ->',LSCHE
         RETURN
      ENDIF
      CALL ECMM(MATABL,'SCHEMAT',LSCHE(IPST))
      IF(IPST.EQ.4)THEN
        CALL LIRREE(XVAL,1,IRET)
          IF(IRET.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
           MOTERR( 1: 8) = 'FLOTTANT'
           CALL ERREUR(38)
           RETURN
          ENDIF
        CALL ECMF(MATABL,'Betat',XVAL)
      ENDIF
 
      CALL ECML(MATABL,'XDIAG',.FALSE.)
 
      ENDIF
 
      DO 51 I=1,NBIC
 
        CALL LIRCHA(CHAI,1,IRET)
           IF(IRET.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT     '
            CALL ERREUR(38)
            RETURN
           ENDIF
 
           CALL OPTLI(IP,LMOTS,CHAI,NBM)
           IF(IP.NE.0)THEN
           write(6,*)' OPTI et ZONE sont des mots cles '
           write(6,*)' Choix mal venu pour un nom de variable'
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT     '
            CALL ERREUR(38)
            RETURN
           ENDIF
 
      WRITE(NOM,FMT='(A3,I1)')'INC',I
      CALL ECMM(MATABL,NOM,CHAI)
 
 51   CONTINUE
 
       IF(IPST.EQ.7.AND.TTRAN)THEN
        DO 52 I=1,NBIC
 
        CALL LIRCHA(CHAI,1,IRET)
           IF(IRET.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
            MOTERR( 1: 8) = 'MOT     '
            CALL ERREUR(38)
            RETURN
           ENDIF
 
        WRITE(NOM,FMT='(A3,I1)')'IMC',I
c     write(6,*)' NOM,CHAI=',NOM,CHAI
        CALL ECMM(MATABL,NOM,CHAI)
 
 52     CONTINUE
       ENDIF
 
 
      ELSEIF(TYPE.EQ.'TABLE')THEN
         CALL LIROBJ(TYPE,MTBLE,0,IRET)
         TYPC=' '
         CALL ACMO(MTBLE,'SOUSTYPE',TYPC,IRET)
         IF(TYPC.EQ.'MOT     ')THEN
          CALL ACMM(MTBLE,'SOUSTYPE',TYPS)
          IF(TYPS.EQ.'DOMAINE')THEN
          KTAB(1)=MTBLE
          KTAB(2)=0
          XEQUA=.FALSE.
          CALL ECML(MTBLE,'XEQUA',XEQUA)
          CALL ECMM(MTBLE,'INEFMD','        ')
          CALL ECMM(KOPT,'INEFMD','        ')
          ELSEIF(TYPS.EQ.'EQEX')THEN
          KTAB(1)=0
          KTAB(2)=MTBLE
          CALL ACML(MTBLE,'XEQUA',XEQUA)
          CALL ACMM(MTBLE,'INEFMD',MNEFMD)
          CALL ECMM(KOPT,'INEFMD',MNEFMD)
          ELSE
          WRITE(IOIMP,*)' On attend une table soustype DOMAINE ou EQEX'
          RETURN
          ENDIF
         ENDIF
      ELSE
          KTAB(1)=-1
          KTAB(2)=0
      ENDIF
 
 6    CONTINUE
 
      IF(KTAB(2).NE.0)THEN
         MTABLE=KTAB(2)
         TYPE='LISTMOTS'
         CALL ACMO(MTABLE,'LISTOPER',TYPE,MLMOT1)
         SEGACT MLMOT1
         NEQUA=MLMOT1.MOTS(/2)
      ELSEIF(KTAB(1).NE.0)THEN
         MTABD=KTAB(1)
 
         IF(MATABL.EQ.0)THEN
         CALL CRTABL(MTABLE)
         CALL ECMM(MTABLE,'SOUSTYPE','EQEX')
         XEQUA=.FALSE.
         CALL ECML(MTABLE,'XEQUA',XEQUA)
         CALL ECMM(MTABLE,'INEFMD','        ')
         CALL ECMM(KOPT,'INEFMD',MNEFMD)
         ELSE
         MTABLE=MATABL
         ENDIF
         CALL CRTABL(MINCO)
         CALL ECMM(MINCO,'SOUSTYPE','INCO')
         NEQUA=0
         CALL ECME(MTABLE,'DISCPRES',0)
         IF(KTAB(1).GT.0)THEN
         IF(MMODEL.EQ.0)CALL ECMO(MTABLE,'DOMAINE','TABLE   ',MTABD)
         CALL ECME(MTABLE,'NAVISTOK',MMODEL)
         ELSE
         CALL ECME(MTABLE,'NAVISTOK',-1)
         ENDIF
C?       IF(MMODEL.NE.0)THEN
C        CALL ECMO(MTABLE,'DOMAINE','MMODEL  ',MMODEL)
C        CALL ECMO(MTABLE,'TDOMAINE','TABLE   ',MTABD)
C?       CALL ECMO(MTABLE,'DOMAINE','TABLE   ',MTABD)
C?       ENDIF
         CALL ECMO(MTABLE,'INCO','TABLE   ',MINCO)
         CALL ECMM(MTABLE,'NOMVI','UN')
         JGN=8
         JGM=0
         SEGINI MLMOT1
         CALL ECMO(MTABLE,'LISTOPER','LISTMOTS',MLMOT1)
         CALL ECME(MTABLE,'ITMA',0)
         CALL ECMF(MTABLE,'ALFA',1.D0)
         CALL ECME(MTABLE,'FIDT',20)
         CALL ECME(MTABLE,'NISTO',20)
         CALL ECME(MTABLE,'NITER',1)
         CALL ECME(MTABLE,'IPP',0)
         CALL ECME(MTABLE,'IMPR',0)
         CALL ECMF(MTABLE,'OMEGA',1.D0)
         CALL ECMF(MTABLE,'EPS',1.D-2)
         CALL ECMF(MTABLE,'TFINAL',1.D30)
         NAT=2
         NSOUPO=0
         SEGINI MCHPOI
         JATTRI(1)=2
         CALL ECMO(MTABLE,'CLIM','CHPOINT',MCHPOI)
 
         CALL CRTABL(MTABT)
         CALL ECMM(MTABT,'SOUSTYPE','PASDETPS')
         CALL ECMO(MTABLE,'PASDETPS','TABLE   ',MTABT)
         IPAT=1
         CALL ECME(MTABT,'NUPASDT',IPAT)
         DT=1.D30
         CALL ECMF(MTABT,'DELTAT',DT)
         CALL ECMF(MTABT,'DELTAT-1',DT)
         TPS=0.D0
         CALL ECMF(MTABT,'TPS',TPS)
         CALL ECMM(MTABT,'OPER','EQEX')
         CALL ECMM(MTABT,'ZONE','EQEX')
         CALL ECMF(MTABT,'DTCONV',0.D0)
         CALL ECMF(MTABT,'DTDIFU',0.D0)
         CALL ECMF(MTABT,'DIAEL',0.D0)
         CALL ECME(MTABT,'NUEL',0)
 
C Définition de la méthode d'inversion et des paramètres
C éventuels associés
         CALL CRTABL(MTINV)
         CALL ECMM(MTINV,'SOUSTYPE','METHINV')
         CALL ECMO(MTABLE,'METHINV ','TABLE   ',MTINV)
C     Méthode d'inversion du système
C                   1 : résolution directe (Choleski)
C                   2 : Gradient Conjugué
C                   3 : Bi-Gradient Conjugué Stabilisé (BiCGSTAB)
C                   4 : BiCGSTAB(2)
         KTYPI=1
         CALL ECME(MTINV,'TYPINV',KTYPI)
C     Niveau d'impression pour la partie résolution itérative
         IMPINV=0
         CALL ECME(MTINV,'IMPINV',IMPINV)
C     Options spécifiques aux méthodes itératives :
C
C        - Pour l'assemblage : type de renumérotation
C        * 'RIEN' : pas de renumérotation
C        * 'SLOA' : algorithme de chez Sloan
C        * 'GIPR' : Gibbs-King (profile reduction)
C        * 'GIBA' : Gibbs-Poole-Stockmeyer (bandwidth reduction)
         CALL ECMM(MTINV,'TYRENU','SLOA')
*         CALL ECMM(MTINV,'TYRENU','RIEN')
C        - Pour l'assemblage : prise en compte des mult.lag
C        * 'RIEN'
C        * 'APR2'
         CALL ECMM(MTINV,'PCMLAG','APR2')
C      - Pour l'assemblage : SCALING (type ENTIER) :
C         Scaling de la matrice :
C         - 0 : pas de scaling
C         - 1 : scaling par les normes L2 des lignes et des colonnes
C         Par défaut : 0
         ISCAL=0
         CALL ECME(MTINV,'SCALING',ISCAL)
C      - Pour l'assemblage : OUBMAT (type ENTIER) :
C         Oublie les matrices élémentaires :
C         - 0 : non
C         - 1 : oui
C         Par défaut : 0
         IOUBL=0
         CALL ECME(MTINV,'OUBMAT',IOUBL)
C        - Champoint d'initialisation de la méthode
C            (i.e. estimation de l'inconnue)
*****         MCHINI=0
         NAT=0
         NSOUPO=0
         SEGINI MCHINI
         SEGDES MCHINI
         CALL ECMO(MTINV,'XINIT','CHPOINT ',MCHINI)
C        - Nombre maxi d'itérations à effectuer
         ITER=2000
         CALL ECME(MTINV,'NITMAX',ITER)
C        - Norme maxi (L2 normé par le second membre) du résidu
         RESID=1.D-10
         CALL ECMF(MTINV,'RESID',RESID)
C        - Type de préconditionnement :
C                   0 : pas de préconditionnement
C                   1 : préconditionnement par la diagonale
C                   2 : préconditionnement D-ILU
C                   3 : préconditionnement ILU(0) (Choleski)
C                   4 : préconditionnement MILU(0) (Choleski modifié)
C                   5 : préconditionnement ILUT (dual truncation)
C                   6 : préconditionnement ILUT2 (une variante du
C                       précédent qui remplit mieux la mémoire et
C                       fonctionne mieux quelquefois)
C                   7 : préconditionnement ILUTP (avec pivoting)
C                   8 : préconditionnement ILUTPG (avec pivoting)
C                       ILUTP version gounand
C                       On traite de manière spéciale les termes
C                       qui sont dans ILU(0)
C                   9 : préconditionnement ILUTPG2 (avec pivoting)
C                       ILUTP version gounand 2
C                       On garde tous les termes qui sont dans ILU(0)
         KPREC=3
         CALL ECME(MTINV,'PRECOND',KPREC)
C        - Pour une méthode ILUT, on a les deux indices suivant :
C          * ILUTLFIL : encombrement maximal (approximatif) du
C                       préconditionneur, par rapport à la matrice.
C          * ILUTDTOL : "drop tolerance" pour le préconditionneur.
C                       i.e. en-dessous de cette valeur relative, les
C                       termes de la factorisation incomplète seront
C                       oubliés.
         XLFIL=2.D0
         CALL ECMF(MTINV,'ILUTLFIL',XLFIL)
* -1. sinon, oubli possible des 0.D0 dans le préconditionneur
         XDTOL=-1.D0
         CALL ECMF(MTINV,'ILUTDTOL',XDTOL)
C        - Pour une méthode ILUTP, on a les deux indices suivant :
C        * ILUTPPIV (type REEL) (compris entre 0.D0 et 1.D0) :
C           0.D0 : on ne pivote pas
C           1.D0 : on pivote tout le temps
C           (recommandation : entre 0.1D0 et 0.01D0)
C           Par défaut : 0.1D0
         XSPIV=0.1D0
         CALL ECMF(MTINV,'ILUTPPIV',XSPIV)
C        - Fréquence de recalcul du préconditionneur en fonction
C          des deux indices de boucle suivant :
C          * indice de boucle sur les pas de temps
C          * indice de boucle sur la boucle d'itérations utilisée
C            pour résoudre les non-linéarités
C        Par défaut, on recalcule tout le temps le préconditionneur
         IFCPRT=1
         IFCPRI=1
         CALL ECME(MTINV,'FCPRECT',IFCPRT)
         CALL ECME(MTINV,'FCPRECI',IFCPRI)
C        - 'Breakdown tolerance' pour les méthodes de type
C           BiCGSTAB. Si un certain produit scalaire de vecteurs
C           "direction" est inférieur à cette tolérance, la
C           méthode s'arrete.
         BRTOL=1.D-40
         CALL ECMF(MTINV,'BCGSBTOL',BRTOL)
C        - Paramètre de relaxation pour le préconditionnement
C          MILU(0) compris entre 0. et 1.
C          S'il est égal à 0, on se ramène à ILU(0)
C          S'il est égal à 1, MILU(0) est dit non relaxé
         RXMILU=1.D0
         CALL ECMF(MTINV,'MILURELX',RXMILU)
C        - Paramètre de redémarrage pour GMRES(m)
         RESTRT=50
         CALL ECME(MTINV,'GMRESTRT',RESTRT)
      ENDIF
 
 1    CONTINUE
      CALL LIRCHA(NOM,0,IRET)
      IF(IRET.EQ.0)GO TO 90
 2    CONTINUE
 
      IF(IDP.NE.0)WRITE(IOIMP,*) ' Directive en cours :',NOM
C                 WRITE(IOIMP,*) ' Directive en cours :',NOM
      CALL OPTLI(IP,LMOTS,NOM,NBM)
C     WRITE(IOIMP,*)' EQEX, IP=',ip,' NOM=',nom
      IF(IP.EQ.0)THEN
         WRITE(IOIMP,*)'Directive : ',NOM
         WRITE(IOIMP,*)'non trouvée dans la liste ->',LMOTS
         RETURN
      ENDIF
 
 
      GO TO (10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,
     &       25,26,27,28,29),IP
 
C ZONE
 
 10   CONTINUE
 
      CALL LITABS(LTAB,KTAB,1,0,IRET)
 
         MMODEL=0
      IF(KTAB(1).EQ.0)THEN
         CALL LIROBJ('MMODEL',MMODEL,0,IRET2)
         IF(IRET2.EQ.0)THEN
         WRITE(IOIMP,*)' On attend un objet TABLE DOMAINE ou MODELE'
         RETURN
         ENDIF
         CALL LEKMOD(MMODEL,MTBLE,INEFMD)
         IF(MTBLE.EQ.0)RETURN
         KTAB(1)=MTBLE
         CALL ECMM(KOPT,'INEFMD',LOPTI(5+INEFMD))
      ENDIF
 
      CALL QUENOM(NOMZ)
      GO TO 1
 
C OPER
 
 11   CONTINUE
 
      CALL LIRCHA(NOMO,1,LNOMO)
      IF(LNOMO.EQ.0)THEN
         WRITE(IOIMP,*)' ON ATTEND LE NOM DE L OPERATEUR'
         RETURN
      ENDIF
 
* ECRITURE DU NOM DE L'OPERATEUR
 
      NEQUA=NEQUA+1
      IF(NEQUA.LT.10)THEN
         LNOMOT=LNOMO+1
         WRITE(MEQUA,FMT='(I1,19A1)')NEQUA,(NOMO(I:I),I=1,LNOMO)
C     WRITE(IOIMP,*)' MEQUA=',MEQUA
      ELSEIF(NEQUA.LT.100.AND.NEQUA.GE.10)THEN
         LNOMOT=LNOMO+2
         WRITE(MEQUA,FMT='(I2,18A1)')NEQUA,(NOMO(I:I),I=1,LNOMO)
      ELSEIF(NEQUA.LT.1000.AND.NEQUA.GE.100)THEN
         LNOMOT=LNOMO+3
         WRITE(MEQUA,FMT='(I3,17A1)')NEQUA,(NOMO(I:I),I=1,LNOMO)
C     WRITE(IOIMP,*)' MEQUA=',MEQUA
      ELSE
         WRITE(IOIMP,*)'PLUS DE 999 OPERATEURS : CAS NON PREVU'
         RETURN
      ENDIF
      JGN=8
      JGM=MLMOT1.MOTS(/2)+1
      SEGADJ MLMOT1
      MLMOT1.MOTS(JGM)=NOMO(1:LNOMO)
C     CALL LENCHA(MEQUA,LC1)
      CALL CRTABL(MTABX)
      CALL ECMM(MTABX,'SOUSTYPE','KIZX')
C     CALL ECMO(MTABLE,MEQUA(1:8),'TABLE',MTABX)
      CALL ECCTAB(MTABLE,'MOT',0,0.D0,MEQUA(1:LNOMOT),.TRUE.,0,
     &     'TABLE',0,0.D0,CHAI,.TRUE.,MTABX)
* ECRITURE DE LA TABLE DE REFERENCE
      CALL ECMO(MTABX,'EQEX','TABLE',MTABLE)
* ECRITURE DU NOM DE LA ZONE
      CALL ECMM(MTABX,'NOMZONE',NOMZ)
* ECRITURE DE MELEMZ
      CALL ECMO(MTABX,'DOMZ','TABLE',KTAB(1))
      IF(MMODEL.NE.0)THEN
      CALL ECMO(MTABX,'DOMZ','MMODEL',MMODEL)
      CALL ECMO(MTABX,'TDOMZ','TABLE',KTAB(1))
      ENDIF
 
      CALL ECMM(MTABX,'NOMOPER',NOMO(1:8))
      CALL ECMO(MTABX,'KOPT','TABLE',KOPT)
* ECRITURE DE LA LISTE DES ARGUMENTS
* 1) on initialise la variable IARG à 0
* 2) chaque fois que l'on trouve un argument pour l'opérateur courant,
*    on incrémente cette variable et on boucle (=> GOTO 110)
* 3) on met à jour 'IARG' dans MTABX dès qu'il n'y a plus d'argument:
*           - soit il n'y a plus d'objet passé à EQEX (=> GOTO 90)
*           - soit on est tombé sur un autre mot-clé  (=> GOTO 2 )
      IARG=0
 
 110  CONTINUE
      CALL QUETYP(MTYP,0,IRET)
C     WRITE(IOIMP,*)' MTYP=',mtyp,' iret=',iret
      IF(IRET.EQ.0)THEN
*       PLUS AUCUN MOT DANS EQEX => ON MET À JOUR 'IARG'
        CALL ECME(MTABX,'IARG',IARG)
        GO TO 90
      ENDIF
 
      IF(MTYP.EQ.'MOT     ')THEN
 
         CALL LIRCHA(NOM,1,IRET)
         CALL OPTLI(IP,LMOTS,NOM,NBM)
         IF(IP.EQ.0)THEN
            IARG=IARG+1
            IF(IARG.GT.9)CALL ARRET(0)
            CHAI=NOM
            WRITE(NOM,FMT='(A3,I1)')'ARG',IARG
            CALL ECMM(MTABX,NOM(1:4),CHAI)
            GO TO 110
         ELSE
*           MOT-CLÉ TROUVÉ => PLUS D'ARGUMENTS => ON MET À JOUR 'IARG'
            CALL ECME(MTABX,'IARG',IARG)
C     WRITE(IOIMP,*)' 1er gt2 nom=',nom
            GO TO 2
         ENDIF
 
      ELSEIF(MTYP.EQ.'FLOTTANT')THEN
         CALL LIRREE(XVAL,1,IRET)
         IARG=IARG+1
         IF(IARG.GT.9)CALL ARRET(0)
         WRITE(NOM,FMT='(A3,I1)')'ARG',IARG
         CALL ECMF(MTABX,NOM(1:4),XVAL)
         GO TO 110
 
      ELSEIF(MTYP.EQ.'ENTIER')THEN
         CALL LIRENT(IENT,1,IRET)
         IARG=IARG+1
         IF(IARG.GT.9)CALL ARRET(0)
         WRITE(NOM,FMT='(A3,I1)')'ARG',IARG
CC       XVAL=DBLE(IENT)
CC       CALL ECMF(MTABX,NOM(1:4),XVAL)
         CALL ECME(MTABX,NOM(1:4),IENT)
         GO TO 110
 
      ELSEIF(MTYP.EQ.'POINT')THEN
         CALL LIROBJ('POINT',IZTAB,1,IRET)
         IARG=IARG+1
         IF(IARG.GT.9)CALL ARRET(0)
         WRITE(NOM,FMT='(A3,I1)')'ARG',IARG
         CALL ECMO(MTABX,NOM(1:4),'POINT',IZTAB)
         GO TO 110
 
      ELSEIF(MTYP.EQ.'LOGIQUE ')THEN
         CALL LIRLOG(LOG1,1,IRET)
         IARG=IARG+1
         IF(IARG.GT.9)CALL ARRET(0)
         WRITE(NOM,FMT='(A3,I1)')'ARG',IARG
         CALL ECML(MTABX,NOM(1:4),LOG1)
         GO TO 110
 
      ELSE
         CALL LIROBJ(MTYP,IZTAB,1,IRET)
         IARG=IARG+1
         IF(IARG.GT.9)RETURN
         WRITE(NOM,FMT='(A3,I1)')'ARG',IARG
         CALL ECMO(MTABX,NOM(1:4),MTYP,IZTAB)
         GO TO 110
 
      ENDIF
 
C INCO
 
 12   CONTINUE
C On crée si ce n'est pas déja fait la liste totale des inconnues
C et on la complète
      TYPE=' '
      CALL ACMO(MTABLE,'LISTINCO',TYPE,MLMOT2)
      IF(TYPE.NE.'LISTMOTS')THEN
         JGN=4
         JGM=0
         NINCT=0
         SEGINI MLMOT2
         CALL ECMO(MTABLE,'LISTINCO','LISTMOTS',MLMOT2)
         SEGDES MLMOT2
      ENDIF
      SEGACT MLMOT2
      NINCT=MLMOT2.MOTS(/2)
      NINCT0=0
      SEGDES MLMOT2
 
C on crée la liste des inconnues associées à l'opérateur
      JGN=4
      JGM=0
      JG =0
      SEGINI MLMOTS,MLENT4
* ECRITURE DE LA LISTE DES INCONNUES
      CALL ECMO(MTABX,'LISTINCO','LISTMOTS',MLMOTS)
      CALL ECMO(MTABX,'NUMEINCO','LISTENTI',MLENT4)
      SEGDES MLMOTS,MLENT4
C on crée la liste des types d'inconnues associés à l'opérateur
      JGN=8
      JGM=0
      SEGINI MLMOT4
* ECRITURE DE LA LISTE DES TYPES D'INCONNUES
      CALL ECMO(MTABX,'TYPEINCO','LISTMOTS',MLMOT4)
      SEGDES MLMOT4
 
 120  CONTINUE
      SEGDES MLMOTS,MLMOT2
      CALL LIRCHA(NOM,0,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
 
      CALL OPTLI(IP,LMOTS,NOM,NBM)
 
      IF(IP.EQ.2.OR.IP.EQ.3)THEN
         WRITE(IOIMP,*)' Il faut recommencer a la directive ZONE '
         RETURN
      ENDIF
C     WRITE(IOIMP,*)' EQEX : ',NOM,IP
 
      IF(IP.EQ.0)THEN
 
         NINCT0=NINCT0+1
         IF(NINCT0.GT.5)THEN
            WRITE(IOIMP,*)' Opérateur EQEX :'
            WRITE(IOIMP,*)' Le nombre d''inconnues semble important ',
     &           NINCT0,' ? ',NOM
         ENDIF
 
         JGN=4
         SEGACT MLMOTS
         JGM=MOTS(/2)+1
         SEGADJ MLMOTS
C      WRITE(IOIMP,*)' EQEX : ',NOM,' NINCT=',ninct
         MOTS(JGM)=NOM
 
C Cas directive EQUA
       IF(XEQUA)THEN
      TYPE=' '
      CALL ACMO(MTABLE,'TYPEINCO',TYPE,MLMOT3)
      IF(TYPE.NE.'LISTMOTS')THEN
      write(6,*)' Petit probleme Non prevu '
      RETURN
      ENDIF
        SEGACT MLMOT2,MLMOT3,MLMOT4,MLENT4
        NBIC=MLMOT2.MOTS(/2)
        DO 122 I=1,NBIC
 
         IF(NOM.EQ.MLMOT2.MOTS(I))THEN
          JGN=8
          JG=JGM
          SEGADJ MLMOT4,MLENT4
          MLMOT4.MOTS(JGM)=MLMOT3.MOTS(I)
          MLENT4.LECT(JGM)=I
          GO TO 123
         ENDIF
 122    CONTINUE
C% L'inconnue : %m1:8 : n'apparait pas dans la liste des inconnues.
        MOTERR( 1: 8) = NOM
        CALL ERREUR(931)
        RETURN
 
 123    CONTINUE
        SEGDES MLMOT2,MLMOT3,MLMOT4,MLENT4
 
 
       ENDIF
C On ecrit aussi directement dans MTABX NBINCO,INC1 INC2
C etc comme pour les arguments
C?    NBINCO=JGM
C?    WRITE(CHAI,FMT='(A3,I1)')'INC',NBINCO
C?    CALL ECMM(MTABX,CHAI(1:4),NOM)
C?    CALL ECME(MTABX,'NBINCO',NBINCO)
 
         SEGACT MLMOT2
         DO 121 I=1,NINCT
C     WRITE(IOIMP,*)' On cherche : ',MLMOT2.MOTS(I),NOM
            IF(MLMOT2.MOTS(I).EQ.NOM)GO TO 120
 121     CONTINUE
         JGM=NINCT+1
         NINCT=NINCT+1
 
         JGN=4
         SEGADJ MLMOT2
         MLMOT2.MOTS(JGM)=NOM
         GO TO 120
      ELSE
         SEGDES MLMOTS,MLMOT2
C     WRITE(IOIMP,*)' 2eme gt2 nom=',nom
         GO TO 2
      ENDIF
 
C CLIM
 
 13   CONTINUE
 
      TYPE=' '
      CALL ACMO(MTABLE,'CLIM',TYPE,MCHPO1)
      IF(TYPE.NE.'CHPOINT')MCHPO1=0
 
      CALL LIRCHA(NOMI,0,IRET)
C     WRITE(IOIMP,*)' MCHPO1=',mchpo1,' NOMI=',nomi
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
      CALL OPTLI(IP,LMOTS,NOMI,NBM)
      IF(IP.NE.0)THEN
         NOM=NOMI
C     WRITE(IOIMP,*)' 3eme gt2 nom=',nom
         GO TO 2
      ENDIF
 
      CALL LENCHA(NOMI,LCI)
 
      CALL LIRMOT(MOIMP,NBL,IP,1)
      IF(IP.EQ.0)THEN
         WRITE(IOIMP,*)' Directive CLIM : '
         WRITE(IOIMP,*)' On attend un mot cle de la liste suivante :'
     $        ,MOIMP
         RETURN
      ENDIF
 
      WRITE(NOML,FMT='(I1,A4)')IP,NOMI
      IF(IP.EQ.4)NOML=NOMI
      CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
 
      CALL ECRCHA('POI1')
      CALL ECROBJ('MAILLAGE',MELEME)
      CALL PRCHAN
      CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
      SEGACT MELEME
      N=NUM(/2)
      NAT=2
      NSOUPO=1
      NC=1
      SEGINI MCHPOI,MSOUPO,MPOVAL
      IFOPOI=IFOUR
      MOCHDE=TITREE
      MTYPOI='CLIM'
      JATTRI(1)=2
      IPCHP(1)=MSOUPO
C     WRITE(IOIMP,*)' EQEX 1 : IPCHP,mchpoi=',IPCHP(1),mchpoi
      NOCOMP(1)=NOML(1:4)
      IGEOC=MELEME
      IPOVAL=MPOVAL
 
      CALL QUETYP(MTYP,1,IRET)
 
      if(mtyp.eq.'TABLE')THEN
         ikkt=0
 1234    continue
         ikkt=ikkt+1
         WRITE(IOIMP,*)' ikkt=',ikkt
         if(ikkt.gt.100)return
         CALL LIROBJ('TABLE',IP,1,IRET)
         WRITE(IOIMP,*)' Petit incident '
         if(ip.ne.ktab(2))then
            WRITE(IOIMP,*)' Gros incident ',ip,ktab
         endif
         CALL QUETYP(MTYP,0,IRET)
         if(mtyp.eq.'TABLE')go to 1234
      endif
 
C     WRITE(IOIMP,*)' MTYP a=',mtyp,iret,' N=',N
 
      IF(MTYP.EQ.'FLOTTANT')THEN
         CALL LIRREE(XVAL,1,IRET)
         CALL INITD(VPOCHA,N,XVAL)
         SEGDES MPOVAL
      ELSEIF(MTYP.EQ.'ENTIER')THEN
         CALL LIRENT(IENT,1,IRET)
         XVAL=DBLE(IENT)
         CALL INITD(VPOCHA,N,XVAL)
         SEGDES MPOVAL
      ELSEIF(MTYP.EQ.'CHPOINT')THEN
         CALL LIROBJ('CHPOINT',MCHPO2,1,IRET)
         SEGACT MCHPO2
         NSP=MCHPO2.IPCHP(/1)
         CALL KRIPAD(MELEME,MLENTI)
         SEGACT MELEME
         DO 3569 L=1,NSP
            MSOUP2=MCHPO2.IPCHP(L)
            SEGACT MSOUP2
            IGEOM=MSOUP2.IGEOC
            MPOVA2=MSOUP2.IPOVAL
            SEGACT IGEOM,MPOVA2
            NBEL=IGEOM.NUM(/2)
            DO 3568 I=1,NBEL
               I1=IGEOM.NUM(1,I)
               II1=LECT(I1)
               IF(II1.EQ.0)GO TO 3568
               VPOCHA(II1,1)=MPOVA2.VPOCHA(I,1)
 3568       CONTINUE
            SEGDES MSOUP2,IGEOM,MPOVA2
 3569    CONTINUE
         SEGSUP MLENTI
         SEGDES MCHPO2
 
      ELSE
         WRITE(IOIMP,*)' TYPE NON ENCORE TRAITE'
         RETURN
      ENDIF
      SEGDES MCHPOI,MSOUPO,MPOVAL
 
      IF(MCHPO1.NE.0)THEN
         CALL ECROBJ('CHPOINT',MCHPOI)
         CALL ECROBJ('CHPOINT',MCHPO1)
         CALL PRFUSE
         CALL LIROBJ('CHPOINT',MCHPOI,1,IRET)
      ENDIF
      CALL ECMO(MTABLE,'CLIM','CHPOINT',MCHPOI)
 
      GO TO 13
 
C ITMA
 
 14   CONTINUE
      CALL LIRENT(IENT,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE ITMA (Nb maximum de pas de temps) :'
         WRITE(IOIMP,*)' On attend un entier '
         RETURN
      ENDIF
      CALL ECME(MTABLE,'ITMA',IENT)
      GO TO 1
 
C ALFA
 
 15   CONTINUE
      CALL LIRREE(XVAL,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE ALFA (Tolerance sur le pas de temps) :'
         WRITE(IOIMP,*)' doit etre compris entre 0 et 1 (1 par defaut)'
         WRITE(IOIMP,*)' On attend un flottant '
         RETURN
      ENDIF
      CALL ECMF(MTABLE,'ALFA',XVAL)
      GO TO 1
 
C DTI
 
 16   CONTINUE
      CALL LIRREE(XVAL,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE DTI (Pas de temps iinitial) :'
         WRITE(IOIMP,*)' On attend un flottant '
         RETURN
      ENDIF
         DT=XVAL
         CALL ECMF(MTABT,'DELTAT',DT)
         CALL ECMF(MTABT,'DELTAT-1',DT)
      GO TO 1
 
C IIMP
 
 17   CONTINUE
 
      TYPE=' '
      CALL ACMO(MTABLE,'IIMP',TYPE,MCHPO1)
      IF(TYPE.NE.'CHPOINT')MCHPO1=0
 
      CALL LIRCHA(NOMI,0,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
      CALL OPTLI(IP,LMOTS,NOMI,NBM)
      IF(IP.NE.0)THEN
         NOM=NOMI
C     WRITE(IOIMP,*)' 4eme gt2 nom=',nom
         GO TO 2
      ENDIF
 
      CALL LENCHA(NOMI,LCI)
 
      CALL LIRMOT(MOIMP,NBL,IP,1)
      IF(IP.EQ.0)THEN
         WRITE(IOIMP,*)' Directive IIMP : '
         WRITE(IOIMP,*)' On attend un mot cle de la liste suivante :'
     $        ,MOIMP
         RETURN
      ENDIF
C     WRITE(IOIMP,*)' MOIMP=',moimp(ip)
 
      WRITE(NOML,FMT='(I1,A4)')IP,NOMI
      IF(IP.EQ.4)NOML=NOMI
 
      CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
 
      CALL ECRCHA('POI1')
      CALL ECROBJ('MAILLAGE',MELEME)
      CALL PRCHAN
      CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
      SEGACT MELEME
      N=NUM(/2)
      NAT=2
      NSOUPO=1
      NC=1
      SEGINI MCHPOI,MSOUPO,MPOVAL
      IFOPOI=IFOUR
      MOCHDE=TITREE
      MTYPOI='IIMP'
      JATTRI(1)=2
      IPCHP(1)=MSOUPO
C     WRITE(IOIMP,*)' EQEX 2 : IPCHP,mchpoi=',IPCHP(1),mchpoi
      NOCOMP(1)=NOML(1:4)
      IGEOC=MELEME
      IPOVAL=MPOVAL
 
      CALL QUETYP(MTYP,1,IRET)
C     WRITE(IOIMP,*)' MTYP=',mtyp
 
      IF(MTYP.EQ.'FLOTTANT')THEN
         CALL LIRREE(XVAL,1,IRET)
         CALL INITD(VPOCHA,N,XVAL)
         SEGDES MPOVAL
      ELSEIF(MTYP.EQ.'ENTIER')THEN
         CALL LIRENT(IENT,1,IRET)
         XVAL=DBLE(IENT)
         CALL INITD(VPOCHA,N,XVAL)
         SEGDES MPOVAL
      ELSE
         WRITE(IOIMP,*)' TYPE NON ENCORE TRAITE'
         RETURN
      ENDIF
      SEGDES MCHPOI,MSOUPO,MPOVAL
 
      IF(MCHPO1.NE.0)THEN
         CALL ECROBJ('CHPOINT',MCHPOI)
         CALL ECROBJ('CHPOINT',MCHPO1)
         CALL PRFUSE
         CALL LIROBJ('CHPOINT',MCHPOI,1,IRET)
      ENDIF
      CALL ECMO(MTABLE,'IIMP','CHPOINT',MCHPOI)
 
      GO TO 17
 
C DUMP
 
 18   CONTINUE
      IDP=1
      GO TO 1
 
C OPTI
 
 19   CONTINUE
C Définition des options par défaut
      CALL CRTABL(KOPT1)
      CALL ECMM(KOPT1,'SOUSTYPE','KOPT')
 
      CALL ACME(KOPT,'IDCEN',IDCEN)
      CALL ACME(KOPT,'RNG  ',KRNG)
      CALL ACME(KOPT,'IKOMP',IKOMP)
      CALL ACME(KOPT,'KMACO',KMACO)
      CALL ACMM(KOPT,'NMACO',NMACO)
      CALL ACME(KOPT,'KIMPL',KIMPL)
      CALL ACME(KOPT,'KFORM',KFORM)
      CALL ACMF(KOPT,'AIMPL',AIMPL)
      CALL ACME(KOPT,'ALE  ',KALE)
      CALL ACME(KOPT,'KMU',KMU)
      CALL ACME(KOPT,'KPOIND',KPOIND)
      CALL ACME(KOPT,'KPOIN',KPOIN)
      CALL ACME(KOPT,'MTRMASS ',MTRMAS)
      CALL ACME(KOPT,'IDEUL   ',IDEUL)
      CALL ACME(KOPT,'ISCHT',ISCHT)
      CALL ACME(KOPT,'IDIV',IDIV)
      CALL ACMF(KOPT,'CMD',CMD)
      CALL ACMF(KOPT,'STABP',STAB)
      CALL ACME(KOPT,'RIGIDITE',IRIG)
      CALL ACME(KOPT,'LIMITE',LIMITE)
      CALL ACMM(KOPT,'INEFMD',MNEFMD)
 
      CALL ECME(KOPT1,'IDCEN',IDCEN)
      CALL ECME(KOPT1,'RNG   ',KRNG)
      CALL ECME(KOPT1,'IKOMP',IKOMP)
      CALL ECME(KOPT1,'KMACO',KMACO)
      CALL ECMM(KOPT1,'NMACO',NMACO)
      CALL ECME(KOPT1,'KIMPL',KIMPL)
      CALL ECME(KOPT1,'KFORM',KFORM)
      IF(KIMPL.EQ.1)AIMPL=1.D0
      IF(KIMPL.EQ.0)AIMPL=0.D0
      CALL ECMF(KOPT1,'AIMPL',AIMPL)
      CALL ECME(KOPT1,'ALE',KALE)
      CALL ECME(KOPT1,'KMU',KMU)
C?    CALL ECME(KOPT1,'KPOIND',99)
      CALL ECME(KOPT1,'KPOIND',KPOIND)
      CALL ECME(KOPT1,'KPOIN',KPOIN)
      CALL ECME(KOPT1,'MTRMASS ',MTRMAS)
      CALL ECME(KOPT1,'IDEUL   ',IDEUL)
      CALL ECME(KOPT1,'ISCHT',ISCHT)
      CALL ECME(KOPT1,'IDIV',IDIV)
      CALL ECMF(KOPT1,'CMD',CMD)
      CALL ECMF(KOPT1,'STABP',STAB)
      CALL ECME(KOPT1,'RIGIDITE',IRIG)
      CALL ECME(KOPT1,'LIMITE',LIMITE)
      CALL ECMM(KOPT1,'INEFMD',MNEFMD)
 
      KOPT=KOPT1
 
 191  CONTINUE
      CALL LIRCHA(NOM,0,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
      CALL OPTLI(IP,LOPTI,NOM,NOPT)
      IF(IP.EQ.0)THEN
         CALL ECRCHA(NOM)
         GO TO 1
      ENDIF
C     write(6,*)' NOM=',NOM,IP
 
      GO TO (1901,1902,1903,1904,1905,1906,1907,1908,1909,1910,
     &     1911,1912,1913,1914,1915,1916,1917,1918,1919,1920,
     &     1921,1922,1923,1924,1925,1926,1927,1928,1929,1930,
     &     1931,1932,1933,1934,1935,1936,1937,1938,1939,1940,
     &     1941,1942,1943,1944,1945,1946,1947,1948,1949,1950,
     &     1951,1952,1953,1954,1955,1956,1957,1958,1959,1960,
     &     1961,1962,1963,1964,1965,1966,1967,1968,1969,1970,
     &     1971,1972,1973,1974,1975                         ),IP
 
C Formulation EFM1
 1901 CALL ECME(KOPT,'KFORM',0)
      GO TO 191
C Formulation EF
 1902 CALL ECME(KOPT,'KFORM',1)
      GO TO 191
C Formulation VF
 1903 CALL ECME(KOPT,'KFORM',2)
      GO TO 191
C Formulation EFMC
 1904 CALL ECME(KOPT,'KFORM',3)
      GO TO 191
 
C Emplacements libres pour une nouvelle formulation
 1905 CONTINUE
      GO TO 191
 
C Formulation EF LINE
 1906 CALL ECMM(KOPT,'INEFMD','LINE')
      GO TO 191
C Formulation EF MACRO
 1907 CALL ECMM(KOPT,'INEFMD','MACRO')
      GO TO 191
C Formulation EF QUAF
 1908 CALL ECMM(KOPT,'INEFMD','QUAF')
      GO TO 191
C Formulation EF LINB
 1909 CALL ECMM(KOPT,'INEFMD','LINB')
      GO TO 191
C Formulation EF LINB
 1910 CALL ECMM(KOPT,'INEFMD','ISOQ')
      GO TO 191
 
C CENTREE
 1911 CALL ECME(KOPT,'IDCEN',1)
      GO TO 191
C SUPGDC
 1912 CALL ECME(KOPT,'IDCEN',2)
      GO TO 191
C SUPG
 1913 CALL ECME(KOPT,'IDCEN',3)
      GO TO 191
C Tenseur visqueux
 1914 CALL ECME(KOPT,'IDCEN',4)
      GO TO 191
C Crank Nicholson généralisé
 1915 CALL ECME(KOPT,'IDCEN',5)
      GO TO 191
C PSI
 1916 CALL ECME(KOPT,'IDCEN',6)
      GO TO 191
C JOHNSON
 1917 CALL ECME(KOPT,'IDCEN',7)
      GO TO 191
C UPWIND
 1918 CALL ECME(KOPT,'IDCEN',8)
      GO TO 191
C GODUNOV
 1919 CALL ECME(KOPT,'IDCEN',9)
      GO TO 191
C VANLEER
 1920 CALL ECME(KOPT,'IDCEN',10)
      GO TO 191
C VLH (VAN LEER - HANEL)
 1921 CALL ECME(KOPT,'IDCEN',11)
      GO TO 191
C HUSVL (VAN LEER + OSHER)
 1922 CALL ECME(KOPT,'IDCEN',12)
      GOTO 191
C HUSVLH (VAN LEER - HANEL + OSHER)
 1923 CALL ECME(KOPT,'IDCEN',13)
      GOTO 191
C AUSM (AUSM+)
 1924 CALL ECME(KOPT,'IDCEN',14)
      GOTO 191
C CG Colella-Glaz
 1925 CALL ECME(KOPT,'IDCEN',15)
      GOTO 191
C VSM  Viscosité de sous-maille
 1926 CALL ECME(KOPT,'IDCEN',16)
      GOTO 191
C VSMCC  Viscosité de sous-maille Capture de choc
 1927 CALL ECME(KOPT,'IDCEN',17)
      GOTO 191
C SUPGDCH
 1928 CALL ECME(KOPT,'IDCEN',18)
      GOTO 191
C SUPGH
 1929 CALL ECME(KOPT,'IDCEN',19)
      GOTO 191
C emplacements libres pour nouveaux schema
 1930 CONTINUE
      GO TO 191
 
C sommet
 1931 CALL ECME(KOPT,'KPOIN',0)
      GO TO 191
C face
 1932 CALL ECME(KOPT,'KPOIN',1)
      GO TO 191
C centre
 1933 CALL ECME(KOPT,'KPOIN',2)
      GO TO 191
C centrep0
 1934 CALL ECME(KOPT,'KPOIN',3)
      GO TO 191
C centrep1
 1935 CALL ECME(KOPT,'KPOIN',4)
      GO TO 191
C msommet
 1936 CALL ECME(KOPT,'KPOIN',5)
      GO TO 191
C Emplacements libres pour de nouveaux points
 1937 CONTINUE
 1938 CONTINUE
 1939 CONTINUE
 1940 CONTINUE
      GO TO 191
 
C Implicite
 1941 CALL ECME(KOPT,'KIMPL',1)
      CALL ECMF(KOPT,'AIMPL',1.D0)
      GO TO 191
C Explicite
 1942 CALL ECME(KOPT,'KIMPL',0)
      CALL ECMF(KOPT,'AIMPL',0.D0)
      GO TO 191
C Semi implicite
 1943 CALL ECME(KOPT,'KIMPL',2)
C?    WRITE(IOIMP,*)' EQEX KIMPL mis a 2 '
      CALL QUETYP(MTYP,0,IRET)
      IF(MTYP.EQ.'FLOTTANT')THEN
         CALL LIRREE(AIMPL,0,IRET)
         CALL ECMF(KOPT,'AIMPL',AIMPL)
      ELSE
         CALL ECMF(KOPT,'AIMPL',0.5D0)
      ENDIF
      GO TO 191
C Schema en temps implicite 2eme ordre BDF2
 1944 CALL ECME(KOPT,'ISCHT',1)
      GO TO 191
C Schema en temps implicite 4eme ordre BDF4
 1945 CALL ECME(KOPT,'ISCHT',2)
      GO TO 191
C Rajout du terme 1/2 T Div U pour stabiliser (par defaut 0)
 1946 CALL ECME(KOPT,'IDIV',1)
      GO TO 191
C Coefficient multiplicateur du decentrement (par defaut 1.)
 1947 CONTINUE
      CALL QUETYP(MTYP,0,IRET)
      IF(MTYP.EQ.'FLOTTANT')THEN
         CALL LIRREE(CMD,0,IRET)
         IF(IRET.EQ.0)THEN
         GO TO 90
         ENDIF
         CALL ECMF(KOPT,'CMD',CMD)
      ELSE
      RETURN
      ENDIF
      GO TO 191
C Format des matrices RIGIDITE IRIG = 1 MATRIK IRIG = 0 defaut
 1948 CONTINUE
      CALL ECME(KOPT,'RIGIDITE',1)
      GO TO 191
C LIMITE Limiteur divers active (Kepsilon ou autre)
 1949 CONTINUE
      CALL ECME(KOPT,'LIMITE',1)
      GO TO 191
C NODIV
 1950 CONTINUE
      CALL ECME(KOPT,'IDIV',0)
      GO TO 191
C Emplacements libres pour de nouveaux Schéma
 
C Formulation conservative
 1951 CALL ECME(KOPT,'IKOMP',1)
C Formulation non conservative
      GO TO 191
 1952 CALL ECME(KOPT,'IKOMP',0)
      GO TO 191
 1953 CALL ECME(KOPT,'IKOMP',2)
      GO TO 191
 1954 CALL ECME(KOPT,'RNG   ',2)
      GO TO 191
 1955 CALL ECME(KOPT,'ALE',1)
      GO TO 191
C Matrice masse pleine
 1956 CALL ECME(KOPT,'MTRMASS ',1)
      GO TO 191
C Matrice masse diagonale
 1957 CALL ECME(KOPT,'MTRMASS ',2)
      GO TO 191
C Matrice masse consistante (Petrov Galerkin) pour le terme source
 1958 CALL ECME(KOPT,'MTRMASS ',3)
      GO TO 191
C Matrice CONSTANTE
 1959 CALL ECME(KOPT,'KMACO',1)
      CALL LIRCHA(CHAI,1,LCHAR)
          IF(LCHAR.EQ.0)THEN
C On ne trouve pas d'objet de type %m1:8
           MOTERR( 1: 8) = 'MOT'
           CALL ERREUR(38)
           RETURN
          ENDIF
      NMACO=CHAI(1:LCHAR)
      CALL ECMM(KOPT,'NMACO',NMACO)
      GO TO 191
C Emplacement libre pour de nouvelles idées (il faudra etre concis)
 1960 CONTINUE
      GO TO 191
 
C Indices  IDEUL
C EULER
 1961 CALL ECME(KOPT,'IDEUL',1)
      GO TO 191
C EULERMS
 1962 CALL ECME(KOPT,'IDEUL',2)
      GO TO 191
C EULERMST
 1963 CALL ECME(KOPT,'IDEUL',3)
      GO TO 191
C Emplacements libres pour de nouveaux Schéma
 1964 CONTINUE
 1965 CONTINUE
      GO TO 191
C Indices  KPOIND
 1966 CONTINUE
      CALL LIRCHA(NOM,0,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
      CALL OPTLI(IP,LOPTI(31),NOM,6)
      IF(IP.EQ.0)THEN
         GO TO 90
      ELSE
         CALL ECME(KOPT,'KPOIND',IP-1)
      ENDIF
      GO TO 191
 
C Emplacements libres
 1967 CONTINUE
      GO TO 191
 
 1968 CONTINUE
C Indice  STABP
      CALL LIRREE(STAB,0,IRET)
      IF(IRET.EQ.0)THEN
         GO TO 90
      ENDIF
         CALL ECMF(KOPT,'STABP',STAB)
      GO TO 191
 
C MUCONS mu constant par élément
 1969 CALL ECME(KOPT,'KMU',0)
      GO TO 191
 
C FTAU mu variable par élément (formulation en grad mu)
 1970 CALL ECME(KOPT,'KMU',1)
      GO TO 191
 
C MUVARI Formulation en Tau
 1971 CALL ECME(KOPT,'KMU',2)
      GO TO 191
 
C Emplacements libres
 1972 CONTINUE
 1973 CONTINUE
 1974 CONTINUE
 1975 CONTINUE
      GO TO 191
 
 
C OPTI
 
 20   CONTINUE
      CALL LIRCHA(NOM,0,IRET)
      IF(IRET.EQ.0)GO TO 90
      CALL ECMM(MTABLE,'NOMVI',NOM)
      GO TO 1
 
C '    '
 
 21   CONTINUE
C     WRITE(IOIMP,*)' nbik=',nbik
      NBIK=NBIK+1
      CALL LIRCHA(NOM,0,IRET)
C     WRITE(IOIMP,*)' NOM=',nom,iret
      IF(IRET.EQ.0)GO TO 90
      TINCD(NBIK)=NOM
      CALL LITABS(LTAB,KTAB,1,1,IRET)
C     WRITE(IOIMP,*)' iret=',iret
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' On attend un objet TABLE DOMAINE'
         RETURN
      ENDIF
      KINCD(NBIK)=KTAB(1)
 
      GO TO 1
 
C TPSI
 
 22   CONTINUE
      CALL LIRREE(XVAL,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE TPSI (Instant initial) :'
         WRITE(IOIMP,*)' On attend un flottant '
         RETURN
      ENDIF
      CALL ECMF(MTABLE,'TPSI',XVAL)
      GO TO 1
 
C TFINAL
 
 23   CONTINUE
      CALL LIRREE(XVAL,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE TFINAL (Temps final) :'
         WRITE(IOIMP,*)' On attend un flottant '
         RETURN
      ENDIF
      CALL ECMF(MTABLE,'TFINAL',XVAL)
      GO TO 1
 
C FIDT
 
 24   CONTINUE
      CALL LIRENT(IENT,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE FIDT (Frequence impression temps) :'
         WRITE(IOIMP,*)' On attend un entier '
         RETURN
      ENDIF
      CALL ECME(MTABLE,'FIDT',IENT)
      GO TO 1
 
C NISTO
 
 25   CONTINUE
      CALL LIRENT(IENT,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE NISTO (Frequence saisie historique) :'
         WRITE(IOIMP,*)' On attend un entier '
         RETURN
      ENDIF
      CALL ECME(MTABLE,'NISTO',IENT)
      GO TO 1
 
C NITER
 
 26   CONTINUE
      CALL LIRENT(IENT,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE NITER (Nombre iterations internes) :'
         WRITE(IOIMP,*)' On attend un entier '
         RETURN
      ENDIF
      CALL ECME(MTABLE,'NITER',IENT)
      GO TO 1
 
C OMEGA
 
 27   CONTINUE
      CALL LIRREE(XVAL,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE OMEGA (Facteur de relaxation) :'
         WRITE(IOIMP,*)' On attend un flottant '
         RETURN
      ENDIF
      CALL ECMF(MTABLE,'OMEGA',XVAL)
      GO TO 1
 
C EPS
 
 28   CONTINUE
      CALL LIRREE(XVAL,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE EPS (Tolerance sur le residu) :'
         WRITE(IOIMP,*)' On attend un flottant '
         RETURN
      ENDIF
      CALL ECMF(MTABLE,'EPS',XVAL)
      GO TO 1
 
C IMPR
 
 29   CONTINUE
      CALL LIRENT(IENT,1,IRET)
      IF(IRET.EQ.0)THEN
         WRITE(IOIMP,*)' MOT CLE IMPR (Niveau d impression) :'
         WRITE(IOIMP,*)' On attend un entier '
         RETURN
      ENDIF
      CALL ECME(MTABLE,'IMPR',IENT)
      GO TO 1
 
 90   CONTINUE
 
      TYPE=' '
      CALL ACMO(MTABLE,'DOMAINE',TYPE,MTABD)
      IF(TYPE.NE.'TABLE')GO TO 900
      TYPE=' '
      CALL ACMO(MTABLE,'LISTINCO',TYPE,MLMOT2)
 
      IF(MLMOT2.NE.0)THEN
 
         SEGACT MLMOT2
         NINCT=MLMOT2.MOTS(/2)
 
         TYPE=' '
         CALL ACMO(MTABLE,'DOMINC',TYPE,MTABI)
         IF(TYPE.EQ.'TABLE')THEN
 
            DO 93 I=1,NINCT
               NOMI=MLMOT2.MOTS(I)
               TYPE=' '
               CALL ACMO(MTABI,NOMI,TYPE,IPT)
               IF(TYPE.NE.'TABLE')CALL ECMO(MTABI,NOMI,'TABLE',MTABD)
 93         CONTINUE
            DO 94 I=1,NBIK
               NOMI=TINCD(NBIK)
               CALL ECMO(MTABI,NOMI,'TABLE',KINCD(NBIK))
 94         CONTINUE
 
         ELSE
 
            CALL CRTABL(MTABI)
            CALL ECMO(MTABLE,'DOMINC','TABLE',MTABI)
            DO 91 I=1,NINCT
               NOMI=MLMOT2.MOTS(I)
               CALL ECMO(MTABI,NOMI,'TABLE',MTABD)
 91         CONTINUE
            DO 92 I=1,NBIK
               NOMI=TINCD(NBIK)
               CALL ECMO(MTABI,NOMI,'TABLE',KINCD(NBIK))
 92         CONTINUE
 
         ENDIF
            SEGDES MLMOT2
      ENDIF
 
 900  CONTINUE
      CALL ECROBJ('TABLE',MTABLE)
C     write(6,*)' RETOUR EQEX '
      RETURN
      END