Télécharger kfm1.eso

Retour à la liste

Numérotation des lignes :

  1. C KFM1 SOURCE CHAT 06/08/24 21:46:35 5529
  2. SUBROUTINE KFM1
  3. C************************************************************************
  4. C
  5. C PROJET : CASTEM 2000
  6. C
  7. C NOM : KFM1
  8. C
  9. C DESCRIPTION : Subroutine appellée par KONV1
  10. C
  11. C Modelisation 2D/3D des equations d'Euler
  12. C
  13. C Free matrix method: calcul d'un stationaire
  14. C avec un jacobi par point, un jacobi par ligne,
  15. C un jacobi avec 2 sweeps
  16. C
  17. C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
  18. C
  19. C AUTEUR : A. BECCANTINI, DRN/DMT/SEMT/LTMF
  20. C
  21. C************************************************************************
  22. C
  23. C APPELES (Calcul) :
  24. C
  25. C************************************************************************
  26. C
  27. C*** SYNTAXE
  28. C
  29. C Discrétisation en VF "cell-centered" des équations d'Euler pour
  30. C un gaz parfait mono-constituent polytropique
  31. C Inconnues: densité, quantité de mouvement, énergie totale par
  32. C unité de volumes (variables conservatives UN)
  33. C
  34. C A(UN) '*' DUN_L = RES(UN) + B(UN) - B(UN_L)
  35. C
  36. C ou
  37. C
  38. C RES(UN) est calculé avec ('KONV' 'VF' 'PERFMONO')
  39. C UN_0 = UN
  40. C UN_L = UN + DUN_L
  41. C
  42. C Ici on calcule DU
  43. C
  44. C DUN = 'KONV' 'VF' 'PMON1FMM' 'JACOBI' LMOT1
  45. C MOD1 MCHPO0 MCHPO1 MCHPO2 MCHPO3 MCHPO4 FLOT0 FLOT1
  46. C ENT1 ('CLIM' LMOT2 MCHPO5 MCHPO7) MCHPO6 ;
  47. C
  48. C
  49. C LMOT1 : Il contient dans l'ordre suivant: le noms de la densité,
  50. C de la qdm, de l'énergie totale par unité de volume
  51. C
  52. C MOD1 : objet modele de type EULER
  53. C
  54. C MCHPO0 : CHPOINT contenant le residu explicite; composantes =
  55. C LMOT1, SPG = (DOMA MOD1 'CENTRE').
  56. C
  57. C MCHPO1 : CHPOINT contenant la masse volumique (en kg/m^3; une
  58. C composante, 'SCAL', SPG =(DOMA MOD1 'CENTRE')).
  59. C
  60. C MCHPO2 : CHPOINT contenant les débits (en kg/s/m^2; deux
  61. C composantes en 2D, 'UX ','UY ', trois composantes
  62. C en 3D, 'UX ','UY ', 'UZ ', SPG =(DOMA MOD1 'CENTRE')).
  63. C
  64. C MCHPO3 : CHPOINT contenant l'énergie totale par unité de volume
  65. C (en J/m^3; une composante, 'SCAL', SPG =(DOMA MOD1 'CENTRE')).
  66. C
  67. C MCHPO4 : CHPOINT contenant le "gamma" du gaz (une composante,
  68. C 'SCAL', SPG =(DOMA MOD1 'CENTRE')).
  69. C
  70. C LMOT2 : Il contient dans l'ordre suivant: le noms de la densité,
  71. C de la vitesse, de la pression.
  72. C
  73. C FLOT0 : pas de temps physique
  74. C
  75. C FLOT1 : facteur de securité (double de la CFL) pour le pas de temps
  76. C dual
  77. C
  78. C ENT1 : iterations dans le jacobi
  79. C
  80. C MCHPO5 : CHPOINT contenant les conditions aux limites
  81. C (composantes en LMOT2).
  82. C
  83. C MCHPO6 : CHPOINT contenant le coeff pour le calcul du rayon
  84. C spectral visc. (une composante,
  85. C 'SCAL', SPG =(DOMA MOD1 'CENTRE')).
  86. C
  87. C MCHPO7 : CHPOINT contenant les types de conditions aux limites
  88. C
  89. C************************************************************************
  90. C
  91. C HISTORIQUE (Anomalies et modifications éventuelles)
  92. C
  93. C HISTORIQUE : crée le 30/04/02
  94. C Janvier 2003: implementation de condition aux limites
  95. C
  96. C HISTORIQUE : rajouts par T. KLOCZKO entre 2003 et 2006
  97. C
  98. C Deux types de méthodes de relaxation disponibles
  99. C
  100. C - Point Jacobi ou PJ (kfm11.eso (2D) et kfm13.eso (3D))
  101. C - Symmetric Gauss-Seidel ou SGS (kfm12.eso (2D) et kfm14.eso (3D))
  102. C
  103. C************************************************************************
  104. C
  105. IMPLICIT INTEGER(I-N)
  106. -INC CCOPTIO
  107. -INC SMLMOTS
  108. -INC SMCHPOI
  109. POINTEUR MLMVIT.MLMOTS
  110. C
  111. C**** Variables de COOPTIO
  112. C
  113. C INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR
  114. C & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ
  115. C & ,IOPER, IOSGB, IOGRA, IOSAU, IORES
  116. C & ,IECHO, IIMPI, IOSPI
  117. C & ,IDIM
  118. C & ,MCOORD
  119. C & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE
  120. C & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU
  121. C & ,NORINC,NORVAL,NORIND,NORVAD
  122. C & ,NUCROU, IPSAUV, IFICLE, IPREFI, IREFOR, ISAFOR
  123. CC
  124. INTEGER ICOEF, IDOMA, IRET, MELEF, MELEC, MELEFL, ICHPSU, INDIC
  125. & , NBCOMP, ICHPDI, IRN, IGN, IRETN, IGAMN, ICHPVO, ILIINC
  126. & , IDU, JGN, JGM, INORM, INEFMD, ICOND, MMODEL
  127. & , ICHLIM, MELLIM, ILIINP, NSOUPO, NJAC, IRES, IPROBL
  128. & , IVCO, MELTFA, ICOEV
  129. C
  130. REAL*8 DCFL, DTPS
  131. C
  132. CHARACTER*4 MOT
  133. CHARACTER*8 TYPE, LJACO(4)
  134. CHARACTER*(40) MESERR
  135. C
  136. DATA LJACO/'PJACO ','LJACOF ','LJACOB ','LJACOFB '/
  137. C
  138. CALL LIRMOT(LJACO,4,ICOEF,1)
  139. IF(IERR .NE. 0)GOTO 9999
  140. C
  141. C**********************************
  142. C**** Lecture de l'objet MODELE ***
  143. C**********************************
  144. C
  145. ICOND = 1
  146. CALL QUETYP(TYPE,ICOND,IRET)
  147.  
  148. IF((IRET.EQ.0).AND.(TYPE.NE.'MMODEL'))THEN
  149. MOTERR(1:40)='MMODEL '
  150. CALL ERREUR(471)
  151. GOTO 9999
  152. ENDIF
  153. CALL LIROBJ('MMODEL',MMODEL,ICOND,IRET)
  154. IF(IERR.NE.0)GOTO 9999
  155. CALL LEKMOD(MMODEL,IDOMA,INEFMD)
  156. IF(IERR.NE.0)GOTO 9999
  157. C
  158. C**** CENTRE, FACE, FACEL, ELTFA
  159. C
  160. CALL LEKTAB(IDOMA,'CENTRE',MELEC)
  161. IF(IERR .NE. 0) GOTO 9999
  162. C
  163. CALL LEKTAB(IDOMA,'FACE',MELEF)
  164. IF(IERR .NE. 0) GOTO 9999
  165. C
  166. CALL LEKTAB(IDOMA,'FACEL',MELEFL)
  167. IF(IERR .NE. 0) GOTO 9999
  168. C
  169. CALL LEKTAB(IDOMA,'ELTFA',MELTFA)
  170. IF(IERR .NE. 0) GOTO 9999
  171. C
  172. C**** Lecture du CHPOINT contenant les surfaces des faces.
  173. C
  174. CALL LEKTAB(IDOMA,'XXSURFAC',ICHPSU)
  175. IF(IERR .NE. 0) GOTO 9999
  176. INDIC = 1
  177. NBCOMP = 1
  178. MOT = 'SCAL'
  179. CALL QUEPOI(ICHPSU, MELEF, INDIC, NBCOMP, MOT)
  180. IF(IERR .NE. 0) GOTO 9999
  181. C
  182. C**** Lecture du CHPOINT contenant les diametres minimums.
  183. C
  184. CALL LEKTAB(IDOMA,'XXDIEMIN',ICHPDI)
  185. IF(IERR .NE. 0) GOTO 9999
  186. INDIC = 1
  187. NBCOMP = 1
  188. MOT = 'SCAL'
  189. CALL QUEPOI(ICHPDI, MELEC, INDIC, NBCOMP, MOT)
  190. IF(IERR .NE. 0) GOTO 9999
  191.  
  192. C
  193. C**** Lecture du CHPOINT contenant les normales aux faces
  194. C
  195. IF(IDIM .EQ. 2)THEN
  196. C Que les normales
  197. CALL LEKTAB(IDOMA,'XXNORMAF',INORM)
  198. IF(IERR .NE. 0) GOTO 9999
  199. JGN = 4
  200. JGM = 2
  201. SEGINI MLMVIT
  202. MLMVIT.MOTS(1) = 'UX '
  203. MLMVIT.MOTS(2) = 'UY '
  204. CALL QUEPO1(INORM, MELEF, MLMVIT)
  205. SEGSUP MLMVIT
  206. C
  207. ELSE
  208. C Les normales et les tangentes
  209. TYPE = ' '
  210. CALL ACMO(IDOMA,'MATROT',TYPE,INORM)
  211. IF (TYPE .NE. 'CHPOINT ') THEN
  212. CALL MATRAN(IDOMA,INORM)
  213. IF(IERR .NE. 0) GOTO 9999
  214. ENDIF
  215. JGN = 4
  216. JGM = 9
  217. SEGINI MLMVIT
  218. MLMVIT.MOTS(1) = 'UX '
  219. MLMVIT.MOTS(2) = 'UY '
  220. MLMVIT.MOTS(3) = 'UZ '
  221. MLMVIT.MOTS(4) = 'RX '
  222. MLMVIT.MOTS(5) = 'RY '
  223. MLMVIT.MOTS(6) = 'RZ '
  224. MLMVIT.MOTS(7) = 'MX '
  225. MLMVIT.MOTS(8) = 'MY '
  226. MLMVIT.MOTS(9) = 'MZ '
  227. CALL QUEPO1(INORM, MELEF, MLMVIT)
  228. SEGSUP MLMVIT
  229. ENDIF
  230. C
  231. C**** Lecture du CHPOINT contenant les volumes
  232. C
  233. CALL LEKTAB(IDOMA,'XXVOLUM',ICHPVO)
  234. IF(IERR .NE. 0) GOTO 9999
  235. INDIC = 1
  236. NBCOMP = 1
  237. MOT = 'SCAL'
  238. CALL QUEPOI(ICHPVO, MELEC, INDIC, NBCOMP, MOT)
  239. IF(IERR .NE. 0) GOTO 9999
  240. C
  241. C********************************
  242. C**** Fin table domaine *********
  243. C********************************
  244. C
  245. C**** La list des inconnues
  246. C
  247. TYPE='LISTMOTS'
  248. CALL LIROBJ(TYPE,ILIINC,1,IRET)
  249. IF(IERR .NE. 0) GOTO 9999
  250. C
  251. C**** On va lire les pointeurs des CHPOINTs
  252. C Lecture du CHPOINT residu
  253. C
  254. TYPE='CHPOINT'
  255. CALL LIROBJ(TYPE,IRES,1,IRET)
  256. IF(IERR.NE.0) GOTO 9999
  257. C
  258. C**** Controle du CHPOINT
  259. C
  260. CALL QUEPO1(IRES, MELEC, ILIINC)
  261. IF(IERR .NE. 0) GOTO 9999
  262. C
  263. C**** On va lire les pointeurs des CHPOINTs
  264. C Lecture du CHPOINT centre densité
  265. C
  266. TYPE='CHPOINT'
  267. CALL LIROBJ(TYPE,IRN,1,IRET)
  268. IF(IERR.NE.0) GOTO 9999
  269. C
  270. C**** Controle du CHPOINT
  271. C
  272. INDIC = 1
  273. NBCOMP = 1
  274. MOT = 'SCAL'
  275. CALL QUEPOI(IRN, MELEC, INDIC, NBCOMP, MOT)
  276. IF(IERR .NE. 0) GOTO 9999
  277. C
  278. C**** Lecture du CHPOINT QDM
  279. C
  280. TYPE='CHPOINT'
  281. CALL LIROBJ(TYPE,IGN,1,IRET)
  282. IF(IERR.NE.0) GOTO 9999
  283. C
  284. C**** Control du CHPOINT
  285. C
  286. INDIC = 1
  287. NBCOMP = IDIM
  288. JGN = 4
  289. JGM = IDIM
  290. SEGINI MLMOTS
  291. MLMOTS.MOTS(1) = 'UX '
  292. MLMOTS.MOTS(2) = 'UY '
  293. IF(IDIM .EQ. 3) MLMOTS.MOTS(3) = 'UZ '
  294. CALL QUEPO1(IGN, MELEC, MLMOTS)
  295. IF(IERR .NE. 0) GOTO 9999
  296. SEGSUP MLMOTS
  297. C
  298. C**** Lecture du CHPOINT centre energie totale
  299. C
  300. TYPE='CHPOINT'
  301. CALL LIROBJ(TYPE,IRETN,1,IRET)
  302. IF(IERR.NE.0) GOTO 9999
  303. C
  304. C**** Controle du CHPOINT
  305. C
  306. INDIC = 1
  307. NBCOMP = 1
  308. MOT = 'SCAL'
  309. CALL QUEPOI(IRETN, MELEC, INDIC, NBCOMP, MOT)
  310. IF(IERR .NE. 0) GOTO 9999
  311. C
  312. C**** Lecture du CHPOINT centre gamma
  313. C
  314. TYPE='CHPOINT'
  315. CALL LIROBJ(TYPE,IGAMN,1,IRET)
  316. IF(IERR.NE.0) GOTO 9999
  317. C
  318. C**** Controle du CHPOINT
  319. C
  320. INDIC = 1
  321. NBCOMP = 1
  322. MOT = 'SCAL'
  323. CALL QUEPOI(IGAMN, MELEC, INDIC, NBCOMP, MOT)
  324. IF(IERR .NE. 0) GOTO 9999
  325. C
  326. C**** Lecture du CHPOINT centre cutoff speed
  327. C
  328. TYPE='CHPOINT'
  329. CALL LIROBJ(TYPE,IVCO,1,IRET)
  330. IF(IERR.NE.0) GOTO 9999
  331. C
  332. C**** Controle du CHPOINT
  333. C
  334. INDIC = 1
  335. NBCOMP = 1
  336. MOT = 'SCAL'
  337. CALL QUEPOI(IVCO, MELEC, INDIC, NBCOMP, MOT)
  338. IF(IERR .NE. 0) GOTO 9999
  339. C
  340. C**** Lecture de DTPS
  341. C
  342. CALL LIRREE(DTPS,1,IRET)
  343. IF(IERR .NE. 0) GOTO 9999
  344. C
  345. C**** Lecture de le double de la CFL pour le temps dual
  346. C
  347. CALL LIRREE(DCFL,1,IRET)
  348. IF(IERR .NE. 0) GOTO 9999
  349. C
  350. C**** Lecture de nombre d'iterations
  351. C
  352. CALL LIRENT(NJAC,1,IRET)
  353. IF(IERR .NE. 0) GOTO 9999
  354. C
  355. C**** Lecture de conditions aux limites
  356. C
  357. IRET=0
  358. CALL LIRCHA(MOT,0,IRET)
  359. IF(IERR .NE. 0) GOTO 9999
  360. IF(IRET .NE. 0)THEN
  361. IF(MOT .EQ. 'CLIM')THEN
  362. C
  363. TYPE='LISTMOTS'
  364. CALL LIROBJ(TYPE,ILIINP,1,IRET)
  365. IF(IERR .NE. 0) GOTO 9999
  366. MLMOT1=ILIINP
  367. C
  368. TYPE='CHPOINT'
  369. CALL LIROBJ(TYPE,ICHLIM,1,IRET)
  370. IF(IERR.NE.0) GOTO 9999
  371. C
  372. MCHPOI = ICHLIM
  373. SEGACT MCHPOI
  374. NSOUPO = MCHPOI.IPCHP(/1)
  375. IF(NSOUPO .EQ. 0) THEN
  376. ICHLIM=0
  377. MELLIM=0
  378. ELSEIF(NSOUPO .GT. 1)THEN
  379. MOTERR(1:8) = 'CHAMPOIN'
  380. C
  381. C**************** Message d'erreur standard
  382. C 132 2
  383. C On veut un objet %m1:8 élémentaire
  384. C
  385. CALL ERREUR(132)
  386. GOTO 9999
  387. ELSE
  388. MSOUPO=MCHPOI.IPCHP(1)
  389. SEGACT MSOUPO
  390. MELLIM=MSOUPO.IGEOC
  391. SEGDES MSOUPO
  392. SEGDES MCHPOI
  393. CALL QUEPO1(ICHLIM, MELLIM, MLMOT1)
  394. IF(IERR.NE.0) GOTO 9999
  395. ENDIF
  396. C
  397. ELSE
  398. CALL REFUS
  399. ENDIF
  400. ELSE
  401. ICHLIM=0
  402. MELLIM=0
  403. ENDIF
  404. C
  405. C**** Le coeff. pour le calcul du rayon spectral visc.
  406. C
  407. C
  408. TYPE='CHPOINT'
  409. CALL LIROBJ(TYPE,ICOEV,1,IRET)
  410. IF(IERR.NE.0) GOTO 9999
  411. C
  412. C**** Controle du CHPOINT
  413. C
  414. INDIC = 1
  415. NBCOMP = 1
  416. MOT = 'SCAL'
  417. CALL QUEPOI(ICOEV, MELEC, INDIC, NBCOMP, MOT)
  418. IF(IERR .NE. 0) GOTO 9999
  419. C
  420. C**** Calcul de DUN
  421. C
  422. TYPE='CHPOINT'
  423. CALL KRCHP1(TYPE, MELEC, IDU, ILIINC)
  424. C
  425. C*** Choix de la procédure de relaxation (PJ ou SGS)
  426. C
  427. IF(ICOEF .EQ. 1)THEN
  428. IF(IDIM .EQ. 2)THEN
  429. C
  430. C*** Procédure de relaxation PJ
  431. C
  432. CALL KFM11(IRES,IRN,IGN,IRETN,IGAMN,IVCO,
  433. & ICHPSU,ICHPDI,ICHPVO,INORM,
  434. & MELEC,MELEF,MELEFL,DTPS,DCFL,
  435. & MELLIM,ICHLIM,NJAC,ICOEV,
  436. & IDU,IPROBL)
  437. IF(IERR .NE. 0)GOTO 9999
  438. ELSE
  439. CALL KFM13(IRES,IRN,IGN,IRETN,IGAMN,IVCO,
  440. & ICHPSU,ICHPDI,ICHPVO,INORM,
  441. & MELEC,MELEF,MELEFL,DTPS,DCFL,
  442. & MELLIM,ICHLIM,NJAC,ICOEV,
  443. & IDU,IPROBL)
  444. IF(IERR .NE. 0)GOTO 9999
  445. c 251 2
  446. c Tentative d'utilisation d'une option non implémentée
  447. c CALL ERREUR(251)
  448. c GOTO 9999
  449. ENDIF
  450. ELSE
  451. IF(IDIM .EQ. 2)THEN
  452. C
  453. C*** Procédure de relaxation SGS
  454. C
  455. CALL KFM12(IRES,IRN,IGN,IRETN,IGAMN,IVCO,
  456. & ICHPSU,ICHPDI,ICHPVO,INORM,
  457. & MELEC,MELEF,MELEFL,MELTFA,DTPS,DCFL,
  458. & MELLIM,ICHLIM,NJAC,ICOEF,ICOEV,
  459. & IDU,IPROBL)
  460. IF(IERR .NE. 0)GOTO 9999
  461. ELSE
  462. CALL KFM14(IRES,IRN,IGN,IRETN,IGAMN,IVCO,
  463. & ICHPSU,ICHPDI,ICHPVO,INORM,
  464. & MELEC,MELEF,MELEFL,MELTFA,DTPS,DCFL,
  465. & MELLIM,ICHLIM,NJAC,ICOEF,ICOEV,
  466. & IDU,IPROBL)
  467. C 251 2
  468. C Tentative d'utilisation d'une option non implémentée
  469. C CALL ERREUR(251)
  470. C GOTO 9999
  471. ENDIF
  472. ENDIF
  473. C
  474. CALL ECRENT(IPROBL)
  475. TYPE='CHPOINT '
  476. CALL ECROBJ(TYPE,IDU)
  477. C
  478. 9999 CONTINUE
  479. RETURN
  480. END
  481.  
  482.  
  483.  
  484.  
  485.  

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