Télécharger pre321.eso

Retour à la liste

pre321
C PRE321    SOURCE    CB215821  20/11/25    13:36:23     10792          
      SUBROUTINE PRE321(LOGTEM,LGAMC,LYC,LSCAC,
     &                  ICEN,IFACE,IFACEL,INORM,
     &                  IROC, IGRROC, IALROC,
     &                  IVITC, IGRVC, IALVC,
     &                  IPC ,IGRPC, IALPC,
     &                  IYC ,IGRYC, IALYC,
     &                  ISCAC ,IGRSC, IALSC,
     &                  IGAMC,
     &                  DELTAT,
     &                  IROF,IVITF,IPF,IYF,ISCAF,
     &                  LOGAN,LOGNEG,LOGBOR,MESERR,VALER,VAL1,VAL2)
C************************************************************************
C
C PROJET            :  CASTEM 2000
C
C NOM               :  PRE321
C
C DESCRIPTION       :  Voir PRE32
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     LGAMC,LYC,LSCAC : LOGICAL : si .TRUE. IGAMC, IYC, ISCAC sont
C                       pointeurs de CHPOINTS
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) Autres pointeurs
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     IYC     : CHPOINT "CENTRE" contenat les fractions massiques ;
C
C     IGRYC   : CHPOINT "CENTRE" contenant les gradient des fr.mass ;
C
C     IALPC   : CHPOINT "CENTRE" contenant les limiteurs des gradients
C               des fr.mass. ;
C
C     ISCAC   : CHPOINT "CENTRE" contenat les scalaires passifs ;
C
C     IGRSC   : CHPOINT "CENTRE" contenant les gradient des scalaires passifs;
C
C     IALPC   : CHPOINT "CENTRE" contenant les limiteurs des gradients
C               des  scalaires passifs ;
C
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     IYF     : MCHAML  defini sur le MELEME de pointeur IFACEL,
C               contenant les  fractions massiques;
C
C     ISCAF   : MCHAML  defini sur le MELEME de pointeur IFACEL,
C               contenant les  scalaires passifs;
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 y 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 21.12.98.
C
C               17.02.2000: transport des scalaires passifs
C
C************************************************************************
C
C
C ATTENTION: Cet programme marche si le MAILLAGE est convex;
C            si non il faut changer l'algoritme 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**** Variables de COOPTIO
C
C      INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
C     &        ,IOTER,   IOLEC,  IOIMP,     IOCAR,  IOACQ
C     &        ,IOPER,   IOSGB,  IOGRA,     IOSAU,  IORES
C     &        ,IECHO,   IIMPI,  IOSPI
C     &        ,IDIM
CC     &        ,MCOORD
C     &        ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
C     &        ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
C     &        ,NORINC,NORVAL,NORIND,NORVAD
C     &        ,NUCROU, IPSAUV
C
C**** Les variables
C
      IMPLICIT INTEGER(I-N)
      INTEGER ICEN, IFACE, IFACEL, INORM,IROC, IGRROC, IALROC
     &       , IVITC, IGRVC, IALVC
     &       , IPC ,IGRPC, IALPC
     &       , IYC, IGRYC, IALYC
     &       , ISCAC, IGRSC, IALSC
     &       , IGAMC
     &       , IROF, IVITF, IPF, IYF, ISCAF
     &       , IGEOM, NFAC, NCEN
     &       , N1PTEL, N1EL, N2PTEL, N2EL, N2, N1, N3, L1, NLCE
     &       , NLCF, NGCF, NGCEG, NLCEG, NGCED, NLCED, NGCF1
     &       , IDIMP1, INDCEL, I1, NESP, NSCA, NMA
      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, UXD, UYD, UND, UTD
     &       , VALCEL, DCEL, ALCEL
     &       , DROX, DROY, DUXX, DUXY, DUYX, DUYY, DPX, DPY
     &       , DRO, DUX, DUY, DP
     &       , ALPHA,  DYMAS, SUMY
      CHARACTER*(40) MESERR
      CHARACTER*(8) TYPE
      LOGICAL LOGAN,LOGNEG, LOGBOR, LOGTEM, LOGI1, LOGI2
     &       ,LGAMC,LYC,LSCAC
C
C**** Les Includes
C
-INC SMCOORD
 
