Sources Esope | Cast3M

pre121
C PRE121    SOURCE    OF166741  24/10/03    21:15:29     12022          
      SUBROUTINE PRE121(LOGTEM,
     &                  ICEN,IFACE,IFACEL,INORM,
     &                  IROC, IGRROC, IALROC,
     &                  IVITC, IGRVC, IALVC,
     &                  IPC ,IGRPC, IALPC,
     &                  IGAMC,
     &                  DELTAT,
     &                  IROF,IVITF,IPF,IGAMF,
     &                  LOGAN,LOGNEG,LOGBOR,MESERR,VALER,VAL1,VAL2)
C************************************************************************
C
C PROJET            :  CASTEM 2000
C
C NOM               :  PRE121
C
C DESCRIPTION       :  Voir PRE12
C
C LANGAGE           :  FORTRAN 77 + ESOPE 2000 (avec extensions CISI)
C
C AUTEUR            :  A. BECCANTINI, DRN/DMT/SEMT/TTMF
C
C************************************************************************
C
C
C APPELES (Outils)  :  KRIPAD, LICHT
C
C APPELES (Calcul)  :  AUCUN
C
C
C************************************************************************
C
C ENTREES
C
C     LOGTEM  : LOGICAL; si .TRUE. 2em ordre en temps
C                        sinon 1er ordre en temps
C
C     1) Pointeurs de MELEMEs et de CHPOINTs de la table DOMAINE
C
C     ICEN    : MELEME de 'POI1' SPG des CENTRES
C
C     IFACE   : MELEME de 'POI1' SPG des FACES
C
C     IFACEL  : MELEME de 'SEG3' avec
C               CENTRE d'Elt "gauche"
C               CENTRE de Face
C               CENTRE d'Elt "droite"
C
C               N.B. = IFACE.NUM(i,1) = IFACEL.NUM(i,2)
C
C     INORM   : CHPOINT des cosinus directeurs de normales aux faces
C
C     2) Pointeurs des  CHPOINTs
C
C     IROC    : CHPOINT "CENTRE" contenant la masse volumique RHO
C
C     IGRROC  : CHPOINT "CENTRE" contenant le gradient de la
C               masse volumique RHO (2 composantes)
C
C     IALROC  : CHPOINT "CENTRE" contenant le limiteur du gradient
C               de la masse volumique
C
C     IVITC   : CHPOINT "CENTRE" contenant la vitesse UX, UY ;
C
C     IGRVC   : CHPOINT "CENTRE" contenant le gradient de la
C               vitesse (4 composantes)
C
C     IALVC  :  CHPOINT "CENTRE" contenant le limiteur du gradient
C               de la vitesse (2 composantes)
C
C     IPC     : CHPOINT "CENTRE" contenat la pression P;
C
C     IGRPC   : CHPOINT "CENTRE" contenant le gradient de la
C               pression (2 composantes)
C
C     IALPC   : CHPOINT "CENTRE" contenant le limiteur du gradient
C               de la pression
C
C     IGAMC   : CHPOINT "CENTRE" contenat le "Gamma" du gaz
C
C     3)
C
C     DELTAT  : REAL*8, encrement en temps pour calculer la prediction
C
C
C SORTIES
C
C
C     IROF    : MCHAML  defini sur le MELEME de pointeur IFACEL,
C               contenant la masse volumique RHO
C
C     IVITF   : MCHAML  defini sur le MELEME de pointeur IFACEL,
C               contenant la vitesse UN, UT dans le repaire local
C               (n,t) et defini sur le MELEME de pointeur IFACE,
C               contenant les cosinus directeurs du repere local
C
C     IPF     : MCHAML  defini sur le MELEME de pointeur IFACEL,
C               contenant la pression P
C
C     IGAMF   : MCHAML  defini sur le MELEME de pointeur IFACEL,
C               contenant le "gamma" du gaz
C
C     LOGAN   : anomalie detectee
C
C     LOGNEG  : (LOGICAL): si .TRUE. une pression ou une densité
C               negative a été detectée -> en interactif le
C               programme s'arrete  en GIBIANE
C               (erreur stocké en MESERR et VALER)
C
C     LOGBOR  : (LOGICAL): si .TRUE. un gamma a ete detecte
C               dehor 1 et 3 (sa valeur stockée en MESERR et VALER;
C               en VAL1 et en VAL2 on stocke 1.0 et 3.0)
C
C     MESERR
C     VALER
C     VAL1,
C     VAL2    : pour les messages d'erreur
C
C************************************************************************
C
C HISTORIQUE (Anomalies et modifications éventuelles)
C
C HISTORIQUE :  Créée le 11.6.98.
C
C************************************************************************
C
C
C ATTENTION: Cet programme marche si le MAILLAGE est convex;
C            si non il faut changer l'argoritme de calcul de
C            l'orientation des normales aux faces.
C
C            La positivité n'est pas controlle parce que c'est déjà fait
C            dans l'operateur PRIM
C
C
C************************************************************************
C
C**** Les variables
C
      IMPLICIT INTEGER(I-N)
 
