Télécharger pre711.eso

Retour à la liste

Numérotation des lignes :

  1. C PRE711 SOURCE BECC 11/05/18 21:15:27 6974
  2. SUBROUTINE PRE711(NESP,MLMESP,
  3. & ICEN,IFACE,IFACEL,INORM,
  4. & IPHI, IGRPHI, ILIPHI,
  5. & IRN1, IGRRN1, ILIRN1,
  6. & IVN1, IGRVN1, ILIVN1,
  7. & IPN1, IGRPN1, ILIPN1,
  8. & IYMA, IGRYMA, ILIYMA,
  9. & IALC, IGRALC, ILIALC,
  10. & IPHIF, IRN1F, IVN1F, IPN1F, IYF, IALF)
  11. C************************************************************************
  12. C
  13. C PROJET : CASTEM 2000
  14. C
  15. C NOM : PRE711
  16. C
  17. C DESCRIPTION : Voir PRE71
  18. C
  19. C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI)
  20. C
  21. C AUTEUR : A. BECCANTINI, DEN/DM2S/SFME/LTMF
  22. C
  23. C************************************************************************
  24. C
  25. C ENTREES
  26. C
  27. C NESP : number of the species involved in the EULER equations.
  28. C
  29. C MLMESP : MLMOTS object; names of the species involved in the
  30. C EULER equations.
  31. C
  32. C 1) Pointeurs de MELEMEs et de CHPOINTs de la table DOMAINE
  33. C
  34. C ICEN : MELEME de 'POI1' SPG des CENTRES
  35. C
  36. C IFACE : MELEME de 'POI1' SPG des FACES
  37. C
  38. C IFACEL : MELEME de 'SEG3' avec
  39. C CENTRE d'Elt "gauche"
  40. C CENTRE de Face
  41. C CENTRE d'Elt "droite"
  42. C
  43. C N.B. = IFACE.NUM(i,1) = IFACEL.NUM(i,2)
  44. C
  45. C INORM : CHPOINT des cosinus directeurs de normales aux faces
  46. C
  47. C 2) Autres pointeurs
  48. C
  49. C IPHI, IGRPHI, ILIPHI,
  50. C CHPOINT "CENTRE" de phi, gradient et limiteur
  51. C
  52. C IRN1, IGRRN1, ILIRN1,
  53. C CHPOINT "CENTRE" de densité, gradient et limiteur
  54. C
  55. C IVN1, IGRVN1, ILIVN1,
  56. C CHPOINT "CENTRE" de vitesse, gradient et limiteur
  57. C
  58. C IPN1, IGRPN1, ILIPN1.
  59. C CHPOINT "CENTRE" de pression, gradients et limiteurs
  60. C
  61. C IYMA, IGRYMA, ILIYMA,
  62. C CHPOINT "CENTRE" de Y, gradient et limiteur
  63. C
  64. C IALC, IGRALC, ILIALC,
  65. C CHPOINT "CENTRE" de ALPHA, gradient et limiteur
  66. C
  67. C SORTIES
  68. C
  69. C IPHIF, IRN1F, IVN1F, IPN1F, IYF, IALF
  70. C MCHAMLs definis sur le MELEME de pointeur IFACEL
  71. C
  72. C************************************************************************
  73. C
  74. C HISTORIQUE (Anomalies et modifications éventuelles)
  75. C
  76. C HISTORIQUE : Crée le 03.12.10.
  77. C
  78. C************************************************************************
  79. C
  80. C
  81. C ATTENTION: Cet programme marche si le MAILLAGE est convex;
  82. C si non il faut changer l'algoritme de calcul de
  83. C l'orientation des normales aux faces.
  84. C
  85. C La positivité n'est pas controlle parce que c'est déjà fait
  86. C dans l'operateur PRIM
  87. C
  88. C
  89. C************************************************************************
  90. C
  91. IMPLICIT INTEGER(I-N)
  92. C
  93. C**** Les variables
  94. C
  95. INTEGER NESP, I1
  96. & , ICEN, IFACE, IFACEL, INORM
  97. & , IPHI, IGRPHI, ILIPHI
  98. & , IRN1, IGRRN1, ILIRN1
  99. & , IVN1, IGRVN1, ILIVN1
  100. & , IPN1, IGRPN1, ILIPN1
  101. & , IYMA, IGRYMA, ILIYMA
  102. & , IALC, IGRALC, ILIALC
  103. & , IPHIF, IRN1F
  104. & , IVN1F, IPN1F
  105. & , IYF, IALF
  106. & , IGEOM, NFAC
  107. & , N1PTEL, N1EL, N2PTEL, N2EL, N2, N1, N3, L1
  108. & , NLCF, NGCF, NGCEG, NLCEG, NGCED, NLCED, NGCF1
  109. & , IDIMP1, INDCEL
  110.  
  111. REAL*8 XG, YG, XC, YC, XD, YD
  112. & ,DXG, DYG, DXD, DYD
  113. & , CNX, CNY, CTX, CTY, ORIENT
  114. & , PHIG, RN1G, PN1G
  115. & , UX1G, UY1G
  116. & , UN1G, UT1G
  117. & , PHID, RN1D, PN1D
  118. & , UX1D, UY1D
  119. & , UN1D, UT1D
  120. & , VALCEL, DCEL, LIMCEL
  121. CHARACTER*(40) MESERR
  122. CHARACTER*(8) TYPE
  123. CHARACTER*(4) CARCEL
  124. LOGICAL LOGI1
  125. C
  126. C**** Segments des fractions massiques gauche et droit
  127. C
  128. INTEGER NS
  129. SEGMENT FRAMAS
  130. REAL*8 FRAMG(NS), FRAMD(NS)
  131. ENDSEGMENT
  132. POINTEUR ALPHA.FRAMAS
  133. C
  134. C**** Les Includes
  135. C
  136. -INC SMCOORD
  137. -INC PPARAM
  138. -INC CCOPTIO
  139. -INC SMCHPOI
  140. POINTEUR MPPHI.MPOVAL, MPGPHI.MPOVAL, MPLPHI.MPOVAL
  141. POINTEUR MPRN1.MPOVAL, MPGRN1.MPOVAL, MPLRN1.MPOVAL
  142. POINTEUR MPVN1.MPOVAL, MPGVN1.MPOVAL, MPLVN1.MPOVAL
  143. POINTEUR MPPN1.MPOVAL, MPGPN1.MPOVAL, MPLPN1.MPOVAL
  144. POINTEUR MPYMA.MPOVAL, MPGYMA.MPOVAL, MPLYMA.MPOVAL
  145. POINTEUR MPALC.MPOVAL, MPGALC.MPOVAL, MPLALC.MPOVAL
  146. POINTEUR MPNORM.MPOVAL
  147. -INC SMCHAML
  148. C Melval des cosinus directeurs
  149. POINTEUR MELVNX.MELVAL, MELVNY.MELVAL,
  150. & MELT1X.MELVAL, MELT1Y.MELVAL
  151. C Melval des vitesses
  152. POINTEUR MEVUN1.MELVAL, MEVUT1.MELVAL
  153. C Melval des densités, pressions, alphas
  154. POINTEUR MELRN1.MELVAL
  155. POINTEUR MELPN1.MELVAL
  156. POINTEUR MELPHI.MELVAL
  157. POINTEUR MCHAMY.MCHAML
  158. POINTEUR MCHAMA.MCHAML
  159. -INC SMLENTI
  160. -INC SMELEME
  161. -INC SMLMOTS
  162. POINTEUR MLMESP.MLMOTS
  163. C
  164. LOGI1 = .FALSE.
  165. C
  166. C**** KRIPAD pour la correspondance global/local de centre
  167. C
  168. CALL KRIPAD(ICEN,MLENT1)
  169. C
  170. C**** MLENTI1 a MCORD.XCOORD(/1)/(IDIM+1) elements
  171. C
  172. C Si i est le numero global d'un noeud de ICEN,
  173. C MLENT1.LECT(i) contient sa position, i.e.
  174. C
  175. C I = numero global du noeud centre
  176. C MLENT1.LECT(i) = numero local du noeud centre
  177. C
  178. C MLENT1 déjà activé, i.e.
  179. C
  180. C SEGACT MLENT1
  181. C
  182. C**** Activation de CHPOINTs
  183. C
  184. C phi + grad + limiteur
  185. C densité + grad + limiteur
  186. C vitesse + grad + limiteur
  187. C pression + grad + limiteur
  188. C ymass + grad + limiteur
  189. C alpha + grad + limiteur
  190. C cosinus directeurs des normales aux surface
  191. C
  192. CALL LICHT(IPHI , MPPHI , TYPE, IGEOM)
  193. C SEGACT MPPHI
  194. CALL LICHT(IGRPHI , MPGPHI , TYPE, IGEOM)
  195. C SEGACT MPGPHI
  196. CALL LICHT(ILIPHI , MPLPHI , TYPE, IGEOM)
  197. C SEGACT MPLPHI
  198. CALL LICHT(IRN1 , MPRN1 , TYPE, IGEOM)
  199. C SEGACT MPRN1
  200. CALL LICHT(IGRRN1 , MPGRN1 , TYPE, IGEOM)
  201. C SEGACT MPGRN1
  202. CALL LICHT(ILIRN1 , MPLRN1 , TYPE, IGEOM)
  203. C SEGACT MPLRN1
  204. CALL LICHT(IVN1 , MPVN1 , TYPE, IGEOM)
  205. C SEGACT MPVN1
  206. CALL LICHT(IGRVN1 , MPGVN1 , TYPE, IGEOM)
  207. C SEGACT MPGVN1
  208. CALL LICHT(ILIVN1 , MPLVN1 , TYPE, IGEOM)
  209. C SEGACT MPLVN1
  210. CALL LICHT(IPN1 , MPPN1 , TYPE, IGEOM)
  211. C SEGACT MPPN1
  212. CALL LICHT(IGRPN1 , MPGPN1 , TYPE, IGEOM)
  213. C SEGACT MPGPN1
  214. CALL LICHT(ILIPN1 , MPLPN1 , TYPE, IGEOM)
  215. C SEGACT MPLPN1
  216. IF (NESP .GE. 1)THEN
  217. CALL LICHT(IYMA , MPYMA , TYPE, IGEOM)
  218. C SEGACT MPYMA
  219. CALL LICHT(IGRYMA , MPGYMA , TYPE, IGEOM)
  220. C SEGACT MPGYMA
  221. CALL LICHT(ILIYMA , MPLYMA , TYPE, IGEOM)
  222. C SEGACT MPLYMA
  223. CALL LICHT(IALC , MPALC , TYPE, IGEOM)
  224. C SEGACT MPALC
  225. CALL LICHT(IGRALC , MPGALC , TYPE, IGEOM)
  226. C SEGACT MPGALC
  227. CALL LICHT(ILIALC , MPLALC , TYPE, IGEOM)
  228. C SEGACT MPLALC
  229. ENDIF
  230. C
  231. C**** Le cosinus directeurs
  232. C
  233. CALL LICHT(INORM , MPNORM , TYPE, IGEOM)
  234. C SEGACT MPNORM
  235. C
  236. C**** Les MPOVAL sont déjà activés i.e.:
  237. C
  238. C
  239. C**** Le MELEME FACEL
  240. C
  241. IPT1 = IFACEL
  242. IPT2 = IFACE
  243. SEGACT IPT1
  244. SEGACT IPT2
  245. NFAC = IPT1.NUM(/2)
  246. C
  247. C**** Creation de MCHAMLs contenant les etats gauche et droite,
  248. C
  249. C i.e.:
  250. C
  251. C vitesse + cosinus directors du repere local
  252. C alpha
  253. C densité
  254. C pression
  255. C
  256. C**********************************************************
  257. C**** Cosinus directors du repere local et vitesse ********
  258. C**********************************************************
  259. C
  260. C Les cosinus directeurs
  261. C
  262. N1 = 2
  263. N3 = 6
  264. L1 = 28
  265. SEGINI MCHEL1
  266. IVN1F = MCHEL1
  267. MCHEL1.TITCHE = 'U '
  268. MCHEL1.IMACHE(1) = IFACE
  269. MCHEL1.IMACHE(2) = IFACEL
  270. MCHEL1.CONCHE(1) = ' (n,t) in (x,y) '
  271. MCHEL1.CONCHE(2) = ' U in (n,t) '
  272. * MCHEL1.NOPHAS(1) = ' '
  273. * MCHEL1.NOPHAS(2) = ' '
  274. C
  275. C**** Valeurs des cosinus definies par respect au repair global, i.e.
  276. C
  277. MCHEL1.INFCHE(1,1) = 2
  278. MCHEL1.INFCHE(1,3) = NIFOUR
  279. MCHEL1.INFCHE(1,4) = 0
  280. MCHEL1.INFCHE(1,5) = 0
  281. MCHEL1.INFCHE(1,6) = 0
  282. MCHEL1.IFOCHE = IFOUR
  283. C
  284. C**** Valeurs de vitesse definies par respect au repair local, i.e.
  285. C
  286. MCHEL1.INFCHE(2,1) = 1
  287. MCHEL1.INFCHE(2,3) = NIFOUR
  288. MCHEL1.INFCHE(2,4) = 0
  289. MCHEL1.INFCHE(2,5) = 0
  290. MCHEL1.INFCHE(2,6) = 0
  291. C
  292. C**** Le cosinus directeurs
  293. C
  294. N1PTEL = 1
  295. N1EL = NFAC
  296. N2PTEL = 0
  297. N2EL = 0
  298. C
  299. C**** MCHAML a N2 composantes:
  300. C
  301. C cosinus directeurs du repere local (n,t1)
  302. C
  303. C IDIM = 2 -> 4 composantes
  304. C
  305. N2 = 4
  306. SEGINI MCHAM1
  307. MCHEL1.ICHAML(1) = MCHAM1
  308. MCHAM1.NOMCHE(1) = 'NX '
  309. MCHAM1.NOMCHE(2) = 'NY '
  310. MCHAM1.NOMCHE(3) = 'TX '
  311. MCHAM1.NOMCHE(4) = 'TY '
  312. MCHAM1.TYPCHE(1) = 'REAL*8 '
  313. MCHAM1.TYPCHE(2) = 'REAL*8 '
  314. MCHAM1.TYPCHE(3) = 'REAL*8 '
  315. MCHAM1.TYPCHE(4) = 'REAL*8 '
  316. SEGINI MELVNX
  317. SEGINI MELVNY
  318. SEGINI MELT1X
  319. SEGINI MELT1Y
  320. MCHAM1.IELVAL(1) = MELVNX
  321. MCHAM1.IELVAL(2) = MELVNY
  322. MCHAM1.IELVAL(3) = MELT1X
  323. MCHAM1.IELVAL(4) = MELT1Y
  324. SEGDES MCHAM1
  325. C
  326. C**** Vitesse
  327. C
  328. N1EL = NFAC
  329. N1PTEL = 3
  330. N2EL = 0
  331. N2PTEL = 0
  332. C
  333. C**** MCHAML a N2 composantes:
  334. C
  335. C IDIM = 2 -> 2 composantes
  336. C
  337. N2 = 2
  338. SEGINI MCHAM1
  339. MCHEL1.ICHAML(2) = MCHAM1
  340. MCHAM1.NOMCHE(1) = 'UN '
  341. MCHAM1.NOMCHE(2) = 'UT '
  342. MCHAM1.TYPCHE(1) = 'REAL*8 '
  343. MCHAM1.TYPCHE(2) = 'REAL*8 '
  344. SEGINI MEVUN1
  345. SEGINI MEVUT1
  346. MCHAM1.IELVAL(1) = MEVUN1
  347. MCHAM1.IELVAL(2) = MEVUT1
  348. SEGDES MCHAM1
  349. C
  350. C**********************************************************
  351. C**** PHI1 ********
  352. C**********************************************************
  353. C
  354. C**** PHI1
  355. C
  356. N1 = 1
  357. N3 = 6
  358. L1 = 15
  359. SEGINI MCHEL2
  360. IPHIF = MCHEL2
  361. MCHEL2.IMACHE(1) = IFACEL
  362. MCHEL2.TITCHE = 'PHI '
  363. MCHEL2.CONCHE(1) = ' '
  364. * MCHEL2.NOPHAS(1) = ' '
  365. C
  366. C**** Valeurs independente du repére, i.e.
  367. C
  368. MCHEL2.INFCHE(1,1) = 0
  369. MCHEL2.INFCHE(1,3) = NIFOUR
  370. MCHEL2.INFCHE(1,4) = 0
  371. MCHEL2.INFCHE(1,5) = 0
  372. MCHEL2.INFCHE(1,6) = 0
  373. MCHEL2.IFOCHE = IFOUR
  374. N2 = 1
  375. SEGINI MCHAM1
  376. MCHEL2.ICHAML(1) = MCHAM1
  377. C We cannot deseactivate MCHEL2 = IPHIF now since we
  378. C use it after...
  379. MCHAM1.NOMCHE(1) = 'SCAL '
  380. MCHAM1.TYPCHE(1) = 'REAL*8 '
  381. SEGINI MELPHI
  382. MCHAM1.IELVAL(1) = MELPHI
  383. SEGDES MCHAM1
  384. C
  385. C**********************************************************
  386. C**** IRN1F and IRN2F ********
  387. C**********************************************************
  388. C
  389. MCHEL1 = IPHIF
  390. SEGINI, MCHEL2 = MCHEL1
  391. IRN1F = MCHEL2
  392. MCHEL2.TITCHE = 'RHO1 '
  393. MCHAM1 = MCHEL1.ICHAML(1)
  394. SEGINI, MCHAM2 = MCHAM1
  395. MCHEL2.ICHAML(1) = MCHAM2
  396. SEGDES MCHEL2
  397. SEGINI MELRN1
  398. MCHAM2.IELVAL(1) = MELRN1
  399. SEGDES MCHAM2
  400. C
  401. C
  402. C
  403. C**********************************************************
  404. C**** IPN1F
  405. C**********************************************************
  406. C
  407. MCHEL1 = IPHIF
  408. SEGINI, MCHEL2 = MCHEL1
  409. IPN1F = MCHEL2
  410. MCHEL2.TITCHE = 'P1 '
  411. MCHAM1 = MCHEL1.ICHAML(1)
  412. SEGINI, MCHAM2 = MCHAM1
  413. MCHEL2.ICHAML(1) = MCHAM2
  414. SEGDES MCHEL2
  415. SEGDES MCHEL1
  416. C We desactivate MCHEL1 = IPHIF now !
  417. SEGINI MELPN1
  418. MCHAM2.IELVAL(1) = MELPN1
  419. SEGDES MCHAM2
  420. C
  421. C write(*,*) 'Qui ci arrivo 1...'
  422. IF (NESP .GE. 1) THEN
  423. C
  424. SEGACT MLMESP
  425. C
  426. C******* YF
  427. C
  428. NS = NESP
  429. SEGINI FRAMAS
  430. MCHEL1 = IRN1F
  431. SEGINI, MCHEL2 = MCHEL1
  432. IYF = MCHEL2
  433. MCHEL2.TITCHE = 'Y '
  434. N2 = NESP
  435. SEGINI MCHAMY
  436. MCHEL2.ICHAML(1) = MCHAMY
  437. SEGDES MCHEL2
  438. N1EL = NFAC
  439. N1PTEL = 3
  440. N2EL = 0
  441. N2PTEL = 0
  442. DO I1 = 1, NESP
  443. SEGINI MELVA1
  444. MCHAMY.IELVAL(I1) = MELVA1
  445. C AB Error. It is corrected by the following line
  446. CARCEL = MLMESP.MOTS(I1)
  447. TYPE = ' '
  448. TYPE(1:4) = CARCEL(1:4)
  449. MCHAMY.NOMCHE(I1) = TYPE
  450. MCHAMY.TYPCHE(I1) = 'REAL*8 '
  451. ENDDO
  452. C
  453. C******* IALF
  454. C
  455. NS = NESP
  456. SEGINI ALPHA
  457. MCHEL1 = IRN1F
  458. SEGINI, MCHEL2 = MCHEL1
  459. IALF = MCHEL2
  460. MCHEL2.TITCHE = 'ALPHA '
  461. N2 = NESP
  462. SEGINI MCHAMA
  463. MCHEL2.ICHAML(1) = MCHAMA
  464. SEGDES MCHEL2
  465. N1EL = NFAC
  466. N1PTEL = 3
  467. N2EL = 0
  468. N2PTEL = 0
  469. DO I1 = 1, NESP
  470. SEGINI MELVA1
  471. MCHAMA.IELVAL(I1) = MELVA1
  472. C AB Error. It is corrected by the following line
  473. CARCEL = MLMESP.MOTS(I1)
  474. TYPE = ' '
  475. TYPE(1:4) = CARCEL(1:4)
  476. MCHAMA.NOMCHE(I1) = TYPE
  477. MCHAMA.TYPCHE(I1) = 'REAL*8 '
  478. ENDDO
  479. C
  480. SEGDES MLMESP
  481. ENDIF
  482. C
  483. C write(*,*) 'Qui ci arrivo 2...'
  484. C
  485. C
  486. C**********************************************************
  487. C**** Boucle sur le faces *********
  488. C**********************************************************
  489. C
  490. IDIMP1 = IDIM + 1
  491. DO NLCF = 1, NFAC
  492. C
  493. C******* NLCF = numero local du centre de face
  494. C NGCF = numero global du centre de face
  495. C NGCEG = numero global du centre ELT "gauche"
  496. C NLCEG = numero local du centre ELT "gauche"
  497. C NGCED = numero global du centre ELT "droite"
  498. C NLCED = numero local du centre ELT "droite"
  499. C
  500. NGCEG = IPT1.NUM(1,NLCF)
  501. NGCF = IPT1.NUM(2,NLCF)
  502. NGCED = IPT1.NUM(3,NLCF)
  503. NLCEG = MLENT1.LECT(NGCEG)
  504. NLCED = MLENT1.LECT(NGCED)
  505. C
  506. C******* TEST: IPT2.NUM(1,NLCF) = IPT1.NUM(2,NLCF)
  507. C
  508. NGCF1 = IPT2.NUM(1,NLCF)
  509. IF( NGCF1 .NE. NGCF) THEN
  510. MESERR(1:40) = 'PRET, subroutine pre611.eso '
  511. WRITE(IOIMP,*) MESERR
  512. CALL ERREUR(5)
  513. GOTO 9999
  514. ENDIF
  515. C
  516. C******* Cosinus directeurs des NORMALES aux faces
  517. C
  518. C On impose que les normales sont direct "Gauche" -> "Centre"
  519. C
  520. INDCEL = (NGCEG-1)*IDIMP1
  521. XG = XCOOR(INDCEL+1)
  522. YG = XCOOR(INDCEL+2)
  523. INDCEL = (NGCF-1)*IDIMP1
  524. XC = XCOOR(INDCEL + 1)
  525. YC = XCOOR(INDCEL + 2)
  526. INDCEL = (NGCED-1)*IDIMP1
  527. XD = XCOOR(INDCEL+1)
  528. YD = XCOOR(INDCEL+2)
  529. DXG = XC - XG
  530. DYG = YC - YG
  531. DXD = XC - XD
  532. DYD = YC - YD
  533. C
  534. C******* On calcule le sign du pruduit scalare
  535. C (Normales de Castem) * (vecteur "gauche" -> "centre")
  536. C
  537. CNX = MPNORM.VPOCHA(NLCF,1)
  538. CNY = MPNORM.VPOCHA(NLCF,2)
  539. ORIENT = CNX * DXG + CNY * DYG
  540. ORIENT = SIGN(1.0D0,ORIENT)
  541. IF(ORIENT .NE. 1.0D0)THEN
  542. MESERR(1:30)=
  543. & 'PRET , subroutine pre611.eso. '
  544. GOTO 9999
  545. ENDIF
  546. CNX = CNX * ORIENT
  547. CNY = CNY * ORIENT
  548. C
  549. C******* Cosinus directeurs de tangent 2D
  550. C
  551. CTX = -1.0D0 * CNY
  552. CTY = CNX
  553. C
  554. C******* Les autres MELVALs
  555. C
  556. C
  557. C******* N.B.: On suppose qu'on a déjà controlle RO, P > 0...
  558. C
  559. C
  560. C******* Etat gauche
  561. C
  562. C PHI
  563. C
  564. VALCEL = MPPHI.VPOCHA(NLCEG, 1)
  565. LIMCEL = MPLPHI.VPOCHA(NLCEG, 1)
  566. DCEL = (MPGPHI.VPOCHA(NLCEG, 1) * DXG) +
  567. & (MPGPHI.VPOCHA(NLCEG, 2) * DYG)
  568. PHIG = VALCEL + (LIMCEL * DCEL)
  569. C write(*,*) valcel, limcel, dcel
  570. C
  571. C
  572. C RN
  573. C
  574. VALCEL = MPRN1.VPOCHA(NLCEG, 1)
  575. LIMCEL = MPLRN1.VPOCHA(NLCEG, 1)
  576. DCEL = (MPGRN1.VPOCHA(NLCEG, 1) * DXG) +
  577. & (MPGRN1.VPOCHA(NLCEG, 2) * DYG)
  578. RN1G = VALCEL + (LIMCEL * DCEL)
  579. C
  580. C
  581. C PN
  582. C
  583. VALCEL = MPPN1.VPOCHA(NLCEG, 1)
  584. LIMCEL = MPLPN1.VPOCHA(NLCEG, 1)
  585. DCEL = (MPGPN1.VPOCHA(NLCEG, 1) * DXG) +
  586. & (MPGPN1.VPOCHA(NLCEG, 2) * DYG)
  587. PN1G = VALCEL + (LIMCEL * DCEL)
  588. C
  589. C VN
  590. C
  591. VALCEL = MPVN1.VPOCHA(NLCEG, 1)
  592. LIMCEL = MPLVN1.VPOCHA(NLCEG, 1)
  593. DCEL = MPGVN1.VPOCHA(NLCEG, 1)*DXG +
  594. & MPGVN1.VPOCHA(NLCEG, 2)*DYG
  595. UX1G = VALCEL + (LIMCEL * DCEL)
  596. C
  597. VALCEL = MPVN1.VPOCHA(NLCEG, 2)
  598. LIMCEL = MPLVN1.VPOCHA(NLCEG, 2)
  599. DCEL = MPGVN1.VPOCHA(NLCEG, 3)*DXG +
  600. & MPGVN1.VPOCHA(NLCEG, 4)*DYG
  601. UY1G = VALCEL + (LIMCEL * DCEL)
  602. C
  603. C YN
  604. C
  605. DO I1 = 1, NESP
  606. INDCEL = 2 * I1 - 1
  607. VALCEL = MPYMA.VPOCHA(NLCEG,I1)
  608. DCEL = MPGYMA.VPOCHA(NLCEG, INDCEL)*DXG +
  609. & MPGYMA.VPOCHA(NLCEG,INDCEL + 1 )*DYG
  610. LIMCEL = MPLYMA.VPOCHA(NLCEG,I1)
  611. FRAMAS.FRAMG(I1) = VALCEL + (LIMCEL * DCEL)
  612. ENDDO
  613. C
  614. C ALPHAN
  615. C
  616. DO I1 = 1, NESP
  617. INDCEL = 2 * I1 - 1
  618. VALCEL = MPALC.VPOCHA(NLCEG,I1)
  619. DCEL = MPGALC.VPOCHA(NLCEG, INDCEL)*DXG +
  620. & MPGALC.VPOCHA(NLCEG,INDCEL + 1 )*DYG
  621. LIMCEL = MPLALC.VPOCHA(NLCEG,I1)
  622. ALPHA.FRAMG(I1) = VALCEL + (LIMCEL * DCEL)
  623. ENDDO
  624. C
  625. C
  626. C******* Si l'on fait pas de prediction, ce n'est pas necessaire de
  627. C controller la positivite' de la pression et de la densité; elle
  628. C est déjà garantie par la proprieté LED de limiteur.
  629. C If we want to check it, just uncomment LOGI1.
  630. C
  631. C LOGI1 = (RN1G .LT. 0.0D0) .OR.
  632. C & (PN1G .LT. 0.0D0)
  633. C
  634. IF ( NGCEG .EQ. NGCED) THEN
  635. C
  636. C********** Cas mur
  637. C
  638. IF(LOGI1)THEN
  639. C
  640. C********** Premier ordre en espace local
  641. C
  642. PHIG = MPPHI.VPOCHA(NLCEG,1)
  643. RN1G = MPRN1.VPOCHA(NLCEG,1)
  644. PN1G = MPPN1.VPOCHA(NLCEG,1)
  645. UX1G = MPVN1.VPOCHA(NLCEG,1)
  646. UY1G = MPVN1.VPOCHA(NLCEG,2)
  647. DO I1 = 1, NESP
  648. FRAMAS.FRAMG(I1) = MPYMA.VPOCHA(NLCEG,I1)
  649. ENDDO
  650. DO I1 = 1, NESP
  651. ALPHA.FRAMG(I1) = MPALC.VPOCHA(NLCEG,I1)
  652. ENDDO
  653. ENDIF
  654. C
  655. UN1G = UX1G * CNX + UY1G * CNY
  656. UT1G = UX1G * CTX + UY1G * CTY
  657. C
  658. C********** Son etat droite
  659. C
  660. PHID = PHIG
  661. RN1D = RN1G
  662. PN1D = PN1G
  663. UN1D = -1.0D0 * UN1G
  664. UT1D = UT1G
  665. DO I1 = 1, NESP
  666. FRAMAS.FRAMD(I1) = FRAMAS.FRAMG(I1)
  667. ENDDO
  668. DO I1 = 1, NESP
  669. ALPHA.FRAMD(I1) = ALPHA.FRAMG(I1)
  670. ENDDO
  671. C
  672. C********** Fin cas mur
  673. C
  674. ELSE
  675. VALCEL = MPPHI.VPOCHA(NLCED, 1)
  676. LIMCEL = MPLPHI.VPOCHA(NLCED, 1)
  677. DCEL = (MPGPHI.VPOCHA(NLCED, 1) * DXD) +
  678. & (MPGPHI.VPOCHA(NLCED, 2) * DYD)
  679. PHID = VALCEL + (LIMCEL * DCEL)
  680. C
  681. C RN
  682. C
  683. VALCEL = MPRN1.VPOCHA(NLCED, 1)
  684. LIMCEL = MPLRN1.VPOCHA(NLCED, 1)
  685. DCEL = (MPGRN1.VPOCHA(NLCED, 1) * DXD) +
  686. & (MPGRN1.VPOCHA(NLCED, 2) * DYD)
  687. RN1D = VALCEL + (LIMCEL * DCEL)
  688. C
  689. C PN
  690. C
  691. VALCEL = MPPN1.VPOCHA(NLCED, 1)
  692. LIMCEL = MPLPN1.VPOCHA(NLCED, 1)
  693. DCEL = (MPGPN1.VPOCHA(NLCED, 1) * DXD) +
  694. & (MPGPN1.VPOCHA(NLCED, 2) * DYD)
  695. PN1D = VALCEL + (LIMCEL * DCEL)
  696. C
  697. C VN
  698. C
  699. VALCEL = MPVN1.VPOCHA(NLCED, 1)
  700. LIMCEL = MPLVN1.VPOCHA(NLCED, 1)
  701. DCEL = MPGVN1.VPOCHA(NLCED, 1)*DXD +
  702. & MPGVN1.VPOCHA(NLCED, 2)*DYD
  703. UX1D = VALCEL + LIMCEL * DCEL
  704. C
  705. VALCEL = MPVN1.VPOCHA(NLCED, 2)
  706. LIMCEL = MPLVN1.VPOCHA(NLCED, 2)
  707. DCEL = MPGVN1.VPOCHA(NLCED, 3)*DXD +
  708. & MPGVN1.VPOCHA(NLCED, 4)*DYD
  709. UY1D = VALCEL + LIMCEL * DCEL
  710. C
  711. C YN
  712. C
  713. DO I1 = 1, NESP
  714. INDCEL = 2 * I1 - 1
  715. VALCEL = MPYMA.VPOCHA(NLCED,I1)
  716. DCEL = MPGYMA.VPOCHA(NLCED, INDCEL)*DXD +
  717. & MPGYMA.VPOCHA(NLCED,INDCEL + 1 )*DYD
  718. LIMCEL = MPLYMA.VPOCHA(NLCED,I1)
  719. FRAMAS.FRAMD(I1) = VALCEL + (LIMCEL * DCEL)
  720. ENDDO
  721. C
  722. C ALPHAN
  723. C
  724. DO I1 = 1, NESP
  725. INDCEL = 2 * I1 - 1
  726. VALCEL = MPALC.VPOCHA(NLCED,I1)
  727. DCEL = MPGALC.VPOCHA(NLCED, INDCEL)*DXD +
  728. & MPGALC.VPOCHA(NLCED,INDCEL + 1 )*DYD
  729. LIMCEL = MPLALC.VPOCHA(NLCED,I1)
  730. ALPHA.FRAMD(I1) = VALCEL + (LIMCEL * DCEL)
  731. ENDDO
  732.  
  733. C
  734. C
  735. C********** Si l'on fait pas de prediction, ce n'est pas necessaire de
  736. C controller la positivite' de la pression et de la densité; elle
  737. C est déjà garantie par la proprieté LED de limiteur.
  738. C If we want to check it, just uncomment LOGI1.
  739. C
  740. C LOGI1 = LOGI1 .OR. (RN1D .LT. 0.0D0)
  741. C $ .OR.(PN1D .LT. 0.0D0)
  742. C
  743. IF(LOGI1)THEN
  744. C
  745. C************* Premier ordre en espace local
  746. C
  747. PHIG = MPPHI.VPOCHA(NLCEG,1)
  748. RN1G = MPRN1.VPOCHA(NLCEG,1)
  749. PN1G = MPPN1.VPOCHA(NLCEG,1)
  750. UX1G = MPVN1.VPOCHA(NLCEG,1)
  751. UY1G = MPVN1.VPOCHA(NLCEG,2)
  752. DO I1 = 1, NESP
  753. FRAMAS.FRAMG(I1) = MPYMA.VPOCHA(NLCEG,I1)
  754. ENDDO
  755. DO I1 = 1, NESP
  756. ALPHA.FRAMG(I1) = MPALC.VPOCHA(NLCEG,I1)
  757. ENDDO
  758. C
  759. PHID = MPPHI.VPOCHA(NLCED,1)
  760. RN1D = MPRN1.VPOCHA(NLCED,1)
  761. PN1D = MPPN1.VPOCHA(NLCED,1)
  762. UX1D = MPVN1.VPOCHA(NLCED,1)
  763. UY1D = MPVN1.VPOCHA(NLCED,2)
  764. DO I1 = 1, NESP
  765. FRAMAS.FRAMD(I1) = MPYMA.VPOCHA(NLCED,I1)
  766. ENDDO
  767. DO I1 = 1, NESP
  768. ALPHA.FRAMD(I1) = MPALC.VPOCHA(NLCED,I1)
  769. ENDDO
  770. ENDIF
  771. C
  772. UN1G = UX1G * CNX + UY1G * CNY
  773. UT1G = UX1G * CTX + UY1G * CTY
  774. C
  775. UN1D = UX1D * CNX + UY1D * CNY
  776. UT1D = UX1D * CTX + UY1D * CTY
  777. C
  778. ENDIF
  779. C
  780. C
  781. C******** Les MELVALs
  782. C
  783. MELPHI.VELCHE(1,NLCF) = PHIG
  784. MELPHI.VELCHE(3,NLCF) = PHID
  785. C
  786. MELRN1.VELCHE(1,NLCF) = RN1G
  787. MELRN1.VELCHE(3,NLCF) = RN1D
  788. C
  789. MELPN1.VELCHE(1,NLCF) = PN1G
  790. MELPN1.VELCHE(3,NLCF) = PN1D
  791. C
  792. MEVUN1.VELCHE(1,NLCF) = UN1G
  793. MEVUN1.VELCHE(3,NLCF) = UN1D
  794. MEVUT1.VELCHE(1,NLCF) = UT1G
  795. MEVUT1.VELCHE(3,NLCF) = UT1D
  796. MELVNX.VELCHE(1,NLCF) = CNX
  797. MELVNY.VELCHE(1,NLCF) = CNY
  798. MELT1X.VELCHE(1,NLCF) = CTX
  799. MELT1Y.VELCHE(1,NLCF) = CTY
  800. C
  801. DO I1 = 1, NESP
  802. MELVA1 = MCHAMY.IELVAL(I1)
  803. MELVA1.VELCHE(1,NLCF) = FRAMAS.FRAMG(I1)
  804. MELVA1.VELCHE(3,NLCF) = FRAMAS.FRAMD(I1)
  805. ENDDO
  806. C
  807. DO I1 = 1, NESP
  808. MELVA1 = MCHAMA.IELVAL(I1)
  809. MELVA1.VELCHE(1,NLCF) = ALPHA.FRAMG(I1)
  810. MELVA1.VELCHE(3,NLCF) = ALPHA.FRAMD(I1)
  811. ENDDO
  812. C
  813. ENDDO
  814. C
  815. C**** Desactivation des SEGMENTs
  816. C
  817. SEGDES IPT1
  818. SEGDES IPT2
  819. C
  820. C MPOVALs
  821. C
  822. SEGDES MPNORM
  823. C
  824. SEGDES MPPHI
  825. SEGDES MPGPHI
  826. SEGDES MPLPHI
  827. C
  828. SEGDES MPRN1
  829. SEGDES MPGRN1
  830. SEGDES MPLRN1
  831. C
  832. SEGDES MPPN1
  833. SEGDES MPGPN1
  834. SEGDES MPLPN1
  835. C
  836. SEGDES MPVN1
  837. SEGDES MPGVN1
  838. SEGDES MPLVN1
  839. C
  840. C MELVALs
  841. C
  842. SEGDES MELVNX
  843. SEGDES MELVNY
  844. SEGDES MELT1X
  845. SEGDES MELT1Y
  846. SEGDES MEVUN1
  847. SEGDES MEVUT1
  848. C
  849. SEGDES MELPHI
  850. C
  851. SEGDES MELRN1
  852. C
  853. SEGDES MELPN1
  854. C
  855. IF (NESP .GE. 1)THEN
  856. SEGDES MPYMA
  857. SEGDES MPGYMA
  858. SEGDES MPLYMA
  859. SEGDES MPALC
  860. SEGDES MPGALC
  861. SEGDES MPLALC
  862. DO I1 = 1, NESP
  863. MELVA1 = MCHAMY.IELVAL(I1)
  864. SEGDES MELVA1
  865. MELVA1 = MCHAMA.IELVAL(I1)
  866. SEGDES MELVA1
  867. ENDDO
  868. SEGDES MCHAMA
  869. SEGDES MCHAMY
  870. SEGSUP FRAMAS
  871. SEGSUP ALPHA
  872. ENDIF
  873. C
  874. C**** Destruction du MELNTI correspondance local/global
  875. C
  876. SEGSUP MLENT1
  877. C
  878. 9999 CONTINUE
  879. C
  880. RETURN
  881. END
  882.  
  883.  
  884.  
  885.  
  886.  
  887.  
  888.  
  889.  
  890.  
  891.  
  892.  
  893.  
  894.  
  895.  

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