-INC PPARAM
-INC CCOPTIO
-INC SMCHPOI
      POINTEUR MPROC.MPOVAL,  MPGRR.MPOVAL,
     &         MPVITC.MPOVAL, MPGRV.MPOVAL,
     &         MPPC.MPOVAL,   MPGRP.MPOVAL,
     &         MPYC.MPOVAL,   MPGRY.MPOVAL,
     &         MPSCAC.MPOVAL, MPGRS.MPOVAL,
     &         MPGAMC.MPOVAL, MPNORM.MPOVAL,
     &         MPROP.MPOVAL,  MPPP.MPOVAL, MPVITP.MPOVAL,
     &         MPYP.MPOVAL,   MPSCAP.MPOVAL
-INC SMCHAML
      POINTEUR MCHAMY.MCHAML, MCHAMS.MCHAML
      POINTEUR MELVNX.MELVAL, MELVNY.MELVAL,
     &         MELT1X.MELVAL, MELT1Y.MELVAL
      POINTEUR MELVUN.MELVAL, MELVUT.MELVAL
      POINTEUR MELRO.MELVAL, MELP.MELVAL
      POINTEUR MELALR.MPOVAL,
     &         MELALV.MPOVAL,
     &         MELALP.MPOVAL,
     &         MELALY.MPOVAL,
     &         MELALS.MPOVAL
-INC SMLENTI
-INC SMELEME
C
C**** Segments des fractions massiques gauche et droit
C
      SEGMENT FRAMAS
         REAL*8 FRAMG(NMA), FRAMD(NMA)
      ENDSEGMENT
      POINTEUR SCALPA.FRAMAS
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     fract.mass + 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)
      IF(LGAMC) 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**** Les 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 (si LGAMC)
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
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 repair 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) = 0
      MCHEL1.IFOCHE = IFOUR
C
C**** Valeurs de vitesse definies par respect au repair 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) = 0
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
      SEGDES MCHAM1
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
      SEGDES MCHEL1
      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
      SEGDES MCHAM1
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) = 0
      MCHEL2.IFOCHE = IFOUR
      N2 = 1
      SEGINI MCHAM1
      MCHEL2.ICHAML(1) = MCHAM1
      SEGDES MCHEL2
      MCHAM1.NOMCHE(1) = 'SCAL    '
      MCHAM1.TYPCHE(1) = 'REAL*8          '
      SEGINI MELRO
      MCHAM1.IELVAL(1) = MELRO
      SEGDES MCHAM1
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
      SEGDES MCHEL2
      SEGINI MELP
      MCHAM2.IELVAL(1) = MELP
      SEGDES MCHAM2
C
C**** Le CHAMELEM des fractions massiques
C
      IF(LYC)THEN
         MCHPO1 = IYC
         SEGACT MCHPO1
         MSOUP1 = MCHPO1.IPCHP(1)
         SEGDES MCHPO1
         SEGACT MSOUP1
         MPYC = MSOUP1.IPOVAL
         SEGACT MPYC
         NESP = MPYC.VPOCHA(/2)
         MCHEL1 = IROF
         SEGINI, MCHEL2 = MCHEL1
         IYF = MCHEL2
         MCHEL2.TITCHE = 'Y            '
         N2 = NESP
         SEGINI MCHAMY
         MCHEL2.ICHAML(1) = MCHAMY
         SEGDES MCHEL2
         N1EL = NFAC
         N1PTEL = 3
         N2EL = 0
         N2PTEL = 0
         DO I1 = 1, NESP
            SEGINI MELVA1
            MCHAMY.IELVAL(I1) = MELVA1
            MCHAMY.NOMCHE(I1) = MSOUP1.NOCOMP(I1)
            MCHAMY.TYPCHE(I1) = 'REAL*8          '
         ENDDO
         SEGDES MSOUP1
C
         CALL LICHT(IGRYC,  MPGRY  , TYPE, IGEOM)
         CALL LICHT(IALYC,  MELALY , TYPE, IGEOM)
         NMA = NESP
         SEGINI FRAMAS
         IF(LOGTEM)THEN
            SEGINI, MPYP = MPYC
         ELSE
            MPYP = MPYC
         ENDIF
      ELSE
         NESP = 0
         IYF = 0
      ENDIF