-INC PPARAM
-INC CCOPTIO
 
      INTEGER ICEN, IFACE, IFACEL, INORM,IROC, IGRROC, IALROC
     &       , IVITC, IGRVC, IALVC
     &       , IPC ,IGRPC, IALPC
     &       , IGAMC
     &       , IROF, IVITF, IPF, IGAMF
     &       , IGEOM, NFAC, NCEN
     &       , N1PTEL, N1EL, N2PTEL, N2EL, N2, N1, N3, L1, NLCE
     &       , NLCF, NGCF, NGCEG, NLCEG, NGCED, NLCED, NGCF1
     &       , IDIMP1, INDCEL
      REAL*8 VALER, VAL1, VAL2, XG, YG, XC, YC, XD, YD, DELTAT
     &        ,DXG, DYG, DXD, DYD
     &       , CNX, CNY, CTX, CTY, ORIENT
     &       , ROG, PG, GAMG, UXG, UYG, UNG, UTG
     &       , ROD, PD, GAMD, UXD, UYD, UND, UTD
     &       , VALCEL, DCEL, ALCEL
     &       , DROX, DROY, DUXX, DUXY, DUYX, DUYY, DPX, DPY
     &       , DRO, DUX, DUY, DP, ALPHA(4)
      CHARACTER*(40) MESERR
      CHARACTER*(8) TYPE
      LOGICAL LOGAN,LOGNEG, LOGBOR, LOGTEM, LOGI1, LOGALP
C
C**** LOGALP = .TRUE. -> Prediction avec limiteur
C
      PARAMETER(LOGALP = .TRUE.)
C
C**** Les Includes
C
-INC SMCOORD
 
-INC SMCHPOI
      POINTEUR MPROC.MPOVAL,  MPGRR.MPOVAL,
     &         MPVITC.MPOVAL, MPGRV.MPOVAL,
     &         MPPC.MPOVAL,   MPGRP.MPOVAL,
     &         MPGAMC.MPOVAL, MPNORM.MPOVAL,
     &         MPROP.MPOVAL,  MPPP.MPOVAL, MPVITP.MPOVAL
-INC SMCHAML
      POINTEUR MELVNX.MELVAL, MELVNY.MELVAL,
     &         MELT1X.MELVAL, MELT1Y.MELVAL
      POINTEUR MELVUN.MELVAL, MELVUT.MELVAL
      POINTEUR MELRO.MELVAL, MELP.MELVAL,
     &         MELGAM.MELVAL
      POINTEUR MELALR.MPOVAL,
     &         MELALV.MPOVAL,
     &         MELALP.MPOVAL
-INC SMLENTI
-INC SMELEME
C
C
C**** Initialisation des parametres d'erreur déjà faite, i.e.
C
C     LOGNEG  = .FALSE.
C     LOGBOR  = .FALSE.
C     MESERR = '                                        '
C     MOTERR(1:40) = MESERR(1:40)
C     VALER = 0.0D0
C     VAL1 = 0.0D0
C     VAL2 = 0.0D0
C
C
C**** KRIPAD pour la correspondance global/local de centre
C
      CALL KRIPAD(ICEN,MLENT1)
