Télécharger konv2f.eso

Retour à la liste

Numérotation des lignes :

  1. C KONV2F SOURCE PV 16/11/17 22:00:14 9180
  2. SUBROUTINE KONV2F()
  3. C************************************************************************
  4. C
  5. C PROJET : CASTEM 2000
  6. C
  7. C NOM : KONV2F
  8. C
  9. C DESCRIPTION : Subroutine appellée par KONV1
  10. C
  11. C Modelling of Two Phase Flow
  12. C (two fluid model, mixtures of water and air)
  13. C
  14. C Calcul du residu / jacobien / DELTAT
  15. C
  16. C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
  17. C
  18. C AUTEUR : A. BECCANTINI, DRN/DMT/SEMT/LTMF
  19. C Modified for two phase flow by
  20. C Jose R. Garcia-Cascales,
  21. C Universidad Politecnica de Cartagena,
  22. C jr.garcia@upct.es
  23. C************************************************************************
  24. C
  25. C APPELES (Calcul) : KONFL8 (2D, two fluid flow)
  26. C KONFL9 (3D, two fluid flow)
  27. C
  28. C************************************************************************
  29. C
  30. C*** SYNTAXE
  31. C FV Discretization of the Two Phase Flow system of equations
  32. C Unknown: mass flow rate, momentum and total energy per unit volume of
  33. C each phase (conservative variables)
  34. C
  35. C
  36. C RMAT1 RCHPO1 RFLOT1 = 'KONV' 'VF' 'PERFMONO' MOT1 MOT2
  37. C TABD LMOT1
  38. C MCHAM1 MCHAM2 MCHAM3 MCHAM4
  39. C MCHAM5 MCHAM6 MCHAM7 MCHAM8
  40. C (MOT3 CHPO1 CHPO2 CHPO3 CHPO5) ;
  41. C
  42. C ENTREES
  43. C
  44. C
  45. C MOT1 : objet de type MOT
  46. C Il vaut 'RESI' si on veut calculer le résidu
  47. C Il vaut 'FLUX' si on veut calculer le flux
  48. C
  49. C MOT2 : objet de type MOT
  50. C Il indique la méthode de décentrement:
  51. C 'AUSMP1' = AUSM +
  52. C 'AUSMP2' = Preconditioned AUSM +
  53. C 'AUSMDV1'= AUSMDV
  54. C 'AUSMDV2'= Preconditioned AUSMDV
  55. C
  56. C LMOT1 : objet de type LISTMOTS
  57. C Noms de composantes du résultat (RCHPO1)
  58. C Il contient dans l'ordre suivant: le noms de la densité,
  59. C de la vitesse, de l'énergie totale par unité de volume
  60. C
  61. C TAB1 : la table domaine
  62. C
  63. C MCHAM1 : MCHAML containing void fraction
  64. C SPG (support géométrique) l'indice 'FACEL' de
  65. C TAB1 (une composante, 'SCAL')
  66. C
  67. C MCHAM2 : MCHAML vapour velocity and director cosines
  68. C du repère locale (n,t) dans le repère
  69. C global (x,y) (dans le cas 2D 6 composantes:
  70. C * 'UN' = vitesse normale (SPG = TAB1 . 'FACEL')
  71. C * 'UT' = vitesse tangentielle (SPG = 'TAB1 . FACEL')
  72. C * 'NX' = n.x (SPG = 'FACE')
  73. C * 'NY' = n.y (SPG = 'FACE')
  74. C * 'TX' = t.x (SPG = 'FACE')
  75. C * 'TY' = t.y (SPG = 'FACE')).
  76. C MCHAM3 : MCHAML liquid velocity and director cosines
  77. C du repère locale (n,t) dans le repère
  78. C global (x,y) (dans le cas 2D 6 composantes:
  79. C * 'UN' = vitesse normale (SPG = TAB1 . 'FACEL')
  80. C * 'UT' = vitesse tangentielle (SPG = 'TAB1 . FACEL')
  81. C * 'NX' = n.x (SPG = 'FACE')
  82. C * 'NY' = n.y (SPG = 'FACE')
  83. C * 'TX' = t.x (SPG = 'FACE')
  84. C * 'TY' = t.y (SPG = 'FACE')).
  85. C
  86. C MCHAM4 : MCHAML (SPG = TAB1 . 'FACEL') containing pressure
  87. C une seule composante, 'SCAL').
  88. C
  89. C MCHAM5 : MCHAML (SPG = TAB1 . 'FACEL') contaning gas/vapour temperature
  90. C (une seule composante, 'SCAL').
  91. C
  92. C MCHAM6 : MCHAML (SPG = TAB1 . 'FACEL') contaning liquid temperature
  93. C (une seule composante, 'SCAL').
  94. C
  95. C MCHAM7 : MCHAML (SPG = TAB1 . 'FACEL') contaning gas/vapour density
  96. C (une seule composante, 'SCAL').
  97. C
  98. C MCHAM8 : MCHAML (SPG = TAB1 . 'FACEL') contaning liquid density
  99. C (une seule composante, 'SCAL').
  100. C
  101. C RCHPO1 : objet de type CHPOINT (composantes = LMOT1)
  102. C Residu si MOT2 = 'RESI' (SPG = TAB1 . 'CENTRE')
  103. C Flux si MOT2 = 'FLUX' (SPG = TAB1 . 'FACE')
  104. C
  105. C RFLOT1 : objet de type FLOTTANT
  106. C Il est le temps caracteristique associé à l'onde la plus
  107. C rapide.
  108. C
  109. C Remarque
  110. C --------
  111. C
  112. C RCHPO1 est égal à:
  113. C * la derivé temporelle des inconnues si l'option 'RESI' est utilisée
  114. C * la projection du flux convectif sur (TAB1 . 'XXNORMAF') si
  115. C l'option 'FLUX' est utilisée
  116. C
  117. C***********************************************************************
  118. C
  119. C************************************************************************
  120. C
  121. C HISTORIQUE (Anomalies et modifications éventuelles)
  122. C
  123. C HISTORIQUE :
  124. C
  125. C************************************************************************
  126. C
  127. IMPLICIT INTEGER(I-N)
  128. IMPLICIT REAL*8(A-H,O-Z)
  129.  
  130. -INC CCOPTIO
  131. -INC SMLMOTS
  132. -INC SMCHPOI
  133. POINTEUR MLMVIT.MLMOTS,MLMOTY.MLMOTS
  134. C
  135. INTEGER IDOMA, MELEMC, MELEMF, MELEFE, ICHPSU, ICHPDI , ICHPVO
  136. & , MELTFA, NBMET, INDMET, IRET, NBJAC
  137. & , IAIF, IALPF, IUVF, IULF, IPF, ITVF, ITLF
  138. & , IRVF, IRLF, IFRMAF
  139. & , ILIINC, NC, ICELL
  140. & , INDIC, NBCOMP, IJACO
  141. & , NKID, NKMT, NMATRI, NRIGE
  142. & , JGM, JGN, INORM, IIMPL
  143. & , ICHFLU, ICHRES
  144.  
  145. C
  146. PARAMETER (NBMET=4,NBJAC=5)
  147. REAL*8 DT
  148. CHARACTER*8 LMETO(NBMET), TYPE, LJACO(NBJAC)
  149. CHARACTER*4 MOT, LFLUX(2)
  150. CHARACTER*(40) MESERR
  151. LOGICAL LOGNC, LOGAN, LOGRES
  152. C**** Depending on if we want to use AUSM+, AUSM+ with "preconditioning",
  153. C AUSMDV or AUSMDV with "preconditioning" we have respectively
  154. DATA LMETO/'AUSMP1','AUSMP2','AUSMDV1','AUSMDV2'/
  155. C**** Option not implemented yet
  156. DATA LJACO/'VLHJACO1','AUSMJAC1','AUSMLMJ1','VLHJACO2',
  157. & 'AUSMJAC2'/
  158. C**** If we want flux or residual
  159. DATA LFLUX/'FLUX','RESI'/
  160. C
  161. C**** Initialisation des variables pour la gestion des erreurs.
  162. C
  163. LOGNC = .FALSE.
  164. LOGAN = .FALSE.
  165. MESERR = ' '
  166. C
  167. C******* Flux ou residu???
  168. C
  169. CALL LIRMOT(LFLUX,2,ICELL,1)
  170. IF(IERR .NE. 0)GOTO 9999
  171. IF(ICELL .EQ. 1)THEN
  172. LOGRES = .FALSE.
  173. ELSEIF(ICELL .EQ. 2)THEN
  174. LOGRES = .TRUE.
  175. ELSE
  176. C
  177. C******** Message d'erreur standard
  178. C 251 2
  179. C Tentative d'utilisation d'une option non implémentée
  180. C
  181. CALL ERREUR(251)
  182. ENDIF
  183. C
  184. C**** Metode utilisée
  185. C
  186. CALL LIRMOT(LMETO,NBMET,INDMET,1)
  187. IF(IERR .NE. 0)GOTO 9999
  188. IF(INDMET .EQ. 0)THEN
  189. C
  190. C******** Message d'erreur standard
  191. C 251 2
  192. C Tentative d'utilisation d'une option non implémentée
  193. C
  194. CALL ERREUR(251)
  195. ENDIF
  196. C
  197. C*******************************
  198. C**** La table domaine *********
  199. C*******************************
  200. C
  201. CALL LIROBJ('TABLE',IDOMA,1,IRET)
  202. IF(IERR .NE. 0)GOTO 9999
  203. C
  204. C**** Centre, FACE, FACEL, ELTFA
  205. C
  206. CALL LEKTAB(IDOMA,'CENTRE',MELEMC)
  207. IF(IERR .NE. 0) GOTO 9999
  208. C
  209. CALL LEKTAB(IDOMA,'FACE',MELEMF)
  210. IF(IERR .NE. 0) GOTO 9999
  211. C
  212. CALL LEKTAB(IDOMA,'FACEL',MELEFE)
  213. IF(IERR .NE. 0) GOTO 9999
  214. C
  215. CALL LEKTAB(IDOMA,'ELTFA',MELTFA)
  216. IF(IERR .NE. 0) GOTO 9999
  217. C
  218. C**** Lecture du CHPOINT contenant les surfaces des faces.
  219. C
  220. CALL LEKTAB(IDOMA,'XXSURFAC',ICHPSU)
  221. IF(IERR .NE. 0) GOTO 9999
  222. INDIC = 1
  223. NBCOMP = 1
  224. MOT = 'SCAL'
  225. CALL QUEPOI(ICHPSU, MELEMF, INDIC, NBCOMP, MOT)
  226. IF(IERR .NE. 0) GOTO 9999
  227. C
  228. C**** Lecture du CHPOINT contenant les diametres minimums.
  229. C
  230. CALL LEKTAB(IDOMA,'XXDIEMIN',ICHPDI)
  231. IF(IERR .NE. 0) GOTO 9999
  232. INDIC = 1
  233. NBCOMP = 1
  234. MOT = 'SCAL'
  235. CALL QUEPOI(ICHPDI, MELEMC, INDIC, NBCOMP, MOT)
  236. IF(IERR .NE. 0) GOTO 9999
  237.  
  238. C
  239. C**** Lecture du CHPOINT contenant les volumes
  240. C
  241. CALL LEKTAB(IDOMA,'XXVOLUM',ICHPVO)
  242. IF(IERR .NE. 0) GOTO 9999
  243. INDIC = 1
  244. NBCOMP = 1
  245. MOT = 'SCAL'
  246. CALL QUEPOI(ICHPVO, MELEMC, INDIC, NBCOMP, MOT)
  247. IF(IERR .NE. 0) GOTO 9999
  248. C
  249. C**** Les normales aux faces
  250. C
  251. IF(IDIM .EQ. 2)THEN
  252. C Que les normales
  253. CALL LEKTAB(IDOMA,'XXNORMAF',INORM)
  254. IF(IERR .NE. 0) GOTO 9999
  255. JGN = 4
  256. JGM = 2
  257. SEGINI MLMVIT
  258. MLMVIT.MOTS(1) = 'UX '
  259. MLMVIT.MOTS(2) = 'UY '
  260. CALL QUEPO1(INORM, MELEMF, MLMVIT)
  261. SEGSUP MLMVIT
  262. IF(IERR .NE. 0) GOTO 9999
  263. C De momento solo 2D
  264. ELSE
  265. C Les normales et les tangentes
  266. TYPE = ' '
  267. CALL ACMO(IDOMA,'MATROT',TYPE,INORM)
  268. IF (TYPE .NE. 'CHPOINT ') THEN
  269. CALL MATRAN(IDOMA,INORM)
  270. IF(IERR .NE. 0) GOTO 9999
  271. ENDIF
  272. JGN = 4
  273. JGM = 9
  274. SEGINI MLMVIT
  275. MLMVIT.MOTS(1) = 'UX '
  276. MLMVIT.MOTS(2) = 'UY '
  277. MLMVIT.MOTS(3) = 'UZ '
  278. MLMVIT.MOTS(4) = 'RX '
  279. MLMVIT.MOTS(5) = 'RY '
  280. MLMVIT.MOTS(6) = 'RZ '
  281. MLMVIT.MOTS(7) = 'MX '
  282. MLMVIT.MOTS(8) = 'MY '
  283. MLMVIT.MOTS(9) = 'MZ '
  284. CALL QUEPO1(INORM, MELEMF, MLMVIT)
  285. SEGSUP MLMVIT
  286. ENDIF
  287. C
  288. C********************************
  289. C**** Fin table domaine *********
  290. C********************************
  291. C
  292. C
  293. C**** On va lire les pointeurs des MCHAMLs
  294. C Reading the MCHAML 'FACEL' containig void fraction
  295. C
  296. TYPE='MCHAML '
  297. CALL LIROBJ(TYPE,IALPF,1,IRET)
  298. IF(IERR.NE.0) GOTO 9999
  299. C
  300. C**** Reading the MCHAML 'FACEL' containing vapour velocity
  301. C
  302. TYPE='MCHAML '
  303. CALL LIROBJ(TYPE,IUVF,1,IRET)
  304. IF(IERR .NE. 0) GOTO 9999
  305. C
  306. C**** Reading the MCHAML 'FACEL' containing liquid velocity
  307. C
  308. TYPE='MCHAML '
  309. CALL LIROBJ(TYPE,IULF,1,IRET)
  310. IF(IERR .NE. 0) GOTO 9999
  311. C
  312. C**** Reading the MCHAML 'FACEL' containing pressure
  313. C
  314. TYPE='MCHAML '
  315. CALL LIROBJ(TYPE,IPF,1,IRET)
  316. IF(IERR .NE. 0) GOTO 9999
  317. C
  318. C**** Reading the MCHAML 'FACEL' containing vapour temperature
  319. C
  320. TYPE='MCHAML '
  321. CALL LIROBJ(TYPE,ITVF,1,IRET)
  322. IF(IERR .NE. 0) GOTO 9999
  323. C
  324. C**** Reading the MCHAML 'FACEL' containing liquid temperature
  325. C
  326. TYPE='MCHAML '
  327. CALL LIROBJ(TYPE,ITLF,1,IRET)
  328. IF(IERR .NE. 0) GOTO 9999
  329. C
  330. C**** Reading the MCHAML 'FACEL' containing vapour density
  331. C
  332. TYPE='MCHAML '
  333. CALL LIROBJ(TYPE,IRVF,1,IRET)
  334. IF(IERR .NE. 0) GOTO 9999
  335. C
  336. C**** Reading the MCHAML 'FACEL' containing liquid density
  337. C
  338. TYPE='MCHAML '
  339. CALL LIROBJ(TYPE,IRLF,1,IRET)
  340. IF(IERR .NE. 0) GOTO 9999
  341.  
  342. C
  343. C**** La list des inconnues
  344. C
  345. TYPE='LISTMOTS'
  346. CALL LIROBJ(TYPE,ILIINC,1,IRET)
  347. IF(IERR .NE. 0) GOTO 9999
  348. MLMOTS = ILIINC
  349. SEGACT MLMOTS
  350. NC = MLMOTS.MOTS(/2)
  351. SEGDES MLMOTS
  352. IF(IDIM .EQ. 3)THEN
  353. C*** In 3D we have 10 equations
  354. NUMEQ = IDIM + 7
  355. ELSE
  356. C*** In 2D we have 8 equations
  357. NUMEQ = IDIM + 6
  358. END IF
  359. IF(NC .NE. NUMEQ)THEN
  360. MOTERR(1:40) = 'LISTINCO = ???'
  361. write(*,*) NC
  362. WRITE(IOIMP,*) MOTERR
  363. C
  364. C******* Message d'erreur standard
  365. C 21 2
  366. C Données incompatibles
  367. C
  368. CALL ERREUR(21)
  369. GOTO 9999
  370. ENDIF
  371. C
  372. C**** Implicite ?
  373. C
  374. CALL LIRMOT(LJACO,NBJAC,IIMPL,0)
  375. IF(IERR .NE. 0)GOTO 9999
  376.  
  377. C For the moment we only have explicit
  378. C
  379. C
  380. IF(IIMPL .EQ. 0)THEN
  381. C
  382. C******** Objet matrik vide
  383. C
  384. NRIGE=7
  385. NMATRI=0
  386. NKID =9
  387. NKMT =7
  388. SEGINI MATRIK
  389. SEGDES MATRIK
  390. IJACO = MATRIK
  391. ELSE
  392. C Tentative d'utilisation d'une option non implémentée
  393. CALL ERREUR(251)
  394. GOTO 9999
  395. ENDIF
  396. C
  397. C**** Creation des flux aux interfaces
  398. C
  399. TYPE = 'CHPOINT '
  400. CALL KRCHP1(TYPE, MELEMF, ICHFLU, MLMOTS)
  401. C
  402. C
  403. C**** Calcul des flux et du pas du temps.
  404. C
  405. IF(IDIM .EQ. 2)THEN
  406. CALL KONF21(INDMET,
  407. & IALPF,IUVF,IULF,IPF,
  408. & ITVF,ITLF,IRVF,IRLF,
  409. & ICHPSU,ICHPDI,
  410. & MELEMC,MELEMF,MELEFE,
  411. & ICHFLU,
  412. & DT,
  413. & LOGNC,LOGAN,MESERR)
  414. ELSEIF(IDIM .EQ. 3)THEN
  415. CALL KONF22(INDMET,
  416. & IALPF,IUVF,IULF,IPF,
  417. & ITVF,ITLF,IRVF,IRLF,
  418. & ICHPSU,ICHPDI,
  419. & MELEMC,MELEMF,MELEFE,
  420. & ICHFLU,DT,
  421. & LOGNC,LOGAN,MESERR)
  422. ELSE
  423. C Erreur anormale.contactez votre support
  424. C Dimension has to be 2 or 3
  425. CALL ERREUR(5)
  426. GOTO 9999
  427. ENDIF
  428. C
  429. IF(LOGAN)THEN
  430. C
  431. C******* Anomalie detectée
  432. C
  433. C
  434. C******* Message d'erreur standard
  435. C -301 0
  436. C %m1:40
  437. C
  438. MOTERR(1:40) = MESERR(1:40)
  439. WRITE(IOIMP,*) MOTERR(1:40)
  440. C
  441. C******* Message d'erreur standard
  442. C 5 3
  443. C Erreur anormale.contactez votre support
  444. C
  445. CALL ERREUR(5)
  446. GOTO 9999
  447. ENDIF
  448. IF(LOGNC)THEN
  449. C
  450. C******* Message d'erreur standard
  451. C -301 0
  452. C %m1:40
  453. C
  454. MOTERR(1:40) = MESERR(1:40)
  455. WRITE(IOIMP,*) MOTERR(1:40)
  456. C
  457. C******* Message d'erreur standard
  458. C 460 2
  459. C Pas de convergence dans les itérations internes
  460. C
  461. CALL ERREUR(460)
  462. GOTO 9999
  463. ENDIF
  464. C
  465. C**** Calcul de residu (si LOGRES = .TRUE.)
  466. C
  467. C LOGRES = .FALSE.
  468.  
  469. IF(LOGRES)THEN
  470. TYPE = 'CHPOINT '
  471. CALL KRCHP1(TYPE, MELEMC, ICHRES, MLMOTS)
  472. C
  473. CALL KONRE8(MELEMC,MELEMF,MELEFE,ICHPVO,
  474. & ICHFLU, ICHRES,
  475. & LOGAN,MESERR)
  476. IF(LOGAN)THEN
  477. C
  478. C******* Anomalie detectée
  479. C
  480. C
  481. C******* Message d'erreur standard
  482. C -301 0
  483. C %m1:40
  484. C
  485. MOTERR(1:40) = MESERR(1:40)
  486. WRITE(IOIMP,*) MOTERR(1:40)
  487. C
  488. C******* Message d'erreur standard
  489. C 5 3
  490. C Erreur anormale.contactez votre support
  491. C
  492. CALL ERREUR(5)
  493. GOTO 9999
  494. ENDIF
  495. ELSE
  496. ICHRES = 0
  497. ENDIF
  498. C
  499. C**** Ecriture des resultats
  500. C
  501. CALL ECRREE(DT)
  502. TYPE = 'CHPOINT '
  503. IF(ICHRES .NE. 0) CALL ECROBJ(TYPE,ICHRES)
  504. IF(ICHFLU .NE. 0) CALL ECROBJ(TYPE,ICHFLU)
  505. TYPE='MATRIK '
  506. C CALL ECROBJ(TYPE,IJACO)
  507. 9999 CONTINUE
  508. RETURN
  509. END
  510.  
  511.  
  512.  
  513.  
  514.  
  515.  
  516.  
  517.  
  518.  
  519.  
  520.  

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