C
C**** Le CHAMELEM des scalaires passifs
C
      IF(LSCAC)THEN
         MCHPO1 = ISCAC
         SEGACT MCHPO1
         MSOUP1 = MCHPO1.IPCHP(1)
         SEGDES MCHPO1
         SEGACT MSOUP1
         MPSCAC = MSOUP1.IPOVAL
         SEGACT MPSCAC
         NSCA = MPSCAC.VPOCHA(/2)
         MCHEL1 = IROF
         SEGINI, MCHEL2 = MCHEL1
         ISCAF = MCHEL2
         MCHEL2.TITCHE = 'SCALPASS     '
         N2 = NSCA
         SEGINI MCHAMS
         MCHEL2.ICHAML(1) = MCHAMS
         SEGDES MCHEL2
         N1EL = NFAC
         N1PTEL = 3
         N2EL = 0
         N2PTEL = 0
         DO I1 = 1, NSCA
            SEGINI MELVA1
            MCHAMS.IELVAL(I1) = MELVA1
            MCHAMS.NOMCHE(I1) = MSOUP1.NOCOMP(I1)
            MCHAMS.TYPCHE(I1) = 'REAL*8          '
         ENDDO
         SEGDES MSOUP1
C
         CALL LICHT(IGRSC,  MPGRS  , TYPE, IGEOM)
         CALL LICHT(IALSC,  MELALS , TYPE, IGEOM)
         NMA = NSCA
         SEGINI SCALPA
         IF(LOGTEM)THEN
            SEGINI, MPSCAP = MPSCAC
         ELSE
            MPSCAP = MPSCAC
         ENDIF
      ELSE
         NSCA = 0
         ISCAF = 0
      ENDIF
C
C**** Donc on a aussi actives le chpoints de fractions massiques
C
C     SEGACT MPYC
C     SEGACT MPGRY
C     SEGACT MPIALY
C
C     SEGACT MPSCAC
C     SEGACT MPGRS
C     SEGACT MPIALS
C
C
C***********************************************************************
C********* PREDICTION **************************************************
C***********************************************************************
C
C**** Prediction avec gradients limités
C
C
      IF(LOGTEM)THEN
C
         IPT3 = ICEN
         SEGACT IPT3
         NCEN = IPT3.NUM(/2)
         DO NLCE = 1, NCEN
            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)*MELALR.VPOCHA(NLCE,1)
            DROY = MPGRR.VPOCHA(NLCE,2)*MELALR.VPOCHA(NLCE,1)
            DUXX = MPGRV.VPOCHA(NLCE,1)*MELALV.VPOCHA(NLCE,1)
            DUXY = MPGRV.VPOCHA(NLCE,2)*MELALV.VPOCHA(NLCE,1)
            DUYX = MPGRV.VPOCHA(NLCE,3)*MELALV.VPOCHA(NLCE,2)
            DUYY = MPGRV.VPOCHA(NLCE,4)* MELALV.VPOCHA(NLCE,2)
            DPX  = MPGRP.VPOCHA(NLCE,1)* MELALP.VPOCHA(NLCE,1)
            DPY  = MPGRP.VPOCHA(NLCE,2)* MELALP.VPOCHA(NLCE,1)
            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
            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
            DO I1 = 1, NESP, 1
               INDCEL = 2 * I1 - 1
               ALPHA = MELALY.VPOCHA(NLCE,I1)
               DYMAS = UXG * MPGRY.VPOCHA(NLCE,INDCEL) * ALPHA +
     &              UYG * MPGRY.VPOCHA(NLCE,INDCEL+1) * ALPHA
               MPYP.VPOCHA(NLCE,I1) = MPYC.VPOCHA(NLCE,I1) -
     &              DELTAT * DYMAS
            ENDDO
            DO I1 = 1, NSCA, 1
               INDCEL = 2 * I1 - 1
               ALPHA = MELALS.VPOCHA(NLCE,I1)
               DYMAS = UXG * MPGRS.VPOCHA(NLCE,INDCEL) * ALPHA +
     &              UYG * MPGRS.VPOCHA(NLCE,INDCEL+1) * ALPHA
               MPSCAP.VPOCHA(NLCE,I1) = MPSCAC.VPOCHA(NLCE,I1) -
     &              DELTAT * DYMAS
            ENDDO
         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 pre321.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 pre321.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              y_i \in (0,1)
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
            LOGI2 = .FALSE.
            SUMY = 0.0D0
            DO I1 = 1, NESP, 1
               INDCEL = 2 * I1 - 1
               VALCEL = MPYP.VPOCHA(NLCEG,I1)
               ALCEL =  MELALY.VPOCHA(NLCEG, I1)
               DCEL = MPGRY.VPOCHA(NLCEG, INDCEL)*DXG +
     &              MPGRY.VPOCHA(NLCEG,INDCEL + 1 )*DYG
               ALCEL = VALCEL + ALCEL * DCEL
               SUMY = SUMY + ALCEL
               LOGI2 = LOGI2 .OR. (ALCEL .LT. 0.0D0)
               FRAMAS.FRAMG(I1) = ALCEL
            ENDDO
            LOGI2 = LOGI2 .OR. (SUMY .GT. 1.0D0)