C
C**** MLENTI1 a MCORD.XCOORD(/1)/(IDIM+1) elements
C
C     Si i est le numero global d'un noeud de ICEN,
C     MLENT1.LECT(i) contient sa position, i.e.
C
C     I              = numero global du noeud centre
C     MLENT1.LECT(i) = numero local du noeud centre
C
C     MLENT1 déjà activé, i.e.
C
C     SEGACT MLENT1
C
C**** Activation de CHPOINTs
C
C     densité  + grad + limiteur
C     vitesse  + grad + limiteur
C     pression + grad + limiteur
C     gamma
C     cosinus directeurs des normales aux surface
C
      CALL LICHT(IROC  , MPROC  , TYPE, IGEOM)
      CALL LICHT(IGRROC, MPGRR  , TYPE, IGEOM)
      CALL LICHT(IVITC,  MPVITC , TYPE, IGEOM)
      CALL LICHT(IGRVC,  MPGRV  , TYPE, IGEOM)
      CALL LICHT(IPC  ,  MPPC   , TYPE, IGEOM)
      CALL LICHT(IGRPC,  MPGRP  , TYPE, IGEOM)
      CALL LICHT(IGAMC,  MPGAMC , TYPE, IGEOM)
      CALL LICHT(INORM,  MPNORM , TYPE, IGEOM)
C
C**** Les MPOVALs 'Prediction'
C
      IF(LOGTEM)THEN
         SEGINI, MPROP = MPROC
         SEGINI, MPPP = MPPC
         SEGINI, MPVITP = MPVITC
      ELSE
         MPROP = MPROC
         MPPP = MPPC
         MPVITP = MPVITC
      ENDIF
C
C**** Les Limiteurs
C
      CALL LICHT(IALROC, MELALR , TYPE, IGEOM)
      CALL LICHT(IALVC,  MELALV , TYPE, IGEOM)
      CALL LICHT(IALPC,  MELALP , TYPE, IGEOM)
C
C
C**** MPOVAL sont déjà activés i.e.:
C
C     SEGACT MPROC
C     SEGACT MPGRR
C     SEGACT MPIALR
C     SEGACT MPVITC
C     SEGACT MPGRV
C     SEGACT MPIALV
C     SEGACT MPPC
C     SEGACT MPGRP
C     SEGACT MPIALP
C     SEGACT MPGAMC
C     SEGACT MPNORM
C
C**** Le MELEME FACEL
C
      IPT1 = IFACEL
      IPT2 = IFACE
      SEGACT IPT1
      SEGACT IPT2
      NFAC = IPT1.NUM(/2)
C
C**** Creation de MCHAMLs contenant les etats gauche et droite,
C
C     i.e.:
C
C     vitesse + cosinus directors du repere local
C     densité
C     pression
C     gamma
C
C**** Cosinus directors du repere local et vitesse
C
C     Les cosinus directeurs
C
      N1 = 2
      N3 = 6
      L1 = 28
      SEGINI MCHEL1
      IVITF = MCHEL1
      MCHEL1.TITCHE = 'U  '
      MCHEL1.IMACHE(1) = IFACE
      MCHEL1.IMACHE(2) = IFACEL
      MCHEL1.CONCHE(1) = ' (n,t) in (x,y) '
      MCHEL1.CONCHE(2) = ' U in (n,t)     '
C
C**** Valeurs des cosinus definies par respect au repere global, i.e.
C
      MCHEL1.INFCHE(1,1) = 2
      MCHEL1.INFCHE(1,3) = NIFOUR
      MCHEL1.INFCHE(1,4) = 0
      MCHEL1.INFCHE(1,5) = 0
      MCHEL1.INFCHE(1,6) = 1
      MCHEL1.IFOCHE = IFOUR
