Télécharger ckon3.eso

Retour à la liste

Numérotation des lignes :

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

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