C
            DO I1 = 1, NSCA, 1
               INDCEL = 2 * I1 - 1
               VALCEL = MPSCAP.VPOCHA(NLCEG,I1)
               ALCEL =  MELALS.VPOCHA(NLCEG, I1)
               DCEL = MPGRS.VPOCHA(NLCEG, INDCEL)*DXG +
     &              MPGRS.VPOCHA(NLCEG,INDCEL + 1 )*DYG
               ALCEL = VALCEL + ALCEL * DCEL
               SCALPA.FRAMG(I1) = ALCEL
            ENDDO
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) .OR. LOGI2
C
            IF(LOGI1)THEN
C
C************* Premier ordre en espace local
C
               ROG = MPROC.VPOCHA(NLCEG,1)
               ROD = ROG
               PG  = MPPC.VPOCHA(NLCEG,1)
               PD = PG
               UNG = MPVITC.VPOCHA(NLCEG,1)*CNX +
     &              MPVITC.VPOCHA(NLCEG,2)*CNY
               UTG = MPVITC.VPOCHA(NLCEG,1)*CTX +
     &              MPVITC.VPOCHA(NLCEG,2)*CTY
               UND = -1.0D0 * UNG
               UTD = UTG
               DO I1 = 1, NESP, 1
                  FRAMAS.FRAMG(I1) = MPYC.VPOCHA(NLCEG,I1)
                  FRAMAS.FRAMD(I1) = FRAMAS.FRAMG(I1)
               ENDDO
               DO I1 = 1, NSCA, 1
                  SCALPA.FRAMG(I1) = MPSCAC.VPOCHA(NLCEG,I1)
                  SCALPA.FRAMD(I1) = SCALPA.FRAMG(I1)
               ENDDO
            ELSE
C
C********** Son etat droite
C
               ROD = ROG
               PD = PG
               UND = -1.0D0 * UNG
               UTD = UTG
               DO I1 = 1, NESP, 1
                  FRAMAS.FRAMD(I1) = FRAMAS.FRAMG(I1)
               ENDDO
               DO I1 = 1, NSCA, 1
                  SCALPA.FRAMD(I1) = SCALPA.FRAMG(I1)
               ENDDO
            ENDIF
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
            LOGI2 = .FALSE.
            SUMY = 0.0D0
            DO I1 = 1, NESP, 1
               INDCEL = 2 * I1 - 1
               VALCEL = MPYP.VPOCHA(NLCEG,I1)
               ALCEL =  MELALY.VPOCHA(NLCEG, I1)
               DCEL = MPGRY.VPOCHA(NLCEG, INDCEL)*DXG +
     &              MPGRY.VPOCHA(NLCEG,INDCEL + 1 )*DYG
               ALCEL = VALCEL + ALCEL * DCEL
               SUMY = SUMY + ALCEL
               LOGI2 = LOGI2 .OR. (ALCEL .LT. 0.0D0)
               FRAMAS.FRAMG(I1) = ALCEL
            ENDDO
            LOGI2 = LOGI2 .OR. (SUMY .GT. 1.0D0)
C
            DO I1 = 1, NSCA, 1
               INDCEL = 2 * I1 - 1
               VALCEL = MPSCAP.VPOCHA(NLCEG,I1)
               ALCEL =  MELALS.VPOCHA(NLCEG, I1)
               DCEL = MPGRS.VPOCHA(NLCEG, INDCEL)*DXG +
     &              MPGRS.VPOCHA(NLCEG,INDCEL + 1 )*DYG
               ALCEL = VALCEL + ALCEL * DCEL
               SCALPA.FRAMG(I1) = ALCEL
            ENDDO
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) .OR. LOGI2
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
               DO I1 = 1, NESP, 1
                  FRAMAS.FRAMG(I1) = MPYC.VPOCHA(NLCEG,I1)
               ENDDO
               DO I1 = 1, NSCA, 1
                  SCALPA.FRAMG(I1) = MPSCAC.VPOCHA(NLCEG,I1)
               ENDDO
            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
            LOGI2 = .FALSE.
            SUMY = 0.0D0
            DO I1 = 1, NESP, 1
               INDCEL = 2 * I1 - 1
               VALCEL = MPYP.VPOCHA(NLCED,I1)
               ALCEL =  MELALY.VPOCHA(NLCED, I1)
               DCEL = MPGRY.VPOCHA(NLCED, INDCEL)*DXD +
     &              MPGRY.VPOCHA(NLCED,INDCEL + 1 )*DYD
               ALCEL = VALCEL + ALCEL * DCEL
               SUMY = SUMY + ALCEL
               LOGI2 = LOGI2 .OR. (ALCEL .LT. 0.0D0)
               FRAMAS.FRAMD(I1) = ALCEL
            ENDDO
            LOGI2 = LOGI2 .OR. (SUMY .GT. 1.0D0)
