1 : C KONFL1 SOURCE BECC 10/05/05 21:15:13 6674
2 : SUBROUTINE KONFL1(LOGME,INDMET,
3 : & IROF,IVITF,IPF,IGAMF,IFRMAF,
4 : & ICHPSU,ICHPDI,
5 : & MELEMC,MELEMF,MELEFE,MELLIM,
6 : & IUINF,IUPRI,
7 : & ICHFLU,DT,
8 : & LOGNC,LOGAN,MESERR)
9 : C************************************************************************
10 : C
11 : C PROJET : CASTEM 2000
12 : C
13 : C NOM : KONFL1
14 : C
15 : C DESCRIPTION : Voir KON11
16 : C
16 : C Cas deux dimensions, gaz "calorically perfect"
17 : C
18 : C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
19 : C
20 : C AUTEUR : A. BECCANTINI, DRN/DMT/SEMT/TTMF
21 : C
22 : C************************************************************************
23 : C
24 : C ENTREES
25 : C
26 : C
27 : C 1) PARAMETRES
28 : C
29 : C LOGME : (LOGICAL); .TRUE. -> MULTI-ESPECES
30 : C .FALSE. -> MONO-ESPECE
31 : C
32 : C INDMET : 1 Godunov
33 : C
34 : C 2 van Leer FVS
35 : C
36 : C 3 van Leer Hanel FVS
37 : C
38 : C 4 HUS (Van Leer FVS + Osher FDS)
39 : C
40 : C 5 HUS (Van Leer - Hanel FVS + Osher FDS)
41 : C
42 : C 6 AUSM+
43 : C
44 : C 7 ROE (GAMMA = 0.5 * (GAMG + GAMD))
45 : C
46 : C 8 SS
47 : C
48 : C 9 AUSM+ Low Mach
49 : C
50 : C 10 Rusanvov
51 : C
52 : C 11 Rusanvov Low Mach
53 : C
54 : C 12 Centrée
55 : C
56 : C 13 RoeLM
57 : C
58 : C 14 HLLC
59 : C
60 : C 15 HLLCLM
61 : C
62 : C 2) Pointeurs des MCHAMLs
63 : C
64 : C IROF : MCHAML sur "FACEL" contenant la masse volumique
65 : C ("gauche" et "droite");
66 : C
67 : C IVITF : MCHAML sur "FACEL" contenant la vitesse dans le repaire
68 : C local (n,t) et les cosinus directeurs des repaire local;
69 : C
70 : C IPF : MCHAML sur "FACEL" contenant la pression;
71 : C
72 : C IGAMF : MCHAML sur "FACEL" contenant le gamma;
73 : C
74 : C IFRAMAF : MCHAML sur "FACEL", contenant les fractions massiques
75 : C si LOGME = .TRUE.;
76 : C LOGME = .FALSE. -> IFRAMAF = 0
77 : C
78 : C
79 : C 3) Pointeurs de CHPOINTs de la table DOMAINE
80 : C
81 : C ICHPSU : CHPOINT "FACE" contenant la surface des faces
82 : C
83 : C ICHPDI : CHPOINT "CENTRE" contenant le diametre minimum
84 : C de chaque element
85 : C
86 : C
87 : C 4) Pointeurs de MELEME
88 : C
89 : C MELEMC : MELEME 'CENTRE' du SPG des CENTRES
90 : C
91 : C MELEMF : MELEME 'FACE' du SPG des FACES
92 : C
93 : C MELEFE : MELEME 'FACEL' du connectivité FACES -> ELEM
94 : C
95 : C MELLIM : MELEME SPG des conditions aux bords
96 : C Sur ces points on ne calcule pas le flux!
97 : C
98 : C 5) Cas Bas MACH
99 : C
100 : C IUINF : CHPOINT, one component, "cut-off velocity"
101 : C 0 if no Bas MACH
102 : C
103 : C IUPRI : CHPOINT, one component, "minimum reference velocity"
104 : C 0 if no BAs Mach
105 : C
106 : C
107 : C SORTIES
108 : C
109 : C ICHFLU : pointeur de CHPOINT "FACE" des flux aux interfaces:
110 : C
111 : C DT : pas de temps pour le respect de la CFL-like condition
112 : C DT < DIAMMIN /2 /max(Lambda_i)
113 : C En maillage regulier cette condition garantie la
114 : C non-interaction des ondes
115 : C
116 : C LOGNC : (LOGICAL): si .TRUE. la methode de Newton-Rapson, utilisée
117 : C dans pour la solution du probleme Riemann exacte ou dans
118 : C l'algorithm HUS, n'a pas bien marchéee; MESERR = 'Goudunov'
119 : C ou 'HUS'.
120 : C
121 : C LOGAN : (LOGICAL): si .TRUE. une anomalie à été detectée
122 : C
123 : C MESERR : pour l'ecriture des messages d'erreurs
124 : C
125 : C************************************************************************
126 : C
127 : C HISTORIQUE (Anomalies et modifications éventuelles)
128 : C
129 : C HISTORIQUE :
130 : C 01/05/2005 par T. KLOCZKO (DEN/DM2S/SFME/LTMF)
131 : C Ajout d'une méthode de décentrement pour les écoulements
132 : C bas-Mach (ROE Low-Mach)
133 : C
134 : C 04/05/2010 par A. BECCANTINI (DEN/DM2S/SFME/LTMF)
135 : C Ajout des méthodes de décentrement HLLC and HLLCLM
136 : C (pour les écoulements bas-Mach).
137 : C
138 : C
139 : C************************************************************************
140 : C
141 : C
142 : C N.B.: On suppose qu'on a déjà controllé RO, P > 0
143 : C GAMMA \in (1,3)
144 : C Y \in (0,1)
145 : C Si non il faut le faire!!!
146 : C
147 : C************************************************************************
148 : C
149 : C
150 : C**** Variables de COOPTIO
151 : C
152 : C INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
153 : C & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ
154 : C & ,IOPER, IOSGB, IOGRA, IOSAU, IORES
155 : C & ,IECHO, IIMPI, IOSPI
156 : C & ,IDIM
157 : C & ,MCOORD
158 : C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
159 : C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
160 : C & ,NORINC,NORVAL,NORIND,NORVAD
161 : C & ,NUCROU, IPSAUV, IFICLE, IPREFI
162 : C & ,IREFOR, ISAFOR
163 : C
164 : C
165 : IMPLICIT INTEGER(I-N)
166 : INTEGER I1
167 : & ,INDMET
168 : & ,IROF,IVITF,IPF,IGAMF,IFRMAF
169 : & ,ICHPSU,ICHPDI,MELEMC,MELEMF,MELEFE
170 : & ,IUINF,IUPRI
171 : & ,IGEOMC,IGEOMF
172 : & ,ICHFLU
173 : & ,NESP, NFAC
174 : & ,NLCF, NGCEG, NGCED, NLCEG, NLCED
175 : & ,NGCF, NLCF1, SPG1, SPG2, NLFL
176 : REAL*8 DT, UNSDT, CELLT
177 : & , ROG, UNG, UTG, PG, GAMG
178 : & , ROD, UND, UTD, PD, GAMD
179 : & , SURF,CNX, CNY, CTX , CTY
180 : & , CELL, DIAMG, DIAMD, DIAM
181 : & , ASON, LAMBDA
182 : & , V_INF
183 : C
184 : LOGICAL LOGME, LOGNC, LOGAN
185 : CHARACTER*(40) MESERR
186 : CHARACTER*(8) TYPE
187 : C
188 : C**** LES INCLUDES
189 : C
190 : -INC CCOPTIO
191 : -INC SMCHAML
192 : POINTEUR MELVNX.MELVAL, MELVNY.MELVAL,
193 : & MELT1X.MELVAL, MELT1Y.MELVAL
194 : POINTEUR MELVUN.MELVAL, MELVUT.MELVAL
195 : POINTEUR MELRO.MELVAL, MELP.MELVAL,
196 : & MELGAM.MELVAL
197 : -INC SMCHPOI
198 : POINTEUR MPOVSU.MPOVAL, MPOVDI.MPOVAL
199 : & , MPOFLU.MPOVAL
200 : & , MPUINF.MPOVAL, MPUPRI.MPOVAL
201 : POINTEUR MCHAMY.MCHAML
202 : -INC SMELEME
203 : POINTEUR MELLIM.MELEME
204 : -INC SMLMOTS
205 : -INC SMLENTI
206 : POINTEUR MLELIM.MLENTI
207 : C
208 : C**** Les fractiones massiques.
209 : C
210 : SEGMENT FRAMAS
211 : REAL*8 YET(NESP)
212 : ENDSEGMENT
213 : POINTEUR FRAMAG.FRAMAS, FRAMAD.FRAMAS
214 : C
215 : C**** Les flux aux interface dans le repaire (n,t)
216 : C
217 : SEGMENT IFLUX
218 : REAL*8 FLUX(NESP+4)
219 : ENDSEGMENT
220 : POINTEUR IFLU1.IFLUX, IFLU2.IFLUX
221 : C
222 : C**** KRIPAD pour la correspondance global/local des conditions limits
223 : C
223 : CALL KRIPAD(MELLIM,MLELIM)
224 : c SEGACT MELLIM
225 : C
226 : C**** Initialisation des MCHAMLs
227 : C
228 : C**** Masse volumique
229 : C
230 : MCHEL1 = IROF
231 : SEGACT MCHEL1
232 : MCHAM1 = MCHEL1.ICHAML(1)
233 : SEGACT MCHAM1
234 : MELRO = MCHAM1.IELVAL(1)
235 : SEGDES MCHEL1
236 : SEGDES MCHAM1
237 : C
238 : C**** Pression
239 : C
240 : MCHEL1 = IPF
241 : SEGACT MCHEL1
242 : MCHAM1 = MCHEL1.ICHAML(1)
243 : SEGACT MCHAM1
244 : MELP = MCHAM1.IELVAL(1)
245 : SEGDES MCHEL1
246 : SEGDES MCHAM1
247 : C
248 : C**** Gamma
249 : C
250 : MCHEL1 = IGAMF
251 : SEGACT MCHEL1
252 : MCHAM1 = MCHEL1.ICHAML(1)
253 : SEGACT MCHAM1
254 : MELGAM = MCHAM1.IELVAL(1)
255 : SEGDES MCHEL1
256 : SEGDES MCHAM1
257 : C
258 : C**** Vitesse et cosinus directeurs du repere (n,t)
259 : C
260 : MCHEL1 = IVITF
261 : SEGACT MCHEL1
262 : C
263 : C**** La vitesse a comme SPG MELEFE
264 : C Le cosinus directeurs ont comme SPG MELEMF
265 : C
266 : C MCHAM1 -> Cosinus directeurs
267 : C MCHAM2 -> Vitesse
268 : C
269 : SPG1 = MCHEL1.IMACHE(1)
270 : SPG2 = MCHEL1.IMACHE(2)
271 : IF((SPG1 .EQ. MELEMF) .AND. (SPG2 .EQ. MELEFE))THEN
272 : MCHAM1 = MCHEL1.ICHAML(1)
273 : MCHAM2 = MCHEL1.ICHAML(2)
274 : ELSEIF((SPG1 .EQ. MELEFE) .AND. (SPG2 .EQ. MELEMF))THEN
275 : MCHAM1 = MCHEL1.ICHAML(2)
276 : MCHAM2 = MCHEL1.ICHAML(1)
277 : ELSE
278 : LOGAN = .TRUE.
279 : GOTO 9999
280 : ENDIF
281 : SEGACT MCHAM1
282 : MELVNX = MCHAM1.IELVAL(1)
283 : MELVNY = MCHAM1.IELVAL(2)
284 : MELT1X = MCHAM1.IELVAL(3)
285 : MELT1Y = MCHAM1.IELVAL(4)
286 : SEGDES MCHAM1
287 : SEGACT MCHAM2
288 : MELVUN = MCHAM2.IELVAL(1)
289 : MELVUT = MCHAM2.IELVAL(2)
290 : SEGDES MCHAM2
291 : SEGDES MCHEL1
292 : C
293 : C**** Fractions massiques
294 : C
295 : IF(LOGME)THEN
296 : MCHEL1 = IFRMAF
297 : SEGACT MCHEL1
298 : MCHAMY = MCHEL1.ICHAML(1)
299 : SEGACT MCHAMY
300 : C
301 : C******* Numero d'especes dans les equations d'Euler
302 : C
303 : NESP = MCHAMY.IELVAL(/1)
304 : DO I1 = 1, NESP
305 : MELVA1 = MCHAMY.IELVAL(I1)
306 : SEGACT MELVA1
307 : ENDDO
308 : SEGINI FRAMAG
309 : SEGINI FRAMAD
310 : SEGDES MCHEL1
311 : ELSE
312 : C
313 : C******* Definition minimale de YET, necessaire pour transmetre YET aux
314 : C subroutines FORTRAN qui calculent les flux
315 : C
316 : NESP = 1
317 : SEGINI FRAMAG
318 : SEGINI FRAMAD
319 : NESP = 0
320 : ENDIF
321 : C
322 : C**** Initialisation des MELEMEs
323 : C
324 : C 'CENTRE', 'FACEL'
325 : C
326 : IPT2 = MELEFE
327 : SEGACT IPT2
328 : NFAC = IPT2.NUM(/2)
329 : C
330 : C**** KRIPAD pour la correspondance global/local de centre
331 : C
331 : CALL KRIPAD(MELEMC,MLENT1)
332 : C
333 : C**** MLENTI1 a MCORD.XCOORD(/1)/(IDIM+1) elements
334 : C
335 : C Si i est le numero global d'un noeud de ICEN,
336 : C MLENT1.LECT(i) contient sa position, i.e.
337 : C
338 : C I = numero global du noeud centre
339 : C MLENT1.LECT(i) = numero local du noeud centre
340 : C
341 : C MLENT1 déjà activé, i.e.
342 : C
343 : C SEGACT MLENT1
344 : C
345 : C
346 : C**** KRIPAD pour la correspondance global/local de 'FACE'
347 : C
347 : CALL KRIPAD(MELEMF,MLENT2)
348 : C SEGACT MELNT2
349 : C
350 : C**** Initialisation de flux
351 : C
352 : SEGINI IFLU1
353 : SEGINI IFLU2
354 : C
355 : C**** IFLU2 = segment de travail en FLUVLH; c'est plus rapide le definir ici
356 : C
357 : C
358 : C**** CHPOINTs de la table DOMAINE
359 : C
359 : CALL LICHT(ICHPSU,MPOVSU,TYPE,IGEOMF)
359 : CALL LICHT(ICHPDI,MPOVDI,TYPE,IGEOMC)
360 : C
361 : C**** LICHT active les MPOVALs en *MOD
362 : C
363 : C i.e.
364 : C
365 : C SEGACT MPOVSU*MOD
366 : C SEGACT MPOVDI*MOD
367 : C
368 : IF(IGEOMF .NE. MELEMF)THEN
369 : MESERR = 'Il ne faut pas jouer avec la console. '
370 : LOGAN = .TRUE.
371 : GOTO 9999
372 : ENDIF
373 : IF(IGEOMC .NE. MELEMC)THEN
374 : MESERR = 'Il ne faut pas jouer avec la console. '
375 : LOGAN = .TRUE.
376 : GOTO 9999
377 : ENDIF
378 : C
379 : C**** Les FLUX aux face
380 : C
381 : C La densité
382 : C
382 : CALL LICHT(ICHFLU,MPOFLU,TYPE,IGEOMF)
383 : C
384 : C SEGACT MPOFLU*MOD
385 : C
386 : IF(IGEOMC .NE. MELEMC)THEN
387 : MESERR = 'Il ne faut pas jouer avec la console. '
388 : LOGAN = .TRUE.
389 : GOTO 9999
390 : ENDIF
391 : C
392 : C**** Activation des MCHAMLs
393 : C
394 : SEGACT MELRO
395 : SEGACT MELP
396 : SEGACT MELGAM
397 : SEGACT MELVUN
398 : SEGACT MELVUT
399 : SEGACT MELVNX
400 : SEGACT MELVNY
401 : SEGACT MELT1X
402 : SEGACT MELT1Y
403 : C
404 : C**** Bas Mach
405 : C
406 : IF(IUINF .NE. 0)THEN
406 : CALL LICHT(IUPRI,MPUPRI,TYPE,IGEOMC)
406 : CALL LICHT(IUINF,MPUINF,TYPE,IGEOMC)
407 : C SEGACT MPUPRI*MOD
408 : C SEGACT MPUINF*MOD
409 : ELSE
410 : MPUPRI=0
411 : MPUINF=0
412 : ENDIF
413 : C
414 : C**** Initialisation de 1/DT
415 : C
416 : UNSDT = 0.0D0
417 : C
418 : C**** BOUCLE SUR FACEL pour le calcul du FLUX
419 : C
420 : DO NLCF = 1, NFAC
421 : C
422 : C******* NLCF = numero local du centre de facel
423 : C NGCF = numero global du centre de facel
424 : C NLCF1 = numero local du centre de face
425 : C NGCEG = numero global du centre ELT "gauche"
426 : C NLCEG = numero local du centre ELT "gauche"
427 : C NGCED = numero global du centre ELT "droite"
428 : C NLCED = numero local du centre ELT "droite"
429 : C
430 : NGCEG = IPT2.NUM(1,NLCF)
431 : NGCED = IPT2.NUM(3,NLCF)
432 : NGCF = IPT2.NUM(2,NLCF)
433 : NLCF1 = MLENT2.LECT(NGCF)
434 : NLCEG = MLENT1.LECT(NGCEG)
435 : NLCED = MLENT1.LECT(NGCED)
436 : NLFL = MLELIM.LECT(NGCF)
437 : C
438 : C******* NLCF != NLCF1 -> l'auteur (MOI) n'a rien compris.
439 : C
440 : IF(NLCF .NE. NLCF1)THEN
441 : MESERR = 'Il ne faut pas jouer avec la console. '
442 : LOGAN = .TRUE.
443 : GOTO 9999
444 : ENDIF
445 : IF(NLFL .EQ. 0)THEN
446 : C
447 : C******* Recuperation des Etats "gauche" et "droite"
448 : C
449 : ROG = MELRO.VELCHE(1,NLCF)
450 : UNG = MELVUN.VELCHE(1,NLCF)
451 : UTG = MELVUT.VELCHE(1,NLCF)
452 : PG = MELP.VELCHE(1,NLCF)
453 : GAMG = MELGAM.VELCHE(1,NLCF)
454 : C
455 : ROD = MELRO.VELCHE(3,NLCF)
456 : UND = MELVUN.VELCHE(3,NLCF)
457 : UTD = MELVUT.VELCHE(3,NLCF)
458 : PD = MELP.VELCHE(3,NLCF)
459 : GAMD = MELGAM.VELCHE(3,NLCF)
460 : C
461 : CNX = MELVNX.VELCHE(1,NLCF)
462 : CNY = MELVNY.VELCHE(1,NLCF)
463 : CTX = MELT1X.VELCHE(1,NLCF)
464 : CTY = MELT1Y.VELCHE(1,NLCF)
465 : C
466 : C******* Le fractiones massiques
467 : C
468 : IF(LOGME)THEN
469 : DO I1 = 1, NESP
470 : MELVA1 = MCHAMY.IELVAL(I1)
471 : FRAMAG.YET(I1) = MELVA1.VELCHE(1,NLCF)
472 : FRAMAD.YET(I1) = MELVA1.VELCHE(3,NLCF)
473 : ENDDO
474 : ENDIF
475 : C
476 : C******* On a defini (ROg,ROUNg,ROUTg,Pg,(Yg)), (ROd,ROUNd,ROUTd,Pd,(Yd))
477 : C et on a déjà verifié ROg, ROd, Pg, Pd > 0 et 0FLURIE(NESP,
488 : & GAMG,ROG,PG,UNG,UTG,
489 : & GAMD,ROD,PD,UND,UTD,
490 : & FRAMAG.YET,FRAMAD.YET,
491 : & IFLU1.FLUX,
492 : & CELLT,
493 : & LOGNC,LOGAN,MESERR)
494 : C
495 : ELSEIF(INDMET .EQ. 2)THEN
496 : C
497 : C******* Van Leer FVS
498 : C
499 : C N.B: FLUXVL en FORTRAN pure
500 : C FRAMAG.YET = table d'un pointeur -> table
501 : C La meme chose pour FRAMAD.YET, IFLU1.FLUX,
502 : C IFLU2.FLUX
503 : C
503 : CALL FLUXVL(NESP,
504 : & GAMG,ROG,PG,UNG,UTG,
505 : & GAMD,ROD,PD,UND,UTD,
506 : & FRAMAG.YET,FRAMAD.YET,
507 : & IFLU1.FLUX,IFLU2.FLUX,
508 : & CELLT)
509 : ELSEIF(INDMET .EQ. 3)THEN
510 : C
511 : C******* Van Leer - Hanel FVS
512 : C
513 : C N.B: FLUVLH en FORTRAN pure
514 : C FRAMAG.YET = table d'un pointeur -> table
515 : C La meme chose pour FRAMAD.YET, IFLU1.FLUX,
516 : C IFLU2.FLUX
517 : C
517 : CALL FLUVLH(NESP,
518 : & GAMG,ROG,PG,UNG,UTG,
519 : & GAMD,ROD,PD,UND,UTD,
520 : & FRAMAG.YET,FRAMAD.YET,
521 : & IFLU1.FLUX,IFLU2.FLUX,
522 : & CELLT)
523 : ELSEIF(INDMET .EQ. 4)THEN
524 : C
525 : C******* HUS (Van Leer FVS + Osher FDS)
526 : C
526 : CALL FLHUS1(NESP,
527 : & GAMG,ROG,PG,UNG,UTG,
528 : & GAMD,ROD,PD,UND,UTD,
529 : & FRAMAG.YET,FRAMAD.YET,
530 : & IFLU1.FLUX,IFLU2.FLUX,
531 : & CELLT,
532 : & LOGNC,MESERR,LOGAN)
533 : ELSEIF(INDMET .EQ. 5)THEN
534 : C
535 : C******* HUS (Van Leer - Hanel FVS + Osher FDS)
536 : C
536 : CALL FLHUS2(NESP,
537 : & GAMG,ROG,PG,UNG,UTG,
538 : & GAMD,ROD,PD,UND,UTD,
539 : & FRAMAG.YET,FRAMAD.YET,
540 : & IFLU1.FLUX,IFLU2.FLUX,
541 : & CELLT,
542 : & LOGNC,MESERR,LOGAN)
543 : ELSEIF(INDMET .EQ. 6)THEN
544 : C
545 : C******* AUSM+ (Liou)
546 : C
546 : CALL FAUSMP(NESP,
547 : & GAMG,ROG,PG,UNG,UTG,
548 : & GAMD,ROD,PD,UND,UTD,
549 : & FRAMAG.YET,FRAMAD.YET,
550 : & IFLU1.FLUX,
551 : & CELLT)
552 : ELSEIF(INDMET .EQ. 7)THEN
553 : C
554 : C******* Roe Scheme
555 : C
555 : CALL FROEMO(NESP,
556 : & GAMG,ROG,PG,UNG,UTG,
557 : & GAMD,ROD,PD,UND,UTD,
558 : & FRAMAG.YET,FRAMAD.YET,
559 : & IFLU1.FLUX,
560 : & CELLT)
561 : ELSEIF(INDMET .EQ. 8)THEN
562 : C
563 : C******* Choc-Choc
564 : C
564 : CALL FCOLCP(NESP,
565 : & GAMG,ROG,PG,UNG,UTG,
566 : & GAMD,ROD,PD,UND,UTD,
567 : & FRAMAG.YET,FRAMAD.YET,
568 : & IFLU1.FLUX,
569 : & CELLT,
570 : & LOGNC,LOGAN,MESERR)
571 : ELSEIF(INDMET .EQ. 9)THEN
572 : C
573 : C******* AUSM+(P) for low Mach number flows (Edwards and Liou)
574 : C
575 : V_INF=MAX(MPUINF.VPOCHA(NLCEG,1),
576 : & MPUINF.VPOCHA(NLCED,1),
577 : & MPUPRI.VPOCHA(NLCEG,1),
578 : & MPUPRI.VPOCHA(NLCED,1))
578 : CALL FAUSM2(NESP,
579 : & GAMG,ROG,PG,UNG,UTG,
580 : & GAMD,ROD,PD,UND,UTD,
581 : & FRAMAG.YET,FRAMAD.YET,
582 : & V_INF,
583 : & IFLU1.FLUX,
584 : & CELLT)
585 : ELSEIF(INDMET .EQ. 10)THEN
586 : C
587 : C******* Rusanov scheme
588 : C
588 : CALL FRUSAN(NESP,
589 : & GAMG,ROG,PG,UNG,UTG,
590 : & GAMD,ROD,PD,UND,UTD,
591 : & FRAMAG.YET,FRAMAD.YET,
592 : & IFLU1.FLUX,
593 : & CELLT)
594 : ELSEIF(INDMET .EQ. 11)THEN
595 : C
596 : C******* Rusanov scheme for Low-Mach number flows
597 : C
598 : V_INF=MAX(MPUINF.VPOCHA(NLCEG,1),
599 : & MPUINF.VPOCHA(NLCED,1),
600 : & MPUPRI.VPOCHA(NLCEG,1),
601 : & MPUPRI.VPOCHA(NLCED,1))
601 : CALL FRUSBM(NESP,
602 : & GAMG,ROG,PG,UNG,UTG,
603 : & GAMD,ROD,PD,UND,UTD,
604 : & FRAMAG.YET,FRAMAD.YET,
605 : & V_INF,
606 : & IFLU1.FLUX,
607 : & CELLT)
608 : ELSEIF(INDMET .EQ. 12)THEN
609 : C
610 : C******* Centered scheme
611 : C
611 : CALL FVCENT(NESP,
612 : & GAMG,ROG,PG,UNG,UTG,
613 : & GAMD,ROD,PD,UND,UTD,
614 : & FRAMAG.YET,FRAMAD.YET,
615 : & IFLU1.FLUX,
616 : & CELLT)
617 : ELSEIF(INDMET .EQ. 13)THEN
618 : C
619 : C******* Roe-Turkel scheme for Low-Mach number flows
620 : C
621 : V_INF=MAX(MPUINF.VPOCHA(NLCEG,1),
622 : & MPUINF.VPOCHA(NLCED,1),
623 : & MPUPRI.VPOCHA(NLCEG,1),
624 : & MPUPRI.VPOCHA(NLCED,1))
624 : CALL FROEBM(NESP,
625 : & GAMG,ROG,PG,UNG,UTG,
626 : & GAMD,ROD,PD,UND,UTD,
627 : & FRAMAG.YET,FRAMAD.YET,
628 : & V_INF,
629 : & IFLU1.FLUX,
630 : & CELLT)
631 : ELSEIF(INDMET .EQ. 14)THEN
632 : C
633 : C******* HLLC
634 : C
634 : CALL FLHLLC(NESP,
635 : & GAMG,ROG,PG,UNG,UTG,
636 : & GAMD,ROD,PD,UND,UTD,
637 : & FRAMAG.YET,FRAMAD.YET,
638 : & IFLU1.FLUX,
639 : & CELLT)
640 : ELSEIF(INDMET .EQ. 15)THEN
641 : C
642 : C******* HLLC-Turkel scheme for Low-Mach number flows
643 : C
644 : V_INF=MAX(MPUINF.VPOCHA(NLCEG,1),
645 : & MPUINF.VPOCHA(NLCED,1),
646 : & MPUPRI.VPOCHA(NLCEG,1),
647 : & MPUPRI.VPOCHA(NLCED,1))
647 : CALL FHLLCL(NESP,
648 : & GAMG,ROG,PG,UNG,UTG,
649 : & GAMD,ROD,PD,UND,UTD,
650 : & FRAMAG.YET,FRAMAD.YET,
651 : & V_INF,
652 : & IFLU1.FLUX,
653 : & CELLT)
654 : ELSEIF(INDMET .EQ. 16)THEN
655 : C
656 : C******* AUSM+up for Low-Mach number flows
657 : C
658 : V_INF=MAX(MPUINF.VPOCHA(NLCEG,1),
659 : & MPUINF.VPOCHA(NLCED,1),
660 : & MPUPRI.VPOCHA(NLCEG,1),
661 : & MPUPRI.VPOCHA(NLCED,1))
661 : CALL FAU_UP(NESP,
662 : & GAMG,ROG,PG,UNG,UTG,
663 : & GAMD,ROD,PD,UND,UTD,
664 : & FRAMAG.YET,FRAMAD.YET,
665 : & V_INF,
666 : & IFLU1.FLUX,
667 : & CELLT)
668 : ENDIF
669 : C
670 : IF(LOGAN) GOTO 9999
671 : IF(LOGNC) GOTO 9999
672 : C
673 : C******* Ecriture des flux
674 : C
675 : C FLUX(1) = RO Un RO Un
676 : C FLUX(2) = RO Un Un + P -> RO Un Ux + P CNX
677 : C FLUX(3) = RO Un Ut -> RO Un Uy + P CNY
678 : C FLUX(4) = RO Un Ht RO Un Ht
679 : C
680 : SURF = MPOVSU.VPOCHA(NLCF,1)
681 : MPOFLU.VPOCHA(NLCF,1) =
682 : & (IFLU1.FLUX(1) * SURF )
683 : MPOFLU.VPOCHA(NLCF,2) =
684 : & ((IFLU1.FLUX(2)*CNX+IFLU1.FLUX(3)*CTX) * SURF)
685 : MPOFLU.VPOCHA(NLCF,3) =
686 : & ((IFLU1.FLUX(2)*CNY+IFLU1.FLUX(3)*CTY) * SURF)
687 : MPOFLU.VPOCHA(NLCF,4) =
688 : & (IFLU1.FLUX(4) * SURF)
689 : IF(LOGME)THEN
690 : DO I1 = 1, NESP
691 : MPOFLU.VPOCHA(NLCF,4+I1)= IFLU1.FLUX(4+I1)* SURF
692 : ENDDO
693 : ENDIF
694 : C
695 : C******* Calcul du pas du temps (CFL)
696 : C
697 : C
698 : C
699 : C****** a) etat a l'interface
700 : C
701 : C CELLT a la dimension du reciproque d'une vitesse
702 : C
703 : DIAMG = MPOVDI.VPOCHA(NLCEG,1)
704 : DIAMD = MPOVDI.VPOCHA(NLCED,1)
705 : DIAM = MIN(DIAMG,DIAMD)
706 : CELL = 1.0D0/(DIAM*CELLT)
707 : IF(CELL .GT. UNSDT)THEN
708 : UNSDT = CELL
709 : ENDIF
710 : C
711 : C****** b) etat gauche
712 : C
713 : ASON = SQRT(GAMG*PG/ROG)
714 : LAMBDA = ABS(UNG) + ASON
715 : CELL = LAMBDA / DIAM
716 : IF(CELL .GT. UNSDT)THEN
717 : UNSDT = CELL
718 : ENDIF
719 : C
720 : C****** C) etat droite
721 : C
722 : ASON = SQRT(GAMD*PD/ROD)
723 : LAMBDA = ABS(UND) + ASON
724 : CELL = LAMBDA / DIAM
725 : IF(CELL .GT. UNSDT)THEN
726 : UNSDT = CELL
727 : ENDIF
728 : ENDIF
729 : C
730 : C**** Fin boucle sur FACEL
731 : C
732 : ENDDO
733 : C
734 : C**** Pas du temps (condition de non interaction en 1D)
735 : C
736 : DT = 0.5D0 / UNSDT
737 : C
738 : C**** Desactivation des segments et
739 : C on detruit les MCHAMLs
740 : C
741 : C
742 : C**** SEGSUP FRAMAG
743 : C SEGSUP FRAMAD
744 : C
745 : C meme si LOGME = .FALSE.
746 : C
747 : SEGSUP FRAMAG
748 : SEGSUP FRAMAD
749 : C
750 : SEGSUP MLENT1
751 : SEGSUP MLENT2
752 : SEGDES IPT2
753 : C
754 : SEGSUP IFLU1
755 : SEGSUP IFLU2
756 : C
757 : SEGDES MPOVSU
758 : SEGDES MPOVDI
759 : C
760 : SEGDES MPOFLU
761 : C
762 : SEGDES MELRO
763 : SEGDES MELP
764 : SEGDES MELGAM
765 : SEGDES MELVUN
766 : SEGDES MELVUT
767 : SEGDES MELVNX
768 : SEGDES MELVNY
769 : SEGDES MELT1X
770 : SEGDES MELT1Y
771 : C
772 : IF(MPUPRI .NE. 0)THEN
773 : SEGDES MPUPRI
774 : SEGDES MPUINF
775 : ENDIF
776 : C
777 : IF(LOGME) THEN
778 : DO I1 = 1, NESP
779 : MELVA1 = MCHAMY.IELVAL(I1)
780 : SEGDES MELVA1
781 : ENDDO
782 : C
783 : SEGDES MCHAMY
784 : ENDIF
785 : SEGSUP MLELIM
786 : CC
787 : 9999 CONTINUE
788 : C
789 : RETURN
790 : END
791 : C
792 :
793 :
794 :
795 :
796 :
797 :
798 :
799 :
800 :
801 :
802 :
803 :
804 :
805 :
© Cast3M 2003 - All rights reserved.
Disclaimer