1 : C KONFL4 SOURCE PV 09/03/12 21:26:31 6325
2 : SUBROUTINE KONFL4(LOGME,LOGSCA,INDMET,NORDP1,
3 : & IROF,IVITF,IPF,IFRMAF,ISCAF,PROPHY,
4 : & ICHPSU,ICHPDI,
5 : & MELEMC,MELEMF,MELEFE,MELLIM,
6 : & ICHFLU,DT,
7 : & LOGNC,LOGAN,MESERR)
8 : C************************************************************************
9 : C
10 : C PROJET : CASTEM 2000
11 : C
12 : C NOM : KONFL4
13 : C
14 : C DESCRIPTION : Voir KONV12
15 : C
16 : C Cas trois dimensions, gaz "thermally perfect"
17 : C
18 : C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
19 : C
20 : C AUTEUR : R. MOREL, A. BECCANTINI, DRN/DMT/SEMT/TTMF
21 : C
22 : C************************************************************************
23 : C
24 : C
25 : C APPELES (Outils
26 : C CASTEM) : KRIPAD, LICHT
27 : C
28 : C APPELES (Calcul) : FVLHT2, FCOLT2
29 : C
30 : C
31 : C************************************************************************
32 : C
33 : C ENTREES
34 : C
35 : C
36 : C 1) PARAMETRES
37 : C
38 : C LOGME : (LOGICAL); .TRUE. -> MULTI-ESPECES
39 : C .FALSE. -> MONO-ESPECE
40 : C
41 : C LOGSCA : (LOGICAL); .TRUE. -> SCALAIRES-PASSIVES
42 : C .FALSE. -> ou non
43 : C
44 : C INDMET : 1 VLH (van Leer Hanel FVS)
45 : C
46 : C 2 Colella-Glaz FDS
47 : C
48 : C NORDP1 : (ordre des polynoms cv(T)) + 1
49 : C
50 : C 2) Pointeurs des MCHAMLs
51 : C
52 : C IROF : MCHAML sur "FACEL" contenant la masse volumique
53 : C ("gauche" et "droite");
54 : C
55 : C IVITF : MCHAML sur "FACEL" contenant la vitesse dans le repaire
56 : C local (n,t) et les cosinus directeurs des repaire local;
57 : C
58 : C IPF : MCHAML sur "FACEL" contenant la pression;
59 : C
60 : C IFRAMAF : MCHAML sur "FACEL", contenant les fractions massiques
61 : C si LOGME = .TRUE.;
62 : C LOGME = .FALSE. -> IFRAMAF = 0
63 : C
64 : C ISCAF : MCHAML sur "FACEL", contenant les scalires passifs a
65 : C transporter (ou ISCAF = 0)
66 : C
67 : C 3) Pointeur sur la table des proprietes de gaz
68 : C
69 : C PROPHY
70 : C
71 : C 4) Pointeurs de CHPOINTs de la table DOMAINE
72 : C
73 : C ICHPSU : CHPOINT "FACE" contenant la surface des faces
74 : C
75 : C ICHPDI : CHPOINT "CENTRE" contenant le diametre minimum
76 : C de chaque element
77 : C
78 : C
79 : C 5) Pointeurs de MELEME de la table DOMAINE
80 : C
81 : C MELEMC : MELEME 'CENTRE' du SPG des CENTRES
82 : C
83 : C MELEMF : MELEME 'FACE' du SPG des FACES
84 : C
85 : C MELEFE : MELEME 'FACEL' du connectivité FACES -> ELEM
86 : C
87 : C MELLIM : MAILLAGE where the flux (or residu) will not be found
88 : C
89 : C SORTIES
90 : C
91 : C ICHFLU : pointeur du CHPOINT flux aux "FACE"s
92 : C
93 : C DT : pas de temps pour le respect de la CFL-like condition
94 : C DT < DIAMEL /2 /max(Lambda_i)
95 : C En maillage regulier cette condition garantie la
96 : C non-interaction des ondes
97 : C
98 : C
99 : C LOGNC : (LOGICAL): si .TRUE. la methode de Newton-Rapson, utilisée
100 : C dans pour la solution du probleme Riemann exacte ou dans
101 : C l'algorithm HUS, n'a pas bien marchéee; MESERR = 'Goudunov'
102 : C ou 'HUS'.
103 : C
104 : C LOGAN : (LOGICAL): si .TRUE. une anomalie à été detectée
105 : C
106 : C MESERR : pour l'ecriture des messages d'erreurs
107 : C
108 : C************************************************************************
109 : C
110 : C HISTORIQUE (Anomalies et modifications éventuelles)
111 : C
112 : C HISTORIQUE : 22.12.98 Creation
113 : C
114 : C 08.02.00 Ajout du flux numerique de Colella-Glaz
115 : C
116 : C 21.02.00 Ajout de transport de scalaires passifs
117 : C
118 : C 03.12.03 Ajout le MAILLIM comme le l'entree
119 : C
120 : C************************************************************************
121 : C
122 : C
123 : C N.B.: On suppose qu'on a déjà controllé RO, P > 0
124 : C GAMMA \in (1,3)
125 : C Y \in (0,1)
126 : C Si non il faut le faire!!!
127 : C
128 : C************************************************************************
129 : C
130 : C
131 : C**** Variables de COOPTIO
132 : C
133 : C INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
134 : C & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ
135 : C & ,IOPER, IOSGB, IOGRA, IOSAU, IORES
136 : C & ,IECHO, IIMPI, IOSPI
137 : C & ,IDIM
138 : CC & ,MCOORD
139 : C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
140 : C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
141 : C & ,NORINC,NORVAL,NORIND,NORVAD
142 : C & ,NUCROU, IPSAUV
143 : C
144 : IMPLICIT INTEGER(I-N)
145 : INTEGER I1, I2
146 : & ,INDMET, NORDP1
147 : & ,IROF,IVITF,IPF,IFRMAF, ISCAF
148 : & ,ICHPSU,ICHPDI,MELEMC,MELEMF,MELEFE
149 : & ,IGEOMC,IGEOMF
150 : & ,ICHFLU
151 : & ,NESP, NFAC, NSCA, NMA, NMA2
152 : & ,NLCF, NGCEG, NGCED, NLCEG, NLCED
153 : & ,NGCF, NLCF1, SPG1, SPG2, NLFL
154 : REAL*8 DT, UNSDT, CELLT
155 : & , ROG, UNG, UTG, UVG, PG, GAMG, TG
156 : & , ROD, UND, UTD, UVD, PD, GAMD, TD
157 : & , SURF,CNX, CNY, CNZ, CTX , CTY, CTZ
158 : & , CVX, CVY, CVZ
159 : & , CELL, DIAMG, DIAMD, DIAM
160 : & , ASON, LAMBDA
161 : & , Y1G, Y1D, RTOTG, RTOTD, CVTOTG, CVTOTD, ETHERG, ETHERD
162 : & , YG, YD, FUNTG, FUNTD, DCV, XST
163 : LOGICAL LOGME, LOGNC, LOGAN, LOGSCA
164 : CHARACTER*(40) MESERR
165 : CHARACTER*(8) TYPE
166 : C
167 : C**** Segment des proprietes du gaz
168 : C
169 : SEGMENT PROPHY
170 : REAL*8 ACV(NORDP1,NESP+1), R(NESP+1), H0K(NESP+1)
171 : & ,ACVTOG(NORDP1), ACVTOD(NORDP1)
172 : ENDSEGMENT
173 : C
174 : C**** Les Includes
175 : C
176 : C
177 : C**** LES INCLUDES
178 : C
179 : -INC CCOPTIO
180 : -INC SMCOORD
181 : -INC SMCHAML
182 : POINTEUR MELVNX.MELVAL, MELVNY.MELVAL,MELVNZ.MELVAL,
183 : & MELT1X.MELVAL, MELT1Y.MELVAL,MELT1Z.MELVAL,
184 : & MELT2X.MELVAL, MELT2Y.MELVAL,MELT2Z.MELVAL
185 : POINTEUR MELVUN.MELVAL, MELVUT.MELVAL, MELVUV.MELVAL
186 : POINTEUR MELRO.MELVAL, MELP.MELVAL
187 : -INC SMCHPOI
188 : POINTEUR MPOVSU.MPOVAL, MPOVDI.MPOVAL
189 : & , MPOFLU.MPOVAL
190 : POINTEUR MCHAMY.MCHAML, MCHAMS.MCHAML
191 : -INC SMELEME
192 : POINTEUR MELLIM.MELEME
193 : -INC SMLMOTS
194 : -INC SMLENTI
195 : POINTEUR MLELIM.MLENTI
196 : C
197 : C**** Les fractiones massiques.
198 : C
199 : SEGMENT FRAMAS
200 : REAL*8 YET(NMA)
201 : ENDSEGMENT
202 : POINTEUR FRAMAG.FRAMAS, FRAMAD.FRAMAS, SCAG.FRAMAS, SCAD.FRAMAS
203 : C
204 : C**** Les flux aux interface dans le repaire (n,t)
205 : C
206 : SEGMENT IFLUX
207 : REAL*8 FLUX(NMA2)
208 : ENDSEGMENT
209 : POINTEUR IFLU1.IFLUX, IFLU2.IFLUX
210 : C--------------------------------------------
210 : CALL KRIPAD(MELLIM,MLELIM)
211 : C--------------------------------------------
212 : C
213 : C**** Initialisation des MCHAMLs
214 : C
215 : C**** Masse volumique
216 : C
217 : MCHEL1 = IROF
218 : SEGACT MCHEL1
219 : MCHAM1 = MCHEL1.ICHAML(1)
220 : SEGACT MCHAM1
221 : MELRO = MCHAM1.IELVAL(1)
222 : SEGDES MCHEL1
223 : SEGDES MCHAM1
224 : C
225 : C**** Pression
226 : C
227 : MCHEL1 = IPF
228 : SEGACT MCHEL1
229 : MCHAM1 = MCHEL1.ICHAML(1)
230 : SEGACT MCHAM1
231 : MELP = MCHAM1.IELVAL(1)
232 : SEGDES MCHEL1
233 : SEGDES MCHAM1
234 : C
235 : C**** Vitesse et cosinus directeurs du repere (n,t)
236 : C
237 : MCHEL1 = IVITF
238 : SEGACT MCHEL1
239 : C
240 : C**** La vitesse a comme SPG MELEFE
241 : C Le cosinus directeurs ont comme SPG MELEMF
242 : C
243 : C MCHAM1 -> Cosinus directeurs
244 : C MCHAM2 -> Vitesse
245 : C
246 : SPG1 = MCHEL1.IMACHE(1)
247 : SPG2 = MCHEL1.IMACHE(2)
248 : IF((SPG1 .EQ. MELEMF) .AND. (SPG2 .EQ. MELEFE))THEN
249 : MCHAM1 = MCHEL1.ICHAML(1)
250 : MCHAM2 = MCHEL1.ICHAML(2)
251 : ELSEIF((SPG1 .EQ. MELEFE) .AND. (SPG2 .EQ. MELEMF))THEN
252 : MCHAM1 = MCHEL1.ICHAML(2)
253 : MCHAM2 = MCHEL1.ICHAML(1)
254 : ELSE
255 : LOGAN = .TRUE.
256 : GOTO 9999
257 : ENDIF
258 : SEGACT MCHAM1
259 : MELVNX = MCHAM1.IELVAL(1)
260 : MELVNY = MCHAM1.IELVAL(2)
261 : MELVNZ = MCHAM1.IELVAL(3)
262 : MELT1X = MCHAM1.IELVAL(4)
263 : MELT1Y = MCHAM1.IELVAL(5)
264 : MELT1Z = MCHAM1.IELVAL(6)
265 : MELT2X = MCHAM1.IELVAL(7)
266 : MELT2Y = MCHAM1.IELVAL(8)
267 : MELT2Z = MCHAM1.IELVAL(9)
268 : SEGDES MCHAM1
269 : SEGACT MCHAM2
270 : MELVUN = MCHAM2.IELVAL(1)
271 : MELVUT = MCHAM2.IELVAL(2)
272 : MELVUV = MCHAM2.IELVAL(3)
273 : SEGDES MCHAM2
274 : SEGDES MCHEL1
275 : C
276 : C**** Fractions massiques
277 : C
278 : IF(LOGME)THEN
279 : MCHEL1 = IFRMAF
280 : SEGACT MCHEL1
281 : MCHAMY = MCHEL1.ICHAML(1)
282 : SEGACT MCHAMY
283 : C
284 : C******* Numero d'especes dans les equations d'Euler
285 : C
286 : NESP = MCHAMY.IELVAL(/1)
287 : DO I1 = 1, NESP
288 : MELVA1 = MCHAMY.IELVAL(I1)
289 : SEGACT MELVA1
290 : ENDDO
291 : NMA = NESP
292 : SEGINI FRAMAG
293 : SEGINI FRAMAD
294 : SEGDES MCHEL1
295 : ELSE
296 : C
297 : C******* Definition minimale de YET, necessaire pour transmetre YET aux
298 : C subroutines FORTRAN qui calculent les flux
299 : C
300 : NMA = 1
301 : SEGINI FRAMAG
302 : SEGINI FRAMAD
303 : NESP = 0
304 : ENDIF
305 : C
306 : C**** Scalaires passifs
307 : C
308 : IF(LOGSCA)THEN
309 : MCHEL1 = ISCAF
310 : SEGACT MCHEL1
311 : MCHAMS = MCHEL1.ICHAML(1)
312 : SEGACT MCHAMS
313 : NSCA = MCHAMS.IELVAL(/1)
314 : DO I1 = 1, NSCA, 1
315 : MELVA1 = MCHAMS.IELVAL(I1)
316 : SEGACT MELVA1
317 : ENDDO
318 : NMA = NSCA
319 : SEGINI SCAG
320 : SEGINI SCAD
321 : SEGDES MCHEL1
322 : ELSE
323 : C
324 : C******* Definition minimale de YET, necessaire pour transmetre YET aux
325 : C subroutines FORTRAN qui calculent les flux
326 : C
327 : NMA = 1
328 : SEGINI SCAG
329 : SEGINI SCAD
330 : NSCA= 0
331 : ENDIF
332 : C
333 : C
334 : C**** Initialisation des MELEMEs
335 : C
336 : C 'CENTRE', 'FACEL'
337 : C
338 : IPT2 = MELEFE
339 : SEGACT IPT2
340 : NFAC = IPT2.NUM(/2)
341 : C
342 : C**** KRIPAD pour la correspondance global/local de centre
343 : C
343 : CALL KRIPAD(MELEMC,MLENT1)
344 : C
345 : C**** MLENTI1 a MCORD.XCOORD(/1)/(IDIM+1) elements
346 : C
347 : C Si i est le numero global d'un noeud de ICEN,
348 : C MLENT1.LECT(i) contient sa position, i.e.
349 : C
350 : C I = numero global du noeud centre
351 : C MLENT1.LECT(i) = numero local du noeud centre
352 : C
353 : C MLENT1 déjà activé, i.e.
354 : C
355 : C SEGACT MLENT1
356 : C
357 : C
358 : C**** KRIPAD pour la correspondance global/local de 'FACE'
359 : C
359 : CALL KRIPAD(MELEMF,MLENT2)
360 : C
361 : C**** Initialisation de flux
362 : C
363 : NMA2 = 5 + NESP + NSCA
364 : SEGINI IFLU1
365 : SEGINI IFLU2
366 : C
367 : C**** IFLU2 = segment de travail en FLUVLH; c'est plus rapide le definir ici
368 : C
369 : C
370 : C**** CHPOINTs de la table DOMAINE
371 : C
371 : CALL LICHT(ICHPSU,MPOVSU,TYPE,IGEOMF)
371 : CALL LICHT(ICHPDI,MPOVDI,TYPE,IGEOMC)
372 : C
373 : C**** LICHT active les MPOVALs en *MOD
374 : C
375 : C i.e.
376 : C
377 : C SEGACT MPOVSU*MOD
378 : C SEGACT MPOVDI*MOD
379 : C
380 : C
381 : C**** Les FLUX aux face
382 : C
382 : CALL LICHT(ICHFLU,MPOFLU,TYPE,IGEOMF)
383 : C
384 : C SEGACT MPOFLU*MOD
385 : C
386 : C**** Activation des MCHAMLs
387 : C
388 : SEGACT MELRO
389 : SEGACT MELP
390 : SEGACT MELVUN
391 : SEGACT MELVUT
392 : SEGACT MELVUV
393 : SEGACT MELVNX
394 : SEGACT MELVNY
395 : SEGACT MELVNZ
396 : SEGACT MELT1X
397 : SEGACT MELT1Y
398 : SEGACT MELT1Z
399 : SEGACT MELT2X
400 : SEGACT MELT2Y
401 : SEGACT MELT2Z
402 : C
403 : C**** Initialisation de 1/DT
404 : C
405 : UNSDT = 0.0D0
406 : C
407 : C**** BOUCLE SUR FACEL pour le calcul du FLUX
408 : C
409 : DO NLCF = 1, NFAC
410 : C
411 : C******* NLCF = numero local du centre de facel
412 : C NGCF = numero global du centre de facel
413 : C NLCF1 = numero local du centre de face
414 : C NGCEG = numero global du centre ELT "gauche"
415 : C NLCEG = numero local du centre ELT "gauche"
416 : C NGCED = numero global du centre ELT "droite"
417 : C NLCED = numero local du centre ELT "droite"
418 : C
419 : NGCEG = IPT2.NUM(1,NLCF)
420 : NGCED = IPT2.NUM(3,NLCF)
421 : NGCF = IPT2.NUM(2,NLCF)
422 : NLCF1 = MLENT2.LECT(NGCF)
423 : NLCEG = MLENT1.LECT(NGCEG)
424 : NLCED = MLENT1.LECT(NGCED)
425 : NLFL = MLELIM.LECT(NGCF)
426 : C
427 : C******* NLCF != NLCF1 -> l'auteur (MOI) n'a rien compris.
428 : C
429 : IF(NLCF .NE. NLCF1)THEN
430 : MESERR = 'Il ne faut pas jouer avec la console. '
431 : LOGAN = .TRUE.
432 : GOTO 9999
433 : ENDIF
434 : IF(NLFL .EQ. 0)THEN
435 : C
436 : C******* Recuperation des Etats "gauche" et "droite"
437 : C
438 : ROG = MELRO.VELCHE(1,NLCF)
439 : UNG = MELVUN.VELCHE(1,NLCF)
440 : UTG = MELVUT.VELCHE(1,NLCF)
441 : UVG = MELVUV.VELCHE(1,NLCF)
442 : PG = MELP.VELCHE(1,NLCF)
443 : C
444 : ROD = MELRO.VELCHE(3,NLCF)
445 : UND = MELVUN.VELCHE(3,NLCF)
446 : UTD = MELVUT.VELCHE(3,NLCF)
447 : UVD = MELVUV.VELCHE(3,NLCF)
448 : PD = MELP.VELCHE(3,NLCF)
449 : C
450 : CNX = MELVNX.VELCHE(1,NLCF)
451 : CNY = MELVNY.VELCHE(1,NLCF)
452 : CNZ = MELVNZ.VELCHE(1,NLCF)
453 : CTX = MELT1X.VELCHE(1,NLCF)
454 : CTY = MELT1Y.VELCHE(1,NLCF)
455 : CTZ = MELT1Z.VELCHE(1,NLCF)
456 : CVX = MELT2X.VELCHE(1,NLCF)
457 : CVY = MELT2Y.VELCHE(1,NLCF)
458 : CVZ = MELT2Z.VELCHE(1,NLCF)
459 : C
460 : C******* Le fractiones massiques, R et ACVTO
461 : C
462 : RTOTG = 0.0D0
463 : RTOTD = 0.0D0
464 : CVTOTG = 0.0D0
465 : CVTOTD = 0.0D0
466 : Y1G = 1.0D0
467 : Y1D = 1.0D0
468 : ETHERG = 0.0D0
469 : ETHERD = 0.0D0
470 : DO I1= 1, NORDP1, 1
471 : PROPHY.ACVTOG(I1) =0.0D0
472 : PROPHY.ACVTOD(I1) =0.0D0
473 : ENDDO
474 : DO I1 = 1, NESP, 1
475 : MELVA1 = MCHAMY.IELVAL(I1)
476 : YG = MELVA1.VELCHE(1,NLCF)
477 : YD = MELVA1.VELCHE(3,NLCF)
478 : Y1G = Y1G - YG
479 : Y1D = Y1D - YD
480 : FRAMAG.YET(I1) = YG
481 : FRAMAD.YET(I1) = YD
482 : RTOTG = RTOTG + YG * PROPHY.R(I1)
483 : RTOTD = RTOTD + YD * PROPHY.R(I1)
484 : DO I2 = 1, NORDP1
485 : PROPHY.ACVTOG(I2) = PROPHY.ACVTOG(I2) +
486 : & YG * PROPHY.ACV(I2,I1)
487 : PROPHY.ACVTOD(I2) = PROPHY.ACVTOD(I2) +
488 : & YD * PROPHY.ACV(I2,I1)
489 : ENDDO
490 : ENDDO
491 : RTOTG = RTOTG + Y1G * PROPHY.R(NESP+1)
492 : RTOTD = RTOTD + Y1D * PROPHY.R(NESP+1)
493 : DO I2 = 1, NORDP1, 1
494 : PROPHY.ACVTOG(I2) = PROPHY.ACVTOG(I2) +
495 : & Y1G * PROPHY.ACV(I2,NESP+1)
496 : PROPHY.ACVTOD(I2) = PROPHY.ACVTOD(I2) +
497 : & Y1D * PROPHY.ACV(I2,NESP+1)
498 : ENDDO
499 : TG = PG / (ROG * RTOTG)
500 : TD = PD / (ROD * RTOTD)
501 : FUNTG = 1.0D0
502 : FUNTD = 1.0D0
503 : CVTOTG = PROPHY.ACVTOG(1)
504 : ETHERG = CVTOTG * TG
505 : CVTOTD = PROPHY.ACVTOD(1)
506 : ETHERD = CVTOTD * TD
507 : DO I1 = 2, NORDP1, 1
508 : FUNTG = FUNTG * TG
509 : FUNTD = FUNTD * TD
510 : DCV = PROPHY.ACVTOG(I1) * FUNTG
511 : CVTOTG = CVTOTG + DCV
512 : ETHERG = ETHERG + DCV * TG / I1
513 : DCV = PROPHY.ACVTOD(I1) * FUNTD
514 : CVTOTD = CVTOTD + DCV
515 : ETHERD = ETHERD + DCV * TD / I1
516 : ENDDO
517 : GAMG = (CVTOTG + RTOTG) /CVTOTG
518 : GAMD = (CVTOTD + RTOTD) /CVTOTD
519 : C
520 : C******* Les scalaires passifs
521 : C
522 : DO I1 = 1, NSCA, 1
523 : MELVA1 = MCHAMS.IELVAL(I1)
524 : SCAG.YET(I1) = MELVA1.VELCHE(1,NLCF)
525 : SCAD.YET(I1) = MELVA1.VELCHE(3,NLCF)
526 : ENDDO
527 : C
528 : C******* On a defini (ROg,ROUNg,ROUTg,Pg,(Yg)), (ROd,ROUNd,ROUTd,Pd,(Yd))
529 : C et on a déjà verifié ROg, ROd, Pg, Pd > 0 et 0FVLHT2(NESP,NSCA,
538 : & GAMG,ROG,PG,UNG,UTG,UVG,ETHERG,
539 : & GAMD,ROD,PD,UND,UTD,UVD,ETHERD,
540 : & FRAMAG.YET,FRAMAD.YET,
541 : & SCAG.YET,SCAD.YET,
542 : & IFLU1.FLUX,IFLU2.FLUX,
543 : & CELLT)
544 : ELSEIF(INDMET .EQ. 2)THEN
545 : C
546 : C******* Colella-Glaz FDS (avec Entropy-fix)
547 : C
548 : XST=0.0D0
548 : CALL FCOLT2(NESP,NSCA,NORDP1,XST,
549 : & GAMG,RTOTG,PROPHY.ACVTOG,ROG,PG,TG,UNG,UTG,UVG,ETHERG,
550 : & GAMD,RTOTD,PROPHY.ACVTOD,ROD,PD,TD,UND,UTD,UVD,ETHERD,
551 : & FRAMAG.YET,FRAMAD.YET,
552 : & SCAG.YET,SCAD.YET,
553 : & IFLU1.FLUX,CELLT,
554 : & LOGNC,MESERR,LOGAN)
555 : ENDIF
556 : C
557 : IF(LOGAN) GOTO 9999
558 : IF(LOGNC) GOTO 9999
559 : C
560 : C******* Ecriture des flux
561 : C
562 : C FLUX(1) = RO Un RO Un
563 : C FLUX(2) = RO Un Un + P -> RO Un Ux + P CNX
564 : C FLUX(3) = RO Un Ut -> RO Un Uy + P CNY
565 : C FLUX(4) = RO Un Et RO Un Et
566 : C
567 : SURF = MPOVSU.VPOCHA(NLCF,1)
568 : MPOFLU.VPOCHA(NLCF,1) =
569 : & (IFLU1.FLUX(1) * SURF )
570 : MPOFLU.VPOCHA(NLCF,2) =
571 : &((IFLU1.FLUX(2)*CNX+IFLU1.FLUX(3)*CTX+IFLU1.FLUX(4)*CVX) * SURF)
572 : MPOFLU.VPOCHA(NLCF,3) =
573 : &((IFLU1.FLUX(2)*CNY+IFLU1.FLUX(3)*CTY+IFLU1.FLUX(4)*CVY) * SURF)
574 : MPOFLU.VPOCHA(NLCF,4) =
575 : &((IFLU1.FLUX(2)*CNZ+IFLU1.FLUX(3)*CTZ+IFLU1.FLUX(4)*CVZ) * SURF)
576 : MPOFLU.VPOCHA(NLCF,5) =
577 : & (IFLU1.FLUX(5) * SURF)
578 : DO I1 = 1, NESP, 1
579 : MPOFLU.VPOCHA(NLCF,5+I1)=IFLU1.FLUX(5+I1)
580 : & * SURF
581 : ENDDO
582 : DO I1 = 1, NSCA, 1
583 : MPOFLU.VPOCHA(NLCF,5+NESP+I1)=IFLU1.FLUX(5+I1+NESP)
584 : & * SURF
585 : ENDDO
586 : C
587 : C******* Calcul du pas du temps (CFL)
588 : C
589 : C****** a) etat a l'interface
590 : C
591 : DIAMG = MPOVDI.VPOCHA(NLCEG,1)
592 : DIAMD = MPOVDI.VPOCHA(NLCED,1)
593 : DIAM = MIN(DIAMG,DIAMD)
594 : CELL = 1.0D0/DIAM/CELLT
595 : IF(CELL .GT. UNSDT)THEN
596 : UNSDT = CELL
597 : ENDIF
598 : C
599 : C****** b) etat gauche
600 : C
601 : ASON = SQRT(GAMG*PG/ROG)
602 : LAMBDA = ABS(UNG) + ASON
603 : CELL = LAMBDA / DIAM
604 : IF(CELL .GT. UNSDT)THEN
605 : UNSDT = CELL
606 : ENDIF
607 : C
608 : C****** C) etat droite
609 : C
610 : ASON = SQRT(GAMD*PD/ROD)
611 : LAMBDA = ABS(UND) + ASON
612 : CELL = LAMBDA / DIAM
613 : IF(CELL .GT. UNSDT)THEN
614 : UNSDT = CELL
615 : ENDIF
616 : ENDIF
617 : C
618 : C
619 : C**** Fin boucle sur FACEL
620 : C
621 : ENDDO
622 : C
623 : C**** Pas du temps (condition de non interaction en 1D)
624 : C
625 : DT = 0.5D0 / UNSDT
626 : C
627 : C**** Desactivation des segments et
628 : C on detruit les MCHAMLs
629 : C
630 : C
631 : C**** SEGSUP FRAMAG
632 : C SEGSUP FRAMAD
633 : C
634 : C meme si LOGME = .FALSE.
635 : C
636 : SEGSUP FRAMAG
637 : SEGSUP FRAMAD
638 : C
639 : SEGSUP MLENT1
640 : SEGDES MLENT2
641 : SEGDES IPT2
642 : C
643 : SEGSUP IFLU1
644 : SEGSUP IFLU2
645 : C
646 : SEGDES MPOVSU
647 : SEGDES MPOVDI
648 : C
649 : SEGDES MPOFLU
650 : C
651 : SEGDES MELRO
652 : SEGDES MELP
653 : SEGDES MELVUN
654 : SEGDES MELVUT
655 : SEGDES MELVNX
656 : SEGDES MELVNY
657 : SEGDES MELT1X
658 : SEGDES MELT1Y
659 : c SEGDES MELLIM
660 : SEGDES MLELIM
661 : C
662 : IF(LOGME) THEN
663 : DO I1 = 1, NESP
664 : MELVA1 = MCHAMY.IELVAL(I1)
665 : SEGDES MELVA1
666 : ENDDO
667 : C
668 : SEGDES MCHAMY
669 : ENDIF
670 : IF(ISCAF .GT. 0) THEN
671 : DO I1 = 1, NSCA, 1
672 : MELVA1 = MCHAMS.IELVAL(I1)
673 : SEGDES MELVA1
674 : ENDDO
675 : C
676 : SEGDES MCHAMS
677 : ENDIF
678 : C
679 : 9999 CONTINUE
680 : C
681 : RETURN
682 : END
683 : C
684 :
685 :
686 :
687 :
688 :
689 :
690 :
691 :
692 :
693 :
694 :
695 :
696 :
697 :
698 :
699 :
700 :
© Cast3M 2003 - All rights reserved.
Disclaimer