C
            DO I1 = 1, NSCA, 1
               INDCEL = 2 * I1 - 1
               VALCEL = MPSCAP.VPOCHA(NLCED,I1)
               ALCEL =  MELALS.VPOCHA(NLCED, I1)
               DCEL = MPGRS.VPOCHA(NLCED, INDCEL)*DXD +
     &              MPGRS.VPOCHA(NLCED,INDCEL + 1 )*DYD
               ALCEL = VALCEL + ALCEL * DCEL
               SCALPA.FRAMD(I1) = ALCEL
            ENDDO
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) .OR. LOGI2
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
               DO I1 = 1, NESP, 1
                  FRAMAS.FRAMD(I1) = MPYC.VPOCHA(NLCED,I1)
               ENDDO
               DO I1 = 1, NSCA, 1
                  SCALPA.FRAMD(I1) = MPSCAC.VPOCHA(NLCED,I1)
               ENDDO
            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
         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
         DO I1 = 1, NESP, 1
            MELVA1 = MCHAMY.IELVAL(I1)
            MELVA1.VELCHE(1,NLCF) = FRAMAS.FRAMG(I1)
            MELVA1.VELCHE(3,NLCF) = FRAMAS.FRAMD(I1)
         ENDDO
         DO I1 = 1, NSCA, 1
            MELVA1 = MCHAMS.IELVAL(I1)
            MELVA1.VELCHE(1,NLCF) = SCALPA.FRAMG(I1)
            MELVA1.VELCHE(3,NLCF) = SCALPA.FRAMD(I1)
         ENDDO
      ENDDO
C
C**** Desactivation des SEGMENTs
C
      SEGDES IPT1
      SEGDES IPT2
C
C**** Le MPOVALs 'Prediction' sont detruits (si existentes)
C
      IF(LOGTEM)THEN
         SEGSUP MPROP
         SEGSUP MPVITP
         SEGSUP MPPP
         IF(LYC) SEGSUP MPYP
         IF(LSCAC) SEGSUP MPSCAP
         SEGDES MPGAMC
      ENDIF
C
      SEGDES MPROC
      SEGDES MPGRR
      SEGDES MELALR
      SEGDES MPVITC
      SEGDES MPGRV
      SEGDES MELALV
      SEGDES MPPC
      SEGDES MPGRP
      SEGDES MELALP
      IF(LYC)THEN
         SEGDES MPYC
         SEGDES MPGRY
         SEGDES MELALY
         DO I1 = 1, NESP, 1
            MELVA1 = MCHAMY.IELVAL(I1)
            SEGDES MELVA1
         ENDDO
         SEGDES MCHAMY
         SEGSUP FRAMAS
      ENDIF
      IF(LSCAC)THEN
         SEGDES MPSCAC
         SEGDES MPGRS
         SEGDES MELALS
         DO I1 = 1, NSCA, 1
            MELVA1 = MCHAMS.IELVAL(I1)
            SEGDES MELVA1
         ENDDO
         SEGDES MCHAMS
         SEGSUP SCALPA
      ENDIF
      SEGDES MPNORM
C
      SEGDES MELRO
      SEGDES MELP
      SEGDES MELVUN
      SEGDES MELVUT
      SEGDES MELVNX
      SEGDES MELVNY
      SEGDES MELT1X
      SEGDES MELT1Y
C
C**** Destruction du MELNTI correspondance local/global
C
      SEGSUP MLENT1
C
 9999 CONTINUE
C
      RETURN
      END
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

© Cast3M 2003 - Tous droits réservés.
Mentions légales