Télécharger pre122.eso

Retour à la liste

Numérotation des lignes :

pre122
  1. C PRE122 SOURCE OF166741 24/10/03 21:15:29 12022
  2. C PRE121 SOURCE BECC 98/08/24 21:17:00 3286
  3. SUBROUTINE PRE122(LOGTEM,
  4. & ICEN,IFACE,IFACEL,INORM,
  5. & IROC, IGRROC, IALROC,
  6. & IVITC, IGRVC, IALVC,
  7. & IPC ,IGRPC, IALPC,
  8. & IGAMC,
  9. & DELTAT,
  10. & IROF,IVITF,IPF,IGAMF,
  11. & LOGAN,LOGNEG,LOGBOR,MESERR,VALER,VAL1,VAL2)
  12. C************************************************************************
  13. C
  14. C PROJET : CASTEM 2000
  15. C
  16. C NOM : PRE121
  17. C
  18. C DESCRIPTION : Voir PRE12
  19. C
  20. C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI)
  21. C
  22. C AUTEUR : A. BECCANTINI, DRN/DMT/SEMT/TTMF
  23. C
  24. C************************************************************************
  25. C
  26. C
  27. C APPELES (Outils) : KRIPAD, LICHT
  28. C
  29. C APPELES (Calcul) : AUCUN
  30. C
  31. C
  32. C************************************************************************
  33. C
  34. C ENTREES
  35. C
  36. C LOGTEM : LOGICAL; si .TRUE. 2em ordre en temps
  37. C sinon 1er ordre en temps
  38. C
  39. C 1) Pointeurs de MELEMEs et de CHPOINTs de la table DOMAINE
  40. C
  41. C ICEN : MELEME de 'POI1' SPG des CENTRES
  42. C
  43. C IFACE : MELEME de 'POI1' SPG des FACES
  44. C
  45. C IFACEL : MELEME de 'SEG3' avec
  46. C CENTRE d'Elt "gauche"
  47. C CENTRE de Face
  48. C CENTRE d'Elt "droite"
  49. C
  50. C N.B. = IFACE.NUM(i,1) = IFACEL.NUM(i,2)
  51. C
  52. C INORM : CHPOINT des cosinus directeurs de normales aux faces
  53. C
  54. C 2) Pointeurs des CHPOINTs
  55. C
  56. C IROC : CHPOINT "CENTRE" contenant la masse volumique RHO
  57. C
  58. C IGRROC : CHPOINT "CENTRE" contenant le gradient de la
  59. C masse volumique RHO (2 composantes)
  60. C
  61. C IALROC : CHPOINT "CENTRE" contenant le limiteur du gradient
  62. C de la masse volumique
  63. C
  64. C IVITC : CHPOINT "CENTRE" contenant la vitesse UX, UY ;
  65. C
  66. C IGRVC : CHPOINT "CENTRE" contenant le gradient de la
  67. C vitesse (4 composantes)
  68. C
  69. C IALVC : CHPOINT "CENTRE" contenant le limiteur du gradient
  70. C de la vitesse (2 composantes)
  71. C
  72. C IPC : CHPOINT "CENTRE" contenat la pression P;
  73. C
  74. C IGRPC : CHPOINT "CENTRE" contenant le gradient de la
  75. C pression (2 composantes)
  76. C
  77. C IALPC : CHPOINT "CENTRE" contenant le limiteur du gradient
  78. C de la pression
  79. C
  80. C IGAMC : CHPOINT "CENTRE" contenat le "Gamma" du gaz
  81. C
  82. C 3)
  83. C
  84. C DELTAT : REAL*8, encrement en temps pour calculer la prediction
  85. C
  86. C
  87. C SORTIES
  88. C
  89. C
  90. C IROF : MCHAML defini sur le MELEME de pointeur IFACEL,
  91. C contenant la masse volumique RHO
  92. C
  93. C IVITF : MCHAML defini sur le MELEME de pointeur IFACEL,
  94. C contenant la vitesse UN, UT dans le repaire local
  95. C (n,t) et defini sur le MELEME de pointeur IFACE,
  96. C contenant les cosinus directeurs du repere local
  97. C
  98. C IPF : MCHAML defini sur le MELEME de pointeur IFACEL,
  99. C contenant la pression P
  100. C
  101. C IGAMF : MCHAML defini sur le MELEME de pointeur IFACEL,
  102. C contenant le "gamma" du gaz
  103. C
  104. C LOGAN : anomalie detectee
  105. C
  106. C LOGNEG : (LOGICAL): si .TRUE. une pression ou une densité
  107. C negative a été detectée -> en interactif le
  108. C programme s'arrete en GIBIANE
  109. C (erreur stocké en MESERR et VALER)
  110. C
  111. C LOGBOR : (LOGICAL): si .TRUE. un gamma a ete detecte
  112. C dehor 1 et 3 (sa valeur stockée en MESERR et VALER;
  113. C en VAL1 et en VAL2 on stocke 1.0 et 3.0)
  114. C
  115. C MESERR
  116. C VALER
  117. C VAL1,
  118. C VAL2 : pour les messages d'erreur
  119. C
  120. C************************************************************************
  121. C
  122. C HISTORIQUE (Anomalies et modifications éventuelles)
  123. C
  124. C HISTORIQUE : Créée le 11.6.98.
  125. C
  126. C************************************************************************
  127. C
  128. C
  129. C ATTENTION: Cet programme marche si le MAILLAGE est convex;
  130. C si non il faut changer l'argoritme de calcul de
  131. C l'orientation des normales aux faces.
  132. C
  133. C La positivité n'est pas controlle parce que c'est déjà fait
  134. C dans l'operateur PRIM
  135. C
  136. C
  137. C************************************************************************
  138. C
  139. C**** Les variables
  140. C
  141. IMPLICIT INTEGER(I-N)
  142. INTEGER ICEN, IFACE, IFACEL, INORM,IROC, IGRROC, IALROC
  143. & , IVITC, IGRVC, IALVC
  144. & , IPC ,IGRPC, IALPC
  145. & , IGAMC
  146. & , IROF, IVITF, IPF, IGAMF
  147. & , IGEOM, NFAC, NCEN
  148. & , N1PTEL, N1EL, N2PTEL, N2EL, N2, N1, N3, L1, NLCE
  149. & , NLCF, NGCF, NGCEG, NLCEG, NGCED, NLCED, NGCF1
  150. & , IDIMP1, INDCEL
  151. REAL*8 VALER, VAL1, VAL2, XG, YG, ZG, XC, YC, ZC, XD, YD, ZD,
  152. & DELTAT
  153. & ,DXG, DYG, DZG, DXD, DYD, DZD,ORIENT
  154. & , CNX, CNY, CNZ, CTX, CTY, CTZ, CVX, CVY, CVZ
  155. & , ROG, PG, GAMG, UXG, UYG, UZG, UNG, UTG, UVG
  156. & , ROD, PD, GAMD, UXD, UYD, UZD, UND, UTD, UVD
  157. & , VALCEL, DCEL, ALCEL
  158. & , DROX, DROY, DROZ, DUXX, DUXY, DUXZ, DUYX, DUYY,
  159. & DUYZ, DUZX, DUZY, DUZZ, DPX, DPY, DPZ
  160. & , DRO, DUX, DUY, DUZ, DP, ALPHA(5)
  161. CHARACTER*(40) MESERR
  162. CHARACTER*(8) TYPE
  163. LOGICAL LOGAN,LOGNEG, LOGBOR, LOGTEM, LOGI1, LOGALP
  164. C
  165. C**** LOGALP = .TRUE. -> Prediction avec limiteur
  166. C
  167. PARAMETER(LOGALP = .TRUE.)
  168. C
  169. C**** Les Includes
  170. C
  171. -INC SMCOORD
  172.  
  173. -INC PPARAM
  174. -INC CCOPTIO
  175. -INC SMCHPOI
  176. POINTEUR MPROC.MPOVAL, MPGRR.MPOVAL,
  177. & MPVITC.MPOVAL, MPGRV.MPOVAL,
  178. & MPPC.MPOVAL, MPGRP.MPOVAL,
  179. & MPGAMC.MPOVAL, MPNORM.MPOVAL,
  180. & MPROP.MPOVAL, MPPP.MPOVAL, MPVITP.MPOVAL
  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. & MELGAM.MELVAL
  188. POINTEUR MELALR.MPOVAL,
  189. & MELALV.MPOVAL,
  190. & MELALP.MPOVAL
  191. -INC SMLENTI
  192. -INC SMELEME
  193. C
  194. C
  195. C**** Initialisation des parametres d'erreur déjà faite, i.e.
  196. C
  197. C LOGNEG = .FALSE.
  198. C LOGBOR = .FALSE.
  199. C MESERR = ' '
  200. C MOTERR(1:40) = MESERR(1:40)
  201. C VALER = 0.0D0
  202. C VAL1 = 0.0D0
  203. C VAL2 = 0.0D0
  204. C
  205. C
  206. C**** KRIPAD pour la correspondance global/local de centre
  207. C
  208. CALL KRIPAD(ICEN,MLENT1)
  209. C
  210. C**** MLENTI1 a MCORD.XCOORD(/1)/(IDIM+1) elements
  211. C
  212. C Si i est le numero global d'un noeud de ICEN,
  213. C MLENT1.LECT(i) contient sa position, i.e.
  214. C
  215. C I = numero global du noeud centre
  216. C MLENT1.LECT(i) = numero local du noeud centre
  217. C
  218. C MLENT1 déjà activé, i.e.
  219. C
  220. C SEGACT MLENT1
  221. C
  222. C**** Activation de CHPOINTs
  223. C
  224. C densité + grad + limiteur
  225. C vitesse + grad + limiteur
  226. C pression + grad + limiteur
  227. C gamma
  228. C cosinus directeurs des normales aux surface
  229. C
  230. CALL LICHT(IROC , MPROC , TYPE, IGEOM)
  231. CALL LICHT(IGRROC, MPGRR , TYPE, IGEOM)
  232. CALL LICHT(IVITC, MPVITC , TYPE, IGEOM)
  233. CALL LICHT(IGRVC, MPGRV , TYPE, IGEOM)
  234. CALL LICHT(IPC , MPPC , TYPE, IGEOM)
  235. CALL LICHT(IGRPC, MPGRP , TYPE, IGEOM)
  236. CALL LICHT(IGAMC, MPGAMC , TYPE, IGEOM)
  237. CALL LICHT(INORM, MPNORM , TYPE, IGEOM)
  238. C
  239. C**** Les MPOVALs 'Prediction'
  240. C
  241. IF(LOGTEM)THEN
  242. SEGINI, MPROP = MPROC
  243. SEGINI, MPPP = MPPC
  244. SEGINI, MPVITP = MPVITC
  245. ELSE
  246. MPROP = MPROC
  247. MPPP = MPPC
  248. MPVITP = MPVITC
  249. ENDIF
  250. C
  251. C**** Les Limiteurs
  252. C
  253. CALL LICHT(IALROC, MELALR , TYPE, IGEOM)
  254. CALL LICHT(IALVC, MELALV , TYPE, IGEOM)
  255. CALL LICHT(IALPC, MELALP , TYPE, IGEOM)
  256. C
  257. C
  258. C**** MPOVAL sont déjà activés i.e.:
  259. C
  260. C SEGACT MPROC
  261. C SEGACT MPGRR
  262. C SEGACT MPIALR
  263. C SEGACT MPVITC
  264. C SEGACT MPGRV
  265. C SEGACT MPIALV
  266. C SEGACT MPPC
  267. C SEGACT MPGRP
  268. C SEGACT MPIALP
  269. C SEGACT MPGAMC
  270. C SEGACT MPNORM
  271. C
  272. C**** Le MELEME FACEL
  273. C
  274. IPT1 = IFACEL
  275. IPT2 = IFACE
  276. SEGACT IPT1
  277. SEGACT IPT2
  278. NFAC = IPT1.NUM(/2)
  279. C
  280. C**** Creation de MCHAMLs contenant les etats gauche et droite,
  281. C
  282. C i.e.:
  283. C
  284. C vitesse + cosinus directors du repere local
  285. C densité
  286. C pression
  287. C gamma
  288. C
  289. C**** Cosinus directors du repere local et vitesse
  290. C
  291. C Les cosinus directeurs
  292. C
  293. N1 = 2
  294. N3 = 6
  295. L1 = 28
  296. SEGINI MCHEL1
  297. IVITF = MCHEL1
  298. MCHEL1.TITCHE = 'U '
  299. MCHEL1.IMACHE(1) = IFACE
  300. MCHEL1.IMACHE(2) = IFACEL
  301. MCHEL1.CONCHE(1) = '(n,t,v)in(x,y,z)'
  302. MCHEL1.CONCHE(2) = ' U in (n,t,v) '
  303. C
  304. C**** Valeurs des cosinus definies par respect au repair global, i.e.
  305. C
  306. MCHEL1.INFCHE(1,1) = 2
  307. MCHEL1.INFCHE(1,3) = NIFOUR
  308. MCHEL1.INFCHE(1,4) = 0
  309. MCHEL1.INFCHE(1,5) = 0
  310. MCHEL1.INFCHE(1,6) = 1
  311. MCHEL1.IFOCHE = IFOUR
  312. C
  313. C**** Valeurs de vitesse definies par respect au repair local, i.e.
  314. C
  315. MCHEL1.INFCHE(2,1) = 1
  316. MCHEL1.INFCHE(2,3) = NIFOUR
  317. MCHEL1.INFCHE(2,4) = 0
  318. MCHEL1.INFCHE(2,5) = 0
  319. MCHEL1.INFCHE(2,6) = 1
  320. C
  321. C**** Le cosinus directeurs
  322. C
  323. N1PTEL = 1
  324. N1EL = NFAC
  325. N2PTEL = 0
  326. N2EL = 0
  327. C
  328. C**** MCHAML a N2 composantes:
  329. C
  330. C cosinus directeurs du repere local (n,t1,t2)
  331. C
  332. C IDIM = 3 -> 9 composantes
  333. C
  334. N2 = 9
  335. SEGINI MCHAM1
  336. MCHEL1.ICHAML(1) = MCHAM1
  337. MCHAM1.NOMCHE(1) = 'NX '
  338. MCHAM1.NOMCHE(2) = 'NY '
  339. MCHAM1.NOMCHE(3) = 'NZ '
  340. MCHAM1.NOMCHE(4) = 'TX '
  341. MCHAM1.NOMCHE(5) = 'TY '
  342. MCHAM1.NOMCHE(6) = 'TZ '
  343. MCHAM1.NOMCHE(7) = 'VX '
  344. MCHAM1.NOMCHE(8) = 'VY '
  345. MCHAM1.NOMCHE(9) = 'VZ '
  346. MCHAM1.TYPCHE(1) = 'REAL*8 '
  347. MCHAM1.TYPCHE(2) = 'REAL*8 '
  348. MCHAM1.TYPCHE(3) = 'REAL*8 '
  349. MCHAM1.TYPCHE(4) = 'REAL*8 '
  350. MCHAM1.TYPCHE(5) = 'REAL*8 '
  351. MCHAM1.TYPCHE(6) = 'REAL*8 '
  352. MCHAM1.TYPCHE(7) = 'REAL*8 '
  353. MCHAM1.TYPCHE(8) = 'REAL*8 '
  354. MCHAM1.TYPCHE(9) = 'REAL*8 '
  355. SEGINI MELVNX
  356. SEGINI MELVNY
  357. SEGINI MELVNZ
  358. SEGINI MELT1X
  359. SEGINI MELT1Y
  360. SEGINI MELT1Z
  361. SEGINI MELT2X
  362. SEGINI MELT2Y
  363. SEGINI MELT2Z
  364. MCHAM1.IELVAL(1) = MELVNX
  365. MCHAM1.IELVAL(2) = MELVNY
  366. MCHAM1.IELVAL(3) = MELVNZ
  367. MCHAM1.IELVAL(4) = MELT1X
  368. MCHAM1.IELVAL(5) = MELT1Y
  369. MCHAM1.IELVAL(6) = MELT1Z
  370. MCHAM1.IELVAL(7) = MELT2X
  371. MCHAM1.IELVAL(8) = MELT2Y
  372. MCHAM1.IELVAL(9) = MELT2Z
  373. SEGDES MCHAM1
  374. C
  375. C**** Vitesse
  376. C
  377. N1EL = NFAC
  378. N1PTEL = 3
  379. N2EL = 0
  380. N2PTEL = 0
  381. C
  382. C**** MCHAML a N2 composantes:
  383. C
  384. C IDIM = 3 -> 3 composantes
  385. C
  386. N2 = 3
  387. SEGINI MCHAM1
  388. MCHEL1.ICHAML(2) = MCHAM1
  389. SEGDES MCHEL1
  390. MCHAM1.NOMCHE(1) = 'UN '
  391. MCHAM1.NOMCHE(2) = 'UT '
  392. MCHAM1.NOMCHE(3) = 'UV '
  393. MCHAM1.TYPCHE(1) = 'REAL*8 '
  394. MCHAM1.TYPCHE(2) = 'REAL*8 '
  395. MCHAM1.TYPCHE(3) = 'REAL*8 '
  396. SEGINI MELVUN
  397. SEGINI MELVUT
  398. SEGINI MELVUV
  399. MCHAM1.IELVAL(1) = MELVUN
  400. MCHAM1.IELVAL(2) = MELVUT
  401. MCHAM1.IELVAL(3) = MELVUV
  402. SEGDES MCHAM1
  403. C
  404. C**** Densite
  405. C
  406. N1 = 1
  407. N3 = 6
  408. L1 = 15
  409. SEGINI MCHEL2
  410. IROF = MCHEL2
  411. MCHEL2.IMACHE(1) = IFACEL
  412. MCHEL2.TITCHE = 'RO '
  413. MCHEL2.CONCHE(1) = ' '
  414. C
  415. C**** Valeurs independente du repére, i.e.
  416. C
  417. MCHEL2.INFCHE(1,1) = 0
  418. MCHEL2.INFCHE(1,3) = NIFOUR
  419. MCHEL2.INFCHE(1,4) = 0
  420. MCHEL2.INFCHE(1,5) = 0
  421. MCHEL2.INFCHE(1,6) = 1
  422. MCHEL2.IFOCHE = IFOUR
  423. N2 = 1
  424. SEGINI MCHAM1
  425. MCHEL2.ICHAML(1) = MCHAM1
  426. SEGDES MCHEL2
  427. MCHAM1.NOMCHE(1) = 'SCAL '
  428. MCHAM1.TYPCHE(1) = 'REAL*8 '
  429. SEGINI MELRO
  430. MCHAM1.IELVAL(1) = MELRO
  431. SEGDES MCHAM1
  432. C
  433. C**** Pression
  434. C
  435. MCHEL1 = IROF
  436. SEGINI, MCHEL2 = MCHEL1
  437. IPF = MCHEL2
  438. MCHEL2.TITCHE = 'P '
  439. C
  440. C**** MCHAM1 = MCHAML de la densite
  441. C
  442. SEGINI, MCHAM2 = MCHAM1
  443. MCHEL2.ICHAML(1) = MCHAM2
  444. SEGDES MCHEL2
  445. SEGINI MELP
  446. MCHAM2.IELVAL(1) = MELP
  447. SEGDES MCHAM2
  448. C
  449. C**** Gamma
  450. C
  451. MCHEL1 = IROF
  452. SEGINI, MCHEL2 = MCHEL1
  453. IGAMF = MCHEL2
  454. MCHEL2.TITCHE = 'GAMMA '
  455. C
  456. C**** MCHAM1 = MCHAML de la densite
  457. C
  458. SEGINI, MCHAM2 = MCHAM1
  459. MCHEL2.ICHAML(1) = MCHAM2
  460. SEGDES MCHEL2
  461. SEGINI MELGAM
  462. MCHAM2.IELVAL(1) = MELGAM
  463. SEGDES MCHAM2
  464. C
  465. C**** Recapitulatif: les MELVALs et les MPOVALs actives
  466. C
  467. C MELVNX, MELVNY, MELVNZ
  468. C MELT1X, MELT1Y, MELT1Z
  469. C MELT2X, MELT2Y, MELT2Z
  470. C
  471. C MELVUN, MELVUT -> vitesse
  472. C
  473. C MELRO -> densite
  474. C
  475. C MELP -> pression
  476. C
  477. C MELGAM -> gamma
  478. C
  479. C MPROC, MPROP -> densite
  480. C
  481. C MPVITC, MPVITP -> vitesse
  482. C
  483. C MPPC, MPPP -> pression
  484. C
  485. C MPGAMC -> gamma
  486. C
  487. C MPNORM -> normales aux faces
  488. C
  489. C
  490. C***********************************************************************
  491. C********* PREDICTION **************************************************
  492. C***********************************************************************
  493. C
  494. IF(LOGTEM)THEN
  495. C
  496. IPT3 = ICEN
  497. SEGACT IPT3
  498. NCEN = IPT3.NUM(/2)
  499. DO NLCE = 1, NCEN
  500. IF(LOGALP)THEN
  501. C
  502. C************* Gradients limités
  503. C
  504. ALPHA(1) = MELALR.VPOCHA(NLCE,1)
  505. ALPHA(2) = MELALV.VPOCHA(NLCE,1)
  506. ALPHA(3) = MELALV.VPOCHA(NLCE,2)
  507. ALPHA(4) = MELALV.VPOCHA(NLCE,3)
  508. ALPHA(5) = MELALP.VPOCHA(NLCE,1)
  509. ELSE
  510. C
  511. C************* Gradients non limités
  512. C
  513. ALPHA(1) = 1.0D0
  514. ALPHA(2) = 1.0D0
  515. ALPHA(3) = 1.0D0
  516. ALPHA(4) = 1.0D0
  517. ALPHA(5) = 1.0D0
  518. ENDIF
  519. ROG = MPROP.VPOCHA(NLCE,1)
  520. UXG = MPVITP.VPOCHA(NLCE,1)
  521. UYG = MPVITP.VPOCHA(NLCE,2)
  522. UZG = MPVITP.VPOCHA(NLCE,3)
  523. PG = MPPP.VPOCHA(NLCE,1)
  524. DROX = MPGRR.VPOCHA(NLCE,1)*ALPHA(1)
  525. DROY = MPGRR.VPOCHA(NLCE,2)*ALPHA(1)
  526. DROZ = MPGRR.VPOCHA(NLCE,3)*ALPHA(1)
  527. DUXX = MPGRV.VPOCHA(NLCE,1)*ALPHA(2)
  528. DUXY = MPGRV.VPOCHA(NLCE,2)*ALPHA(2)
  529. DUXZ = MPGRV.VPOCHA(NLCE,3)*ALPHA(2)
  530. DUYX = MPGRV.VPOCHA(NLCE,4)*ALPHA(3)
  531. DUYY = MPGRV.VPOCHA(NLCE,5)*ALPHA(3)
  532. DUYZ = MPGRV.VPOCHA(NLCE,6)*ALPHA(3)
  533. DUZX = MPGRV.VPOCHA(NLCE,7)*ALPHA(4)
  534. DUZY = MPGRV.VPOCHA(NLCE,8)*ALPHA(4)
  535. DUZZ = MPGRV.VPOCHA(NLCE,9)*ALPHA(4)
  536. DPX = MPGRP.VPOCHA(NLCE,1)*ALPHA(5)
  537. DPY = MPGRP.VPOCHA(NLCE,2)*ALPHA(5)
  538. DPZ = MPGRP.VPOCHA(NLCE,3)*ALPHA(5)
  539. GAMG = MPGAMC.VPOCHA(NLCE,1)
  540. DRO = UXG * DROX + ROG * ( DUXX + DUYY + DUZZ)
  541. & + UYG * DROY + UZG * DROZ
  542. DUX = UXG * DUXX + DPX / ROG + UYG * DUXY
  543. & + UZG * DUXZ
  544. DUY = UXG * DUYX + UYG * DUYY + DPY / ROG
  545. & + UZG * DUYZ
  546. DUZ = UXG * DUZX + UYG * DUZY + UZG * DUZZ
  547. & + DPZ / ROG
  548. DP = GAMG * PG * (DUXX + DUYY + DUZZ)
  549. & + UXG * DPX + UYG * DPY + UZG * DPZ
  550. C
  551. C************* Prediction avec/sans gradients limités
  552. C
  553. MPROP.VPOCHA(NLCE,1) = ROG - DELTAT * DRO
  554. MPVITP.VPOCHA(NLCE,1) = UXG - DELTAT * DUX
  555. MPVITP.VPOCHA(NLCE,2) = UYG - DELTAT * DUY
  556. MPVITP.VPOCHA(NLCE,3) = UZG - DELTAT * DUZ
  557. MPPP.VPOCHA(NLCE,1) = PG - DELTAT * DP
  558. ENDDO
  559. C
  560. ENDIF
  561. C
  562. C
  563. C***********************************************************************
  564. C********* CORRECTION **************************************************
  565. C***********************************************************************
  566. C
  567. C**** Boucle sur le faces
  568. C
  569. IDIMP1 = IDIM + 1
  570. DO NLCF = 1, NFAC
  571. C
  572. C******* NLCF = numero local du centre de face
  573. C NGCF = numero global du centre de face
  574. C NGCEG = numero global du centre ELT "gauche"
  575. C NLCEG = numero local du centre ELT "gauche"
  576. C NGCED = numero global du centre ELT "droite"
  577. C NLCED = numero local du centre ELT "droite"
  578. C
  579. NGCEG = IPT1.NUM(1,NLCF)
  580. NGCF = IPT1.NUM(2,NLCF)
  581. NGCED = IPT1.NUM(3,NLCF)
  582. NLCEG = MLENT1.LECT(NGCEG)
  583. NLCED = MLENT1.LECT(NGCED)
  584. C
  585. C******* TEST: IPT2.NUM(1,NLCF) = IPT1.NUM(2,NLCF)
  586. C
  587. NGCF1 = IPT2.NUM(1,NLCF)
  588. IF( NGCF1 .NE. NGCF) THEN
  589. LOGAN = .TRUE.
  590. MESERR(1:40) = 'PRET, subroutine pre122.eso '
  591. GOTO 9999
  592. ENDIF
  593. C
  594. C******* Cosinus directeurs des NORMALES aux faces
  595. C
  596. C On impose que les normales sont direct "Gauche" -> "Centre"
  597. C
  598. INDCEL = (NGCEG-1)*IDIMP1
  599. XG = XCOOR(INDCEL+1)
  600. YG = XCOOR(INDCEL+2)
  601. ZG = XCOOR(INDCEL+3)
  602. INDCEL = (NGCF-1)*IDIMP1
  603. XC = XCOOR(INDCEL + 1)
  604. YC = XCOOR(INDCEL + 2)
  605. ZC = XCOOR(INDCEL+3)
  606. INDCEL = (NGCED-1)*IDIMP1
  607. XD = XCOOR(INDCEL+1)
  608. YD = XCOOR(INDCEL+2)
  609. ZD = XCOOR(INDCEL+3)
  610. DXG = XC - XG
  611. DYG = YC - YG
  612. DZG = ZC - ZG
  613. DXD = XC - XD
  614. DYD = YC - YD
  615. DZD = ZC - ZD
  616. C
  617. C******* On calcule le sign du pruduit scalare
  618. C (Normales de Castem) * (vecteur "gauche" -> "centre")
  619. C
  620. CNX = MPNORM.VPOCHA(NLCF,7)
  621. CNY = MPNORM.VPOCHA(NLCF,8)
  622. CNZ = MPNORM.VPOCHA(NLCF,9)
  623. ORIENT = CNX * DXG + CNY * DYG + CNZ * DZG
  624. ORIENT = SIGN(1.0D0,ORIENT)
  625. IF(ORIENT .NE. 1.0D0)THEN
  626. LOGAN = .TRUE.
  627. MESERR(1:30)=
  628. & 'PRET , subroutine pre122.eso. '
  629. GOTO 9999
  630. ENDIF
  631. CNX = CNX * ORIENT
  632. CNY = CNY * ORIENT
  633. CNZ = CNZ * ORIENT
  634. C
  635. C********** Cosinus directeurs de tangente 1
  636. C
  637. CTX = MPNORM.VPOCHA(NLCF,1) * ORIENT
  638. CTY = MPNORM.VPOCHA(NLCF,2) * ORIENT
  639. CTZ = MPNORM.VPOCHA(NLCF,3) * ORIENT
  640. C
  641. C********** Cosinus directeurs de tangente 2
  642. C
  643. CVX = MPNORM.VPOCHA(NLCF,4) * ORIENT
  644. CVY = MPNORM.VPOCHA(NLCF,5) * ORIENT
  645. CVZ = MPNORM.VPOCHA(NLCF,6) * ORIENT
  646. C
  647. C
  648. C******* Les autres MELVALs
  649. C
  650. C
  651. C******* N.B.: On suppose qu'on a déjà controlle RO, P > 0
  652. C GAMMA \in (1,3)
  653. C Si non il faut le faire, en utilisant LOGBOR,
  654. C LOGNEG, VALER, VAL1, VAL2
  655. C
  656. C
  657. C
  658. C******* NGCEG = NGCED -> Mur
  659. C
  660. IF( NGCEG .EQ. NGCED)THEN
  661. C
  662. C********** Sur le mur on fait de reconstruction sur l'etat gauche
  663. C
  664. C
  665. C********** Etat gauche
  666. C
  667. VALCEL = MPROP.VPOCHA(NLCEG, 1)
  668. ALCEL = MELALR.VPOCHA(NLCEG, 1)
  669. DCEL = MPGRR.VPOCHA(NLCEG, 1)*DXG +
  670. & MPGRR.VPOCHA(NLCEG, 2)*DYG +
  671. & MPGRR.VPOCHA(NLCEG, 3)*DZG
  672. ROG = VALCEL + ALCEL * DCEL
  673. C
  674. VALCEL = MPPP.VPOCHA(NLCEG, 1)
  675. ALCEL = MELALP.VPOCHA(NLCEG, 1)
  676. DCEL = MPGRP.VPOCHA(NLCEG, 1)*DXG +
  677. & MPGRP.VPOCHA(NLCEG, 2)*DYG +
  678. & MPGRP.VPOCHA(NLCEG, 3)*DZG
  679. PG = VALCEL + ALCEL * DCEL
  680. C
  681. GAMG = MPGAMC.VPOCHA(NLCEG, 1)
  682. C
  683. VALCEL = MPVITP.VPOCHA(NLCEG, 1)
  684. ALCEL = MELALV.VPOCHA(NLCEG, 1)
  685. DCEL = MPGRV.VPOCHA(NLCEG, 1)*DXG +
  686. & MPGRV.VPOCHA(NLCEG, 2)*DYG +
  687. & MPGRV.VPOCHA(NLCEG, 3)*DZG
  688. UXG = VALCEL + ALCEL * DCEL
  689. C
  690. VALCEL = MPVITP.VPOCHA(NLCEG, 2)
  691. ALCEL = MELALV.VPOCHA(NLCEG, 2)
  692. DCEL = MPGRV.VPOCHA(NLCEG, 4)*DXG +
  693. & MPGRV.VPOCHA(NLCEG, 5)*DYG +
  694. & MPGRV.VPOCHA(NLCEG, 6)*DZG
  695. UYG = VALCEL + ALCEL * DCEL
  696. C
  697. VALCEL = MPVITP.VPOCHA(NLCEG, 3)
  698. ALCEL = MELALV.VPOCHA(NLCEG, 3)
  699. DCEL = MPGRV.VPOCHA(NLCEG, 7)*DXG +
  700. & MPGRV.VPOCHA(NLCEG, 8)*DYG +
  701. & MPGRV.VPOCHA(NLCEG, 9)*DZG
  702. UZG = VALCEL + ALCEL * DCEL
  703. C
  704. UNG = UXG * CNX + UYG * CNY + UZG * CNZ
  705. UTG = UXG * CTX + UYG * CTY + UZG * CTZ
  706. UVG = UXG * CVX + UYG * CVY + UZG * CVZ
  707. C
  708. C********** Si l'on fait pas de prediction, ce n'est pas necessaire de
  709. C controller la positivite' de la pression et de la densité; elle
  710. C est déjà garantie par la proprieté LED de limiteur; mais, vu
  711. C que le limiteur n'est pas calculé ici, mais dans un autre
  712. C operateur, on le fait
  713. C
  714. LOGI1 = (PG .LT. 0.0D0) .OR. (ROG .LT. 0.0D0)
  715. C
  716. IF(LOGI1)THEN
  717. C
  718. C************* Premier ordre en espace local
  719. C
  720. ROG = MPROC.VPOCHA(NLCEG,1)
  721. PG = MPPC.VPOCHA(NLCEG,1)
  722. UNG = MPVITC.VPOCHA(NLCEG,1)*CNX +
  723. & MPVITC.VPOCHA(NLCEG,2)*CNY +
  724. & MPVITC.VPOCHA(NLCEG,3)*CNZ
  725. UTG = MPVITC.VPOCHA(NLCEG,1)*CTX +
  726. & MPVITC.VPOCHA(NLCEG,2)*CTY +
  727. & MPVITC.VPOCHA(NLCEG,3)*CTZ
  728. UVG = MPVITC.VPOCHA(NLCEG,1)*CVX +
  729. & MPVITC.VPOCHA(NLCEG,2)*CVY +
  730. & MPVITC.VPOCHA(NLCEG,3)*CVZ
  731. ENDIF
  732. C
  733. C********** Son etat droite
  734. C
  735. ROD = ROG
  736. PD = PG
  737. GAMD = GAMG
  738. UND = -1.0D0 * UNG
  739. UTD = UTG
  740. UVD = UVG
  741. C
  742. C************* Fin cas mur
  743. C
  744. ELSE
  745. C
  746. C************* Etat gauche
  747. C
  748. VALCEL = MPROP.VPOCHA(NLCEG, 1)
  749. ALCEL = MELALR.VPOCHA(NLCEG, 1)
  750. DCEL = MPGRR.VPOCHA(NLCEG, 1)*DXG +
  751. & MPGRR.VPOCHA(NLCEG, 2)*DYG +
  752. & MPGRR.VPOCHA(NLCEG, 3)*DZG
  753. ROG = VALCEL + ALCEL * DCEL
  754. C
  755. VALCEL = MPPP.VPOCHA(NLCEG, 1)
  756. ALCEL = MELALP.VPOCHA(NLCEG, 1)
  757. DCEL = MPGRP.VPOCHA(NLCEG, 1)*DXG +
  758. & MPGRP.VPOCHA(NLCEG, 2)*DYG +
  759. & MPGRP.VPOCHA(NLCEG, 3)*DZG
  760. PG = VALCEL + ALCEL * DCEL
  761. C
  762. GAMG = MPGAMC.VPOCHA(NLCEG, 1)
  763. C
  764. VALCEL = MPVITP.VPOCHA(NLCEG, 1)
  765. ALCEL = MELALV.VPOCHA(NLCEG, 1)
  766. DCEL = MPGRV.VPOCHA(NLCEG, 1)*DXG +
  767. & MPGRV.VPOCHA(NLCEG, 2)*DYG +
  768. & MPGRV.VPOCHA(NLCEG, 3)*DZG
  769. UXG = VALCEL + ALCEL * DCEL
  770. C
  771. VALCEL = MPVITP.VPOCHA(NLCEG, 2)
  772. ALCEL = MELALV.VPOCHA(NLCEG, 2)
  773. DCEL = MPGRV.VPOCHA(NLCEG, 4)*DXG +
  774. & MPGRV.VPOCHA(NLCEG, 5)*DYG +
  775. & MPGRV.VPOCHA(NLCEG, 6)*DZG
  776. UYG = VALCEL + ALCEL * DCEL
  777. C
  778. VALCEL = MPVITP.VPOCHA(NLCEG, 3)
  779. ALCEL = MELALV.VPOCHA(NLCEG, 3)
  780. DCEL = MPGRV.VPOCHA(NLCEG, 7)*DXG +
  781. & MPGRV.VPOCHA(NLCEG, 8)*DYG +
  782. & MPGRV.VPOCHA(NLCEG, 9)*DZG
  783. UZG = VALCEL + ALCEL * DCEL
  784. C
  785. UNG = UXG * CNX + UYG * CNY + UZG * CNZ
  786. UTG = UXG * CTX + UYG * CTY + UZG * CTZ
  787. UVG = UXG * CVX + UYG * CVY + UZG * CVZ
  788. C
  789. C********** Positivite
  790. C
  791. LOGI1 = (PG .LT. 0.0D0) .OR. (ROG .LT. 0.0D0)
  792. C
  793. IF(LOGI1)THEN
  794. C
  795. C************* Premier ordre en espace local
  796. C
  797. ROG = MPROC.VPOCHA(NLCEG,1)
  798. PG = MPPC.VPOCHA(NLCEG,1)
  799. UNG = MPVITC.VPOCHA(NLCEG,1)*CNX +
  800. & MPVITC.VPOCHA(NLCEG,2)*CNY +
  801. & MPVITC.VPOCHA(NLCEG,3)*CNZ
  802. UTG = MPVITC.VPOCHA(NLCEG,1)*CTX +
  803. & MPVITC.VPOCHA(NLCEG,2)*CTY +
  804. & MPVITC.VPOCHA(NLCEG,3)*CTZ
  805. UVG = MPVITC.VPOCHA(NLCEG,1)*CVX +
  806. & MPVITC.VPOCHA(NLCEG,2)*CVY +
  807. & MPVITC.VPOCHA(NLCEG,3)*CVZ
  808. ENDIF
  809. C
  810. C********** Etat droite
  811. C
  812. VALCEL = MPROP.VPOCHA(NLCED, 1)
  813. ALCEL = MELALR.VPOCHA(NLCED, 1)
  814. DCEL = MPGRR.VPOCHA(NLCED, 1)*DXD +
  815. & MPGRR.VPOCHA(NLCED, 2)*DYD +
  816. & MPGRR.VPOCHA(NLCED, 3)*DZD
  817. ROD = VALCEL + ALCEL * DCEL
  818. C
  819. VALCEL = MPPP.VPOCHA(NLCED, 1)
  820. ALCEL = MELALP.VPOCHA(NLCED, 1)
  821. DCEL = MPGRP.VPOCHA(NLCED, 1)*DXD +
  822. & MPGRP.VPOCHA(NLCED, 2)*DYD +
  823. & MPGRP.VPOCHA(NLCED, 3)*DZD
  824. PD = VALCEL + ALCEL * DCEL
  825. C
  826. GAMD = MPGAMC.VPOCHA(NLCED, 1)
  827. C
  828. VALCEL = MPVITP.VPOCHA(NLCED, 1)
  829. ALCEL = MELALV.VPOCHA(NLCED, 1)
  830. DCEL = MPGRV.VPOCHA(NLCED, 1)*DXD +
  831. & MPGRV.VPOCHA(NLCED, 2)*DYD +
  832. & MPGRV.VPOCHA(NLCED, 3)*DZD
  833. UXD = VALCEL + ALCEL * DCEL
  834. C
  835. VALCEL = MPVITP.VPOCHA(NLCED, 2)
  836. ALCEL = MELALV.VPOCHA(NLCED, 2)
  837. DCEL = MPGRV.VPOCHA(NLCED, 4)*DXD +
  838. & MPGRV.VPOCHA(NLCED, 5)*DYD +
  839. & MPGRV.VPOCHA(NLCED, 6)*DZD
  840. UYD = VALCEL + ALCEL * DCEL
  841. C
  842. VALCEL = MPVITP.VPOCHA(NLCED, 3)
  843. ALCEL = MELALV.VPOCHA(NLCED, 3)
  844. DCEL = MPGRV.VPOCHA(NLCED, 7)*DXD +
  845. & MPGRV.VPOCHA(NLCED, 8)*DYD +
  846. & MPGRV.VPOCHA(NLCED, 9)*DZD
  847. UZD = VALCEL + ALCEL * DCEL
  848. C
  849. UND = UXD * CNX + UYD * CNY + UZD * CNZ
  850. UTD = UXD * CTX + UYD * CTY + UZD * CTZ
  851. UVD = UXD * CVX + UYD * CVY + UZD * CVZ
  852. C
  853. C********** Positivite
  854. C
  855. LOGI1 = (PD .LT. 0.0D0) .OR. (ROD .LT. 0.0D0)
  856. C
  857. IF(LOGI1)THEN
  858. C
  859. C************* Premier ordre en espace local
  860. C
  861. ROD = MPROC.VPOCHA(NLCED,1)
  862. PD = MPPC.VPOCHA(NLCED,1)
  863. UND = MPVITC.VPOCHA(NLCED,1)*CNX +
  864. & MPVITC.VPOCHA(NLCED,2)*CNY +
  865. & MPVITC.VPOCHA(NLCED,3)*CNZ
  866. UTD = MPVITC.VPOCHA(NLCED,1)*CTX +
  867. & MPVITC.VPOCHA(NLCED,2)*CTY +
  868. & MPVITC.VPOCHA(NLCED,3)*CTZ
  869. UVD = MPVITC.VPOCHA(NLCED,1)*CVX +
  870. & MPVITC.VPOCHA(NLCED,2)*CVY +
  871. & MPVITC.VPOCHA(NLCED,3)*CVZ
  872. ENDIF
  873. ENDIF
  874. C
  875. C******** Les MELVALs
  876. C
  877. C
  878. MELRO.VELCHE(1,NLCF) = ROG
  879. MELRO.VELCHE(3,NLCF) = ROD
  880. MELP.VELCHE(1,NLCF) = PG
  881. MELP.VELCHE(3,NLCF) = PD
  882. MELGAM.VELCHE(1,NLCF) = GAMG
  883. MELGAM.VELCHE(3,NLCF) = GAMD
  884. MELVUN.VELCHE(1,NLCF) = UNG
  885. MELVUN.VELCHE(3,NLCF) = UND
  886. MELVUT.VELCHE(1,NLCF) = UTG
  887. MELVUT.VELCHE(3,NLCF) = UTD
  888. MELVUV.VELCHE(1,NLCF) = UVG
  889. MELVUV.VELCHE(3,NLCF) = UVD
  890. MELVNX.VELCHE(1,NLCF) = CNX
  891. MELVNY.VELCHE(1,NLCF) = CNY
  892. MELVNZ.VELCHE(1,NLCF) = CNZ
  893. MELT1X.VELCHE(1,NLCF) = CTX
  894. MELT1Y.VELCHE(1,NLCF) = CTY
  895. MELT1Z.VELCHE(1,NLCF) = CTZ
  896. MELT2X.VELCHE(1,NLCF) = CVX
  897. MELT2Y.VELCHE(1,NLCF) = CVY
  898. MELT2Z.VELCHE(1,NLCF) = CVZ
  899. ENDDO
  900. C
  901. C**** Desactivation des SEGMENTs
  902. C
  903. SEGDES IPT1
  904. SEGDES IPT2
  905. C
  906. C**** Le MPOVALs 'Prediction' sont detruits (si existentes)
  907. C
  908. IF(LOGTEM)THEN
  909. SEGSUP MPROP
  910. SEGSUP MPVITP
  911. SEGSUP MPPP
  912. ENDIF
  913. C
  914. SEGDES MPROC
  915. SEGDES MPGRR
  916. SEGDES MELALR
  917. SEGDES MPVITC
  918. SEGDES MPGRV
  919. SEGDES MELALV
  920. SEGDES MPPC
  921. SEGDES MPGRP
  922. SEGDES MELALP
  923. SEGDES MPGAMC
  924. SEGDES MPNORM
  925. C
  926. SEGDES MELRO
  927. SEGDES MELP
  928. SEGDES MELGAM
  929. SEGDES MELVUN
  930. SEGDES MELVUT
  931. SEGDES MELVUV
  932. SEGDES MELVNX
  933. SEGDES MELVNY
  934. SEGDES MELVNZ
  935. SEGDES MELT1X
  936. SEGDES MELT1Y
  937. SEGDES MELT1Z
  938. SEGDES MELT2X
  939. SEGDES MELT2Y
  940. SEGDES MELT2Z
  941. C
  942. C
  943. C**** Destruction du MELNTI correspondance local/global
  944. C
  945. SEGSUP MLENT1
  946. C
  947. 9999 CONTINUE
  948. C
  949. RETURN
  950. END
  951.  
  952.  
  953.  
  954.  
  955.  
  956.  
  957.  
  958.  
  959.  
  960.  
  961.  
  962.  
  963.  

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