konfl1
C KONFL1 SOURCE CB215821 20/11/25 13:32:01 10792 & IROF,IVITF,IPF,IGAMF,IFRMAF, & ICHPSU,ICHPDI, & MELEMC,MELEMF,MELEFE,MELLIM, & IUINF,IUPRI, & ICHFLU,DT, & LOGNC,LOGAN,MESERR) C************************************************************************ C C PROJET : CASTEM 2000 C C NOM : KONFL1 C C DESCRIPTION : Voir KON11 C C Cas deux dimensions, gaz "calorically perfect" C C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI) C C AUTEUR : A. BECCANTINI, DRN/DMT/SEMT/TTMF C C************************************************************************ C C ENTREES C C C 1) PARAMETRES C C LOGME : (LOGICAL); .TRUE. -> MULTI-ESPECES C .FALSE. -> MONO-ESPECE C C INDMET : 1 Godunov C C 2 van Leer FVS C C 3 van Leer Hanel FVS C C 4 HUS (Van Leer FVS + Osher FDS) C C 5 HUS (Van Leer - Hanel FVS + Osher FDS) C C 6 AUSM+ C C 7 ROE (GAMMA = 0.5 * (GAMG + GAMD)) C C 8 SS C C 9 AUSM+ Low Mach C C 10 Rusanvov C C 11 Rusanvov Low Mach C C 12 Centrée C C 13 RoeLM C C 14 HLLC C C 15 HLLCLM C C 2) Pointeurs des MCHAMLs C C IROF : MCHAML sur "FACEL" contenant la masse volumique C ("gauche" et "droite"); C C IVITF : MCHAML sur "FACEL" contenant la vitesse dans le repaire C local (n,t) et les cosinus directeurs des repaire local; C C IPF : MCHAML sur "FACEL" contenant la pression; C C IGAMF : MCHAML sur "FACEL" contenant le gamma; C C IFRAMAF : MCHAML sur "FACEL", contenant les fractions massiques C si LOGME = .TRUE.; C LOGME = .FALSE. -> IFRAMAF = 0 C C C 3) Pointeurs de CHPOINTs de la table DOMAINE C C ICHPSU : CHPOINT "FACE" contenant la surface des faces C C ICHPDI : CHPOINT "CENTRE" contenant le diametre minimum C de chaque element C C C 4) Pointeurs de MELEME C C MELEMC : MELEME 'CENTRE' du SPG des CENTRES C C MELEMF : MELEME 'FACE' du SPG des FACES C C MELEFE : MELEME 'FACEL' du connectivité FACES -> ELEM C C MELLIM : MELEME SPG des conditions aux bords C Sur ces points on ne calcule pas le flux! C C 5) Cas Bas MACH C C IUINF : CHPOINT, one component, "cut-off velocity" C 0 if no Bas MACH C C IUPRI : CHPOINT, one component, "minimum reference velocity" C 0 if no BAs Mach C C C SORTIES C C ICHFLU : pointeur de CHPOINT "FACE" des flux aux interfaces: C C DT : pas de temps pour le respect de la CFL-like condition C DT < DIAMMIN /2 /max(Lambda_i) C En maillage regulier cette condition garantie la C non-interaction des ondes C C LOGNC : (LOGICAL): si .TRUE. la methode de Newton-Rapson, utilisée C dans pour la solution du probleme Riemann exacte ou dans C l'algorithm HUS, n'a pas bien marchéee; MESERR = 'Goudunov' C ou 'HUS'. C C LOGAN : (LOGICAL): si .TRUE. une anomalie à été detectée C C MESERR : pour l'ecriture des messages d'erreurs C C************************************************************************ C C HISTORIQUE (Anomalies et modifications éventuelles) C C HISTORIQUE : C 01/05/2005 par T. KLOCZKO (DEN/DM2S/SFME/LTMF) C Ajout d'une méthode de décentrement pour les écoulements C bas-Mach (ROE Low-Mach) C C 04/05/2010 par A. BECCANTINI (DEN/DM2S/SFME/LTMF) C Ajout des méthodes de décentrement HLLC and HLLCLM C (pour les écoulements bas-Mach). C C C************************************************************************ C C C N.B.: On suppose qu'on a déjà controllé RO, P > 0 C GAMMA \in (1,3) C Y \in (0,1) C Si non il faut le faire!!! 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 C & ,MCOORD C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU C & ,NORINC,NORVAL,NORIND,NORVAD C & ,NUCROU, IPSAUV, IFICLE, IPREFI C & ,IREFOR, ISAFOR C C IMPLICIT INTEGER(I-N) INTEGER I1 & ,INDMET & ,IROF,IVITF,IPF,IGAMF,IFRMAF & ,ICHPSU,ICHPDI,MELEMC,MELEMF,MELEFE & ,IUINF,IUPRI & ,IGEOMC,IGEOMF & ,ICHFLU & ,NESP, NFAC & ,NLCF, NGCEG, NGCED, NLCEG, NLCED & ,NGCF, NLCF1, SPG1, SPG2, NLFL REAL*8 DT, UNSDT, CELLT & , ROG, UNG, UTG, PG, GAMG & , ROD, UND, UTD, PD, GAMD & , SURF,CNX, CNY, CTX , CTY & , CELL, DIAMG, DIAMD, DIAM & , ASON, LAMBDA & , V_INF C LOGICAL LOGME, LOGNC, LOGAN CHARACTER*(40) MESERR CHARACTER*(8) TYPE C C**** LES INCLUDES C -INC PPARAM -INC CCOPTIO -INC SMCHAML POINTEUR MELVNX.MELVAL, MELVNY.MELVAL, & MELT1X.MELVAL, MELT1Y.MELVAL POINTEUR MELVUN.MELVAL, MELVUT.MELVAL POINTEUR MELRO.MELVAL, MELP.MELVAL, & MELGAM.MELVAL -INC SMCHPOI POINTEUR MPOVSU.MPOVAL, MPOVDI.MPOVAL & , MPOFLU.MPOVAL & , MPUINF.MPOVAL, MPUPRI.MPOVAL POINTEUR MCHAMY.MCHAML -INC SMELEME POINTEUR MELLIM.MELEME -INC SMLMOTS -INC SMLENTI POINTEUR MLELIM.MLENTI C C**** Les fractiones massiques. C SEGMENT FRAMAS REAL*8 YET(NESP) ENDSEGMENT POINTEUR FRAMAG.FRAMAS, FRAMAD.FRAMAS C C**** Les flux aux interface dans le repaire (n,t) C SEGMENT IFLUX ENDSEGMENT POINTEUR IFLU1.IFLUX, IFLU2.IFLUX C C**** KRIPAD pour la correspondance global/local des conditions limits C c SEGACT MELLIM C C**** Initialisation des MCHAMLs C C**** Masse volumique C MCHEL1 = IROF SEGACT MCHEL1 MCHAM1 = MCHEL1.ICHAML(1) SEGACT MCHAM1 MELRO = MCHAM1.IELVAL(1) SEGDES MCHEL1 SEGDES MCHAM1 C C**** Pression C MCHEL1 = IPF SEGACT MCHEL1 MCHAM1 = MCHEL1.ICHAML(1) SEGACT MCHAM1 MELP = MCHAM1.IELVAL(1) SEGDES MCHEL1 SEGDES MCHAM1 C C**** Gamma C MCHEL1 = IGAMF SEGACT MCHEL1 MCHAM1 = MCHEL1.ICHAML(1) SEGACT MCHAM1 MELGAM = MCHAM1.IELVAL(1) SEGDES MCHEL1 SEGDES MCHAM1 C C**** Vitesse et cosinus directeurs du repere (n,t) C MCHEL1 = IVITF SEGACT MCHEL1 C C**** La vitesse a comme SPG MELEFE C Le cosinus directeurs ont comme SPG MELEMF C C MCHAM1 -> Cosinus directeurs C MCHAM2 -> Vitesse C SPG1 = MCHEL1.IMACHE(1) SPG2 = MCHEL1.IMACHE(2) IF((SPG1 .EQ. MELEMF) .AND. (SPG2 .EQ. MELEFE))THEN MCHAM1 = MCHEL1.ICHAML(1) MCHAM2 = MCHEL1.ICHAML(2) ELSEIF((SPG1 .EQ. MELEFE) .AND. (SPG2 .EQ. MELEMF))THEN MCHAM1 = MCHEL1.ICHAML(2) MCHAM2 = MCHEL1.ICHAML(1) ELSE LOGAN = .TRUE. GOTO 9999 ENDIF SEGACT MCHAM1 MELVNX = MCHAM1.IELVAL(1) MELVNY = MCHAM1.IELVAL(2) MELT1X = MCHAM1.IELVAL(3) MELT1Y = MCHAM1.IELVAL(4) SEGDES MCHAM1 SEGACT MCHAM2 MELVUN = MCHAM2.IELVAL(1) MELVUT = MCHAM2.IELVAL(2) SEGDES MCHAM2 SEGDES MCHEL1 C C**** Fractions massiques C IF(LOGME)THEN MCHEL1 = IFRMAF SEGACT MCHEL1 MCHAMY = MCHEL1.ICHAML(1) SEGACT MCHAMY C C******* Numero d'especes dans les equations d'Euler C NESP = MCHAMY.IELVAL(/1) DO I1 = 1, NESP MELVA1 = MCHAMY.IELVAL(I1) SEGACT MELVA1 ENDDO SEGINI FRAMAG SEGINI FRAMAD SEGDES MCHEL1 ELSE C C******* Definition minimale de YET, necessaire pour transmetre YET aux C subroutines FORTRAN qui calculent les flux C NESP = 1 SEGINI FRAMAG SEGINI FRAMAD NESP = 0 ENDIF C C**** Initialisation des MELEMEs C C 'CENTRE', 'FACEL' C IPT2 = MELEFE SEGACT IPT2 NFAC = IPT2.NUM(/2) C C**** KRIPAD pour la correspondance global/local de centre C 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 C**** KRIPAD pour la correspondance global/local de 'FACE' C C SEGACT MELNT2 C C**** Initialisation de flux C SEGINI IFLU1 SEGINI IFLU2 C C**** IFLU2 = segment de travail en FLUVLH; c'est plus rapide le definir ici C C C**** CHPOINTs de la table DOMAINE C C C**** LICHT active les MPOVALs en *MOD C C i.e. C C SEGACT MPOVSU*MOD C SEGACT MPOVDI*MOD C IF(IGEOMF .NE. MELEMF)THEN MESERR = 'Il ne faut pas jouer avec la console. ' LOGAN = .TRUE. GOTO 9999 ENDIF IF(IGEOMC .NE. MELEMC)THEN MESERR = 'Il ne faut pas jouer avec la console. ' LOGAN = .TRUE. GOTO 9999 ENDIF C C**** Les FLUX aux face C C La densité C C C SEGACT MPOFLU*MOD C IF(IGEOMC .NE. MELEMC)THEN MESERR = 'Il ne faut pas jouer avec la console. ' LOGAN = .TRUE. GOTO 9999 ENDIF C C**** Activation des MCHAMLs C SEGACT MELRO SEGACT MELP SEGACT MELGAM SEGACT MELVUN SEGACT MELVUT SEGACT MELVNX SEGACT MELVNY SEGACT MELT1X SEGACT MELT1Y C C**** Bas Mach C IF(IUINF .NE. 0)THEN C SEGACT MPUPRI*MOD C SEGACT MPUINF*MOD ELSE MPUPRI=0 MPUINF=0 ENDIF C C**** Initialisation de 1/DT C UNSDT = 0.0D0 C C**** BOUCLE SUR FACEL pour le calcul du FLUX C DO NLCF = 1, NFAC C C******* NLCF = numero local du centre de facel C NGCF = numero global du centre de facel C NLCF1 = numero local 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 = IPT2.NUM(1,NLCF) NGCED = IPT2.NUM(3,NLCF) NGCF = IPT2.NUM(2,NLCF) NLCF1 = MLENT2.LECT(NGCF) NLCEG = MLENT1.LECT(NGCEG) NLCED = MLENT1.LECT(NGCED) NLFL = MLELIM.LECT(NGCF) C C******* NLCF != NLCF1 -> l'auteur (MOI) n'a rien compris. C IF(NLCF .NE. NLCF1)THEN MESERR = 'Il ne faut pas jouer avec la console. ' LOGAN = .TRUE. GOTO 9999 ENDIF IF(NLFL .EQ. 0)THEN C C******* Recuperation des Etats "gauche" et "droite" C ROG = MELRO.VELCHE(1,NLCF) UNG = MELVUN.VELCHE(1,NLCF) UTG = MELVUT.VELCHE(1,NLCF) PG = MELP.VELCHE(1,NLCF) GAMG = MELGAM.VELCHE(1,NLCF) C ROD = MELRO.VELCHE(3,NLCF) UND = MELVUN.VELCHE(3,NLCF) UTD = MELVUT.VELCHE(3,NLCF) PD = MELP.VELCHE(3,NLCF) GAMD = MELGAM.VELCHE(3,NLCF) C CNX = MELVNX.VELCHE(1,NLCF) CNY = MELVNY.VELCHE(1,NLCF) CTX = MELT1X.VELCHE(1,NLCF) CTY = MELT1Y.VELCHE(1,NLCF) C C******* Le fractiones massiques C IF(LOGME)THEN DO I1 = 1, NESP MELVA1 = MCHAMY.IELVAL(I1) FRAMAG.YET(I1) = MELVA1.VELCHE(1,NLCF) FRAMAD.YET(I1) = MELVA1.VELCHE(3,NLCF) ENDDO ENDIF C C******* On a defini (ROg,ROUNg,ROUTg,Pg,(Yg)), (ROd,ROUNd,ROUTd,Pd,(Yd)) C et on a déjà verifié ROg, ROd, Pg, Pd > 0 et 0<Y_i<1 C C C******* Calcul du flux aux interfaces C IF(INDMET .EQ. 1)THEN C C******* GODUNOV C FLURIE en FORTRAN STANDARD C CELLT = C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT, & LOGNC,LOGAN,MESERR) C ELSEIF(INDMET .EQ. 2)THEN C C******* Van Leer FVS C C N.B: FLUXVL en FORTRAN pure C FRAMAG.YET = table d'un pointeur -> table C La meme chose pour FRAMAD.YET, IFLU1.FLUX, C IFLU2.FLUX C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & CELLT) ELSEIF(INDMET .EQ. 3)THEN C C******* Van Leer - Hanel FVS C C N.B: FLUVLH en FORTRAN pure C FRAMAG.YET = table d'un pointeur -> table C La meme chose pour FRAMAD.YET, IFLU1.FLUX, C IFLU2.FLUX C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & CELLT) ELSEIF(INDMET .EQ. 4)THEN C C******* HUS (Van Leer FVS + Osher FDS) C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & CELLT, & LOGNC,MESERR,LOGAN) ELSEIF(INDMET .EQ. 5)THEN C C******* HUS (Van Leer - Hanel FVS + Osher FDS) C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & CELLT, & LOGNC,MESERR,LOGAN) ELSEIF(INDMET .EQ. 6)THEN C C******* AUSM+ (Liou) C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 7)THEN C C******* Roe Scheme C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 8)THEN C C******* Choc-Choc C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT, & LOGNC,LOGAN,MESERR) ELSEIF(INDMET .EQ. 9)THEN C C******* AUSM+(P) for low Mach number flows (Edwards and Liou) C V_INF=MAX(MPUINF.VPOCHA(NLCEG,1), & MPUINF.VPOCHA(NLCED,1), & MPUPRI.VPOCHA(NLCEG,1), & MPUPRI.VPOCHA(NLCED,1)) & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & V_INF, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 10)THEN C C******* Rusanov scheme C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 11)THEN C C******* Rusanov scheme for Low-Mach number flows C V_INF=MAX(MPUINF.VPOCHA(NLCEG,1), & MPUINF.VPOCHA(NLCED,1), & MPUPRI.VPOCHA(NLCEG,1), & MPUPRI.VPOCHA(NLCED,1)) & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & V_INF, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 12)THEN C C******* Centered scheme C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 13)THEN C C******* Roe-Turkel scheme for Low-Mach number flows C V_INF=MAX(MPUINF.VPOCHA(NLCEG,1), & MPUINF.VPOCHA(NLCED,1), & MPUPRI.VPOCHA(NLCEG,1), & MPUPRI.VPOCHA(NLCED,1)) & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & V_INF, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 14)THEN C C******* HLLC C & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 15)THEN C C******* HLLC-Turkel scheme for Low-Mach number flows C V_INF=MAX(MPUINF.VPOCHA(NLCEG,1), & MPUINF.VPOCHA(NLCED,1), & MPUPRI.VPOCHA(NLCEG,1), & MPUPRI.VPOCHA(NLCED,1)) & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & V_INF, & IFLU1.FLUX, & CELLT) ELSEIF(INDMET .EQ. 16)THEN C C******* AUSM+up for Low-Mach number flows C V_INF=MAX(MPUINF.VPOCHA(NLCEG,1), & MPUINF.VPOCHA(NLCED,1), & MPUPRI.VPOCHA(NLCEG,1), & MPUPRI.VPOCHA(NLCED,1)) & GAMG,ROG,PG,UNG,UTG, & GAMD,ROD,PD,UND,UTD, & FRAMAG.YET,FRAMAD.YET, & V_INF, & IFLU1.FLUX, & CELLT) ENDIF C IF(LOGAN) GOTO 9999 IF(LOGNC) GOTO 9999 C C******* Ecriture des flux C C FLUX(1) = RO Un RO Un C FLUX(2) = RO Un Un + P -> RO Un Ux + P CNX C FLUX(3) = RO Un Ut -> RO Un Uy + P CNY C FLUX(4) = RO Un Ht RO Un Ht C SURF = MPOVSU.VPOCHA(NLCF,1) MPOFLU.VPOCHA(NLCF,1) = MPOFLU.VPOCHA(NLCF,2) = MPOFLU.VPOCHA(NLCF,3) = MPOFLU.VPOCHA(NLCF,4) = IF(LOGME)THEN DO I1 = 1, NESP ENDDO ENDIF C C******* Calcul du pas du temps (CFL) C C C C****** a) etat a l'interface C C CELLT a la dimension du reciproque d'une vitesse C DIAMG = MPOVDI.VPOCHA(NLCEG,1) DIAMD = MPOVDI.VPOCHA(NLCED,1) DIAM = MIN(DIAMG,DIAMD) CELL = 1.0D0/(DIAM*CELLT) IF(CELL .GT. UNSDT)THEN UNSDT = CELL ENDIF C C****** b) etat gauche C ASON = SQRT(GAMG*PG/ROG) LAMBDA = ABS(UNG) + ASON CELL = LAMBDA / DIAM IF(CELL .GT. UNSDT)THEN UNSDT = CELL ENDIF C C****** C) etat droite C ASON = SQRT(GAMD*PD/ROD) LAMBDA = ABS(UND) + ASON CELL = LAMBDA / DIAM IF(CELL .GT. UNSDT)THEN UNSDT = CELL ENDIF ENDIF C C**** Fin boucle sur FACEL C ENDDO C C**** Pas du temps (condition de non interaction en 1D) C DT = 0.5D0 / UNSDT C C**** Desactivation des segments et C on detruit les MCHAMLs C C C**** SEGSUP FRAMAG C SEGSUP FRAMAD C C meme si LOGME = .FALSE. C SEGSUP FRAMAG SEGSUP FRAMAD C SEGSUP MLENT1 SEGSUP MLENT2 SEGDES IPT2 C SEGSUP IFLU1 SEGSUP IFLU2 C SEGDES MPOVSU SEGDES MPOVDI C SEGDES MPOFLU C SEGDES MELRO SEGDES MELP SEGDES MELGAM SEGDES MELVUN SEGDES MELVUT SEGDES MELVNX SEGDES MELVNY SEGDES MELT1X SEGDES MELT1Y C IF(MPUPRI .NE. 0)THEN SEGDES MPUPRI SEGDES MPUINF ENDIF C IF(LOGME) THEN DO I1 = 1, NESP MELVA1 = MCHAMY.IELVAL(I1) SEGDES MELVA1 ENDDO C SEGDES MCHAMY ENDIF SEGSUP MLELIM CC 9999 CONTINUE C RETURN END C
© Cast3M 2003 - Tous droits réservés.
Mentions légales