C
C**** Valeurs de vitesse definies par respect au repere local, i.e.
C
      MCHEL1.INFCHE(2,1) = 1
      MCHEL1.INFCHE(2,3) = NIFOUR
      MCHEL1.INFCHE(2,4) = 0
      MCHEL1.INFCHE(2,5) = 0
      MCHEL1.INFCHE(2,6) = 1
C
C**** Le cosinus directeurs
C
      N1PTEL = 1
      N1EL = NFAC
      N2PTEL = 0
      N2EL =   0
C
C**** MCHAML a N2 composantes:
C
C     cosinus directeurs du repere local (n,t1)
C
C     IDIM = 2 -> 4 composantes
C
      N2 = 4
      SEGINI MCHAM1
      MCHEL1.ICHAML(1) = MCHAM1
      MCHAM1.NOMCHE(1) = 'NX      '
      MCHAM1.NOMCHE(2) = 'NY      '
      MCHAM1.NOMCHE(3) = 'TX      '
      MCHAM1.NOMCHE(4) = 'TY      '
      MCHAM1.TYPCHE(1) = 'REAL*8          '
      MCHAM1.TYPCHE(2) = 'REAL*8          '
      MCHAM1.TYPCHE(3) = 'REAL*8          '
      MCHAM1.TYPCHE(4) = 'REAL*8          '
      SEGINI MELVNX
      SEGINI MELVNY
      SEGINI MELT1X
      SEGINI MELT1Y
      MCHAM1.IELVAL(1) = MELVNX
      MCHAM1.IELVAL(2) = MELVNY
      MCHAM1.IELVAL(3) = MELT1X
      MCHAM1.IELVAL(4) = MELT1Y
C
C**** Vitesse
C
      N1EL = NFAC
      N1PTEL = 3
      N2EL = 0
      N2PTEL = 0
C
C**** MCHAML a N2 composantes:
C
C     IDIM = 2 -> 2 composantes
C
      N2 = 2
      SEGINI MCHAM1
      MCHEL1.ICHAML(2) = MCHAM1
      MCHAM1.NOMCHE(1) = 'UN      '
      MCHAM1.NOMCHE(2) = 'UT      '
      MCHAM1.TYPCHE(1) = 'REAL*8          '
      MCHAM1.TYPCHE(2) = 'REAL*8          '
      SEGINI MELVUN
      SEGINI MELVUT
      MCHAM1.IELVAL(1) = MELVUN
      MCHAM1.IELVAL(2) = MELVUT
C
C**** Densite
C
      N1 = 1
      N3 = 6
      L1 = 15
      SEGINI MCHEL2
      IROF = MCHEL2
      MCHEL2.IMACHE(1) = IFACEL
      MCHEL2.TITCHE = 'RO             '
      MCHEL2.CONCHE(1) = '                '
C
C**** Valeurs independente du repére, i.e.
C
      MCHEL2.INFCHE(1,1) = 0
      MCHEL2.INFCHE(1,3) = NIFOUR
      MCHEL2.INFCHE(1,4) = 0
      MCHEL2.INFCHE(1,5) = 0
      MCHEL2.INFCHE(1,6) = 1
      MCHEL2.IFOCHE = IFOUR
      N2 = 1
      SEGINI MCHAM1
      MCHEL2.ICHAML(1) = MCHAM1
      MCHAM1.NOMCHE(1) = 'SCAL    '
      MCHAM1.TYPCHE(1) = 'REAL*8          '
      SEGINI MELRO
      MCHAM1.IELVAL(1) = MELRO
C
C**** Pression
C
      MCHEL1 = IROF
      SEGINI, MCHEL2 = MCHEL1
      IPF = MCHEL2
      MCHEL2.TITCHE = 'P            '
C
C**** MCHAM1 = MCHAML de la densite
C
      SEGINI, MCHAM2 = MCHAM1
      MCHEL2.ICHAML(1) = MCHAM2
      SEGINI MELP
      MCHAM2.IELVAL(1) = MELP
C
C**** Gamma
C
      MCHEL1 = IROF
      SEGINI, MCHEL2 = MCHEL1
      IGAMF = MCHEL2
      MCHEL2.TITCHE = 'GAMMA        '
C
C**** MCHAM1 = MCHAML de la densite
C
      SEGINI, MCHAM2 = MCHAM1
      MCHEL2.ICHAML(1) = MCHAM2
      SEGINI MELGAM
      MCHAM2.IELVAL(1) = MELGAM
C
C**** Recapitulatif: les MELVALs et les MPOVALs actives
C
C     MELVNX, MELVNY
C     MELT1X, MELT1Y
C
C     MELVUN, MELVUT -> vitesse
C
C     MELRO  -> densite
C
C     MELP   -> pression
C
C     MELGAM -> gamma
C
C     MPROC, MPROP  -> densite
C
C     MPVITC, MPVITP -> vitesse
C
C     MPPC, MPPP   -> pression
C
C     MPGAMC -> gamma
C
C     MPNORM -> normales aux faces
C
C
C***********************************************************************
C********* PREDICTION **************************************************
C***********************************************************************
C
      IF(LOGTEM)THEN
C
         IPT3 = ICEN
         SEGACT IPT3
         NCEN = IPT3.NUM(/2)
         DO NLCE = 1, NCEN
            IF(LOGALP)THEN
C
C************* Gradients limités
C
               ALPHA(1) = MELALR.VPOCHA(NLCE,1)
               ALPHA(2) = MELALV.VPOCHA(NLCE,1)
               ALPHA(3) = MELALV.VPOCHA(NLCE,2)
               ALPHA(4) = MELALP.VPOCHA(NLCE,1)
            ELSE
C
C************* Gradients non limités
C
               ALPHA(1) = 1.0D0
               ALPHA(2) = 1.0D0
               ALPHA(3) = 1.0D0
               ALPHA(4) = 1.0D0
            ENDIF
            ROG = MPROP.VPOCHA(NLCE,1)
            UXG = MPVITP.VPOCHA(NLCE,1)
            UYG = MPVITP.VPOCHA(NLCE,2)
            PG = MPPP.VPOCHA(NLCE,1)
            DROX = MPGRR.VPOCHA(NLCE,1)*ALPHA(1)
            DROY = MPGRR.VPOCHA(NLCE,2)*ALPHA(1)
            DUXX = MPGRV.VPOCHA(NLCE,1)*ALPHA(2)
            DUXY = MPGRV.VPOCHA(NLCE,2)*ALPHA(2)
            DUYX = MPGRV.VPOCHA(NLCE,3)*ALPHA(3)
            DUYY = MPGRV.VPOCHA(NLCE,4)*ALPHA(3)
            DPX  = MPGRP.VPOCHA(NLCE,1)*ALPHA(4)
            DPY  = MPGRP.VPOCHA(NLCE,2)*ALPHA(4)
            GAMG = MPGAMC.VPOCHA(NLCE,1)
            DRO  = UXG * DROX + ROG * ( DUXX + DUYY )
     &           + UYG * DROY
            DUX  = UXG * DUXX + DPX / ROG + UYG * DUXY
            DUY  = UXG * DUYX + UYG * DUYY + DPY / ROG
            DP   = GAMG * PG * (DUXX + DUYY)
     &           + UXG * DPX + UYG * DPY
C
C************* Prediction avec/sans gradients limités
C
            MPROP.VPOCHA(NLCE,1) = ROG - DELTAT * DRO
            MPVITP.VPOCHA(NLCE,1) = UXG - DELTAT * DUX
            MPVITP.VPOCHA(NLCE,2) = UYG - DELTAT * DUY
            MPPP.VPOCHA(NLCE,1) = PG - DELTAT * DP
         ENDDO
C
      ENDIF
C
C
C***********************************************************************
C********* CORRECTION **************************************************
C***********************************************************************
C
C**** Boucle sur le faces
C
      IDIMP1 = IDIM + 1
      DO NLCF = 1, NFAC
C
C******* NLCF  = numero local du centre de face
C        NGCF  = numero global du centre de face
C        NGCEG = numero global du centre ELT "gauche"
C        NLCEG = numero local du centre ELT "gauche"
C        NGCED = numero global du centre ELT "droite"
C        NLCED = numero local du centre ELT "droite"
C
         NGCEG = IPT1.NUM(1,NLCF)
         NGCF  = IPT1.NUM(2,NLCF)
         NGCED = IPT1.NUM(3,NLCF)
         NLCEG = MLENT1.LECT(NGCEG)
         NLCED = MLENT1.LECT(NGCED)
C
C******* TEST: IPT2.NUM(1,NLCF) = IPT1.NUM(2,NLCF)
C
         NGCF1 = IPT2.NUM(1,NLCF)
         IF( NGCF1 .NE. NGCF) THEN
            LOGAN = .TRUE.
            MESERR(1:40) = 'PRET, subroutine pre121.eso             '
            GOTO 9999
         ENDIF
C
C******* Cosinus directeurs des NORMALES aux faces
C
C        On impose que les normales sont direct "Gauche" -> "Centre"
C
         INDCEL = (NGCEG-1)*IDIMP1
         XG = XCOOR(INDCEL+1)
         YG = XCOOR(INDCEL+2)
         INDCEL = (NGCF-1)*IDIMP1
         XC = XCOOR(INDCEL + 1)
         YC = XCOOR(INDCEL + 2)
         INDCEL = (NGCED-1)*IDIMP1
         XD = XCOOR(INDCEL+1)
         YD = XCOOR(INDCEL+2)
         DXG = XC - XG
         DYG = YC - YG
         DXD = XC - XD
         DYD = YC - YD
C
C******* On calcule le sign du pruduit scalare
C        (Normales de Castem) * (vecteur "gauche" -> "centre")
C
         CNX = MPNORM.VPOCHA(NLCF,1)
         CNY = MPNORM.VPOCHA(NLCF,2)
         ORIENT = CNX * DXG + CNY * DYG
         ORIENT = SIGN(1.0D0,ORIENT)
         IF(ORIENT .NE. 1.0D0)THEN
            LOGAN = .TRUE.
            MESERR(1:30)=
     &           'PRET , subroutine pre121.eso. '
            GOTO 9999
         ENDIF
         CNX = CNX * ORIENT
         CNY = CNY * ORIENT
C
C********** Cosinus directeurs de tangent 2D
C
         CTX = -1.0D0 * CNY
         CTY = CNX
C
C
C******* Les autres MELVALs
C
C
C******* N.B.: On suppose qu'on a déjà controlle RO, P > 0
C              GAMMA \in (1,3)
C              Si non il faut le faire, en utilisant LOGBOR,
C              LOGNEG, VALER, VAL1, VAL2
C
C
C
C******* NGCEG = NGCED -> Mur
C
         IF( NGCEG .EQ. NGCED)THEN
C
C********** Sur le mur on fait de reconstruction sur l'etat gauche
C
C
C********** Etat gauche
C
            VALCEL = MPROP.VPOCHA(NLCEG, 1)
            ALCEL = MELALR.VPOCHA(NLCEG, 1)
            DCEL = MPGRR.VPOCHA(NLCEG, 1)*DXG +
     &             MPGRR.VPOCHA(NLCEG, 2)*DYG
            ROG = VALCEL + ALCEL * DCEL
C
            VALCEL = MPPP.VPOCHA(NLCEG, 1)
            ALCEL =  MELALP.VPOCHA(NLCEG, 1)
            DCEL = MPGRP.VPOCHA(NLCEG, 1)*DXG +
     &             MPGRP.VPOCHA(NLCEG, 2)*DYG
            PG = VALCEL + ALCEL * DCEL
C
            GAMG = MPGAMC.VPOCHA(NLCEG, 1)
C
            VALCEL = MPVITP.VPOCHA(NLCEG, 1)
            ALCEL = MELALV.VPOCHA(NLCEG, 1)
            DCEL = MPGRV.VPOCHA(NLCEG, 1)*DXG +
     &             MPGRV.VPOCHA(NLCEG, 2)*DYG
            UXG = VALCEL + ALCEL * DCEL
C
            VALCEL = MPVITP.VPOCHA(NLCEG, 2)
            ALCEL = MELALV.VPOCHA(NLCEG, 2)
            DCEL = MPGRV.VPOCHA(NLCEG, 3)*DXG +
     &             MPGRV.VPOCHA(NLCEG, 4)*DYG
            UYG = VALCEL + ALCEL * DCEL
C
            UNG = UXG * CNX + UYG * CNY
            UTG = UXG * CTX + UYG * CTY
C
C********** Si l'on fait pas de prediction, ce n'est pas necessaire de
C           controller la positivite' de la pression et de la densité; elle
C           est déjà garantie par la proprieté LED de limiteur; mais, vu
C           que le limiteur n'est pas calculé ici, mais dans un autre
C           operateur, on le fait
C
            LOGI1 = (PG .LT. 0.0D0) .OR. (ROG .LT. 0.0D0)
C
            IF(LOGI1)THEN
C
C************* Premier ordre en espace local
C
               ROG = MPROC.VPOCHA(NLCEG,1)
               PG  = MPPC.VPOCHA(NLCEG,1)
               UNG = MPVITC.VPOCHA(NLCEG,1)*CNX +
     &              MPVITC.VPOCHA(NLCEG,2)*CNY
               UTG = MPVITC.VPOCHA(NLCEG,1)*CTX +
     &              MPVITC.VPOCHA(NLCEG,2)*CTY
            ENDIF
C
C********** Son etat droite
C
            ROD = ROG
            PD = PG
            GAMD = GAMG
            UND = -1.0D0 * UNG
            UTD = UTG
C
C************* Fin cas mur
C
         ELSE
C
C************* Etat gauche
C
            VALCEL = MPROP.VPOCHA(NLCEG, 1)
            ALCEL = MELALR.VPOCHA(NLCEG, 1)
            DCEL = MPGRR.VPOCHA(NLCEG, 1)*DXG +
     &             MPGRR.VPOCHA(NLCEG, 2)*DYG
            ROG = VALCEL + ALCEL * DCEL
C
            VALCEL = MPPP.VPOCHA(NLCEG, 1)
            ALCEL =  MELALP.VPOCHA(NLCEG, 1)
            DCEL = MPGRP.VPOCHA(NLCEG, 1)*DXG +
     &             MPGRP.VPOCHA(NLCEG, 2)*DYG
            PG = VALCEL + ALCEL * DCEL
C
            GAMG = MPGAMC.VPOCHA(NLCEG, 1)
C
            VALCEL = MPVITP.VPOCHA(NLCEG, 1)
            ALCEL = MELALV.VPOCHA(NLCEG, 1)
            DCEL = MPGRV.VPOCHA(NLCEG, 1)*DXG +
     &             MPGRV.VPOCHA(NLCEG, 2)*DYG
            UXG = VALCEL + ALCEL * DCEL
C
            VALCEL = MPVITP.VPOCHA(NLCEG, 2)
            ALCEL = MELALV.VPOCHA(NLCEG, 2)
            DCEL = MPGRV.VPOCHA(NLCEG, 3)*DXG +
     &             MPGRV.VPOCHA(NLCEG, 4)*DYG
            UYG = VALCEL + ALCEL * DCEL
C
            UNG = UXG * CNX + UYG * CNY
            UTG = UXG * CTX + UYG * CTY
C
C********** Positivite
C
            LOGI1 = (PG .LT. 0.0D0) .OR. (ROG .LT. 0.0D0)
C
            IF(LOGI1)THEN
C
C************* Premier ordre en espace local
C
               ROG = MPROC.VPOCHA(NLCEG,1)
               PG  = MPPC.VPOCHA(NLCEG,1)
               UNG = MPVITC.VPOCHA(NLCEG,1)*CNX +
     &              MPVITC.VPOCHA(NLCEG,2)*CNY
               UTG = MPVITC.VPOCHA(NLCEG,1)*CTX +
     &              MPVITC.VPOCHA(NLCEG,2)*CTY
            ENDIF
C
C********** Etat droite
C
            VALCEL = MPROP.VPOCHA(NLCED, 1)
            ALCEL = MELALR.VPOCHA(NLCED, 1)
            DCEL = MPGRR.VPOCHA(NLCED, 1)*DXD +
     &             MPGRR.VPOCHA(NLCED, 2)*DYD
            ROD = VALCEL + ALCEL * DCEL
C
            VALCEL = MPPP.VPOCHA(NLCED, 1)
            ALCEL =  MELALP.VPOCHA(NLCED, 1)
            DCEL = MPGRP.VPOCHA(NLCED, 1)*DXD +
     &             MPGRP.VPOCHA(NLCED, 2)*DYD
            PD = VALCEL + ALCEL * DCEL
C
            GAMD = MPGAMC.VPOCHA(NLCED, 1)
C
            VALCEL = MPVITP.VPOCHA(NLCED, 1)
            ALCEL = MELALV.VPOCHA(NLCED, 1)
            DCEL = MPGRV.VPOCHA(NLCED, 1)*DXD +
     &             MPGRV.VPOCHA(NLCED, 2)*DYD
            UXD = VALCEL + ALCEL * DCEL
C
            VALCEL = MPVITP.VPOCHA(NLCED, 2)
            ALCEL = MELALV.VPOCHA(NLCED, 2)
            DCEL = MPGRV.VPOCHA(NLCED, 3)*DXD +
     &             MPGRV.VPOCHA(NLCED, 4)*DYD
            UYD = VALCEL + ALCEL * DCEL
C
            UND = UXD * CNX + UYD * CNY
            UTD = UXD * CTX + UYD * CTY
C
C********** Positivite
C
            LOGI1 = (PD .LT. 0.0D0) .OR. (ROD .LT. 0.0D0)
C
            IF(LOGI1)THEN
C
C************* Premier ordre en espace local
C
               ROD = MPROC.VPOCHA(NLCED,1)
               PD  = MPPC.VPOCHA(NLCED,1)
               UND = MPVITC.VPOCHA(NLCED,1)*CNX +
     &              MPVITC.VPOCHA(NLCED,2)*CNY
               UTD = MPVITC.VPOCHA(NLCED,1)*CTX +
     &              MPVITC.VPOCHA(NLCED,2)*CTY
            ENDIF
         ENDIF
C
C******** Les MELVALs
C
         MELRO.VELCHE(1,NLCF)  = ROG
         MELRO.VELCHE(3,NLCF)  = ROD
         MELP.VELCHE(1,NLCF)   = PG
         MELP.VELCHE(3,NLCF)   = PD
         MELGAM.VELCHE(1,NLCF) = GAMG
         MELGAM.VELCHE(3,NLCF) = GAMD
         MELVUN.VELCHE(1,NLCF) = UNG
         MELVUN.VELCHE(3,NLCF) = UND
         MELVUT.VELCHE(1,NLCF) = UTG
         MELVUT.VELCHE(3,NLCF) = UTD
         MELVNX.VELCHE(1,NLCF) = CNX
         MELVNY.VELCHE(1,NLCF) = CNY
         MELT1X.VELCHE(1,NLCF) = CTX
         MELT1Y.VELCHE(1,NLCF) = CTY
C
      ENDDO
C
C**** Le MPOVALs 'Prediction' sont detruits (si existentes)
C
      IF(LOGTEM)THEN
         SEGSUP MPROP
         SEGSUP MPVITP
         SEGSUP MPPP
      ENDIF
C**** Destruction du MLENTI correspondance local/global
C
      SEGSUP MLENT1
C
 9999 CONTINUE
 
      END