Télécharger la plaquette Cast3M


fr en

Télécharger kfmm1.eso

Retour à la liste

Numérotation des lignes : 

kfmm1

  1. C KFMM1     SOURCE    OF166741  24/12/13    21:16:09     12097          
  2.       SUBROUTINE KFMM1(IRN,IGN,IRETN,IGAMN,
  3.      &     ICHPSU,ICHPDI,ICHPVO,INORM,
  4.      &     MELEMC,MELEMF,MELEFE,DCFL,
  5.      &     MELLIM,ICHLIM,
  6.      &     ICHRES)
  7. C************************************************************************
  8. C
  9. C PROJET            :  CASTEM 2000
  10. C
  11. C NOM               :  KFMM1
  12. C
  13. C DESCRIPTION       :  Voir KFMM
  14. C
  15. C                      Cas deux dimensions, gaz "calorically perfect"
  16. C
  17. C LANGAGE           :  FORTRAN 77 + ESOPE 2000 (avec estensions CISI)
  18. C
  19. C AUTEUR            :  A. BECCANTINI, DRN/DMT/SEMT/TTMF
  20. C
  21. C************************************************************************
  22. C ENTREES
  23. C
  24. C
  25. C  1) Pointeurs des CHPOINTs
  26. C
  27. C     IRN     : CHPOINT 'CENTRE' contenant la masse volumique
  28. C
  29. C     IGN     : CHPOINT 'CENTRE' contenant la q.d.m.
  30. C
  31. C     IRETN   : CHPOINT 'CENTRE' contenant l'energie totale
  32. C
  33. C     IGAMN   : CHPOINT 'CENTRE' contenant le gamma
  34. C
  35. C     ICHLIM  : CHPOINT conditions aux bords
  36. C
  37. C  3) Pointeurs de CHPOINTs de la table DOMAINE
  38. C
  39. C     ICHPSU  : CHPOINT "FACE" contenant la surface des faces
  40. C
  41. C     ICHPDI  : CHPOINT "CENTRE" contenant le diametre minimum
  42. C               de chaque element
  43. C
  44. C     ICHPVO  : CHPOINT "CENTRE" contenant le volume
  45. C               de chaque element
  46. C
  47. C     INORM   : CHPOINT "CENTRE" contenant le normales aux faces
  48. C
  49. C
  50. C  4) Pointeurs de MELEME de la table DOMAINE
  51. C
  52. C     MELEMC  : MELEME 'CENTRE' du SPG des CENTRES
  53. C
  54. C     MELEMF  : MELEME 'FACE' du SPG des FACES
  55. C
  56. C     MELEFE  : MELEME 'FACEL' du connectivité FACES -> ELEM
  57. C
  58. C     MELLIM  : MELEME SPG conditions aux bords
  59. C
  60. C  5)
  61. C
  62. C     DCFL = le double de la CFL
  63. C
  64. C
  65. C SORTIES
  66. C
  67. C     ICHRES  : resultat
  68. C
  69. C************************************************************************
  70. C
  71. C HISTORIQUE (Anomalies et modifications éventuelles)
  72. C
  73. C HISTORIQUE : Avril 2002  : creation
  74. C              Janvier 2003: implementation de condition aux limites
  75. C
  76. C************************************************************************
  77. C
  78. C
  79. C N.B.: On suppose qu'on a déjà controllé RO, P > 0
  80. C                                         GAMMA \in (1,3)
  81. C                                         Y \in (0,1)
  82. C       Si non il faut le faire!!!
  83. C
  84. C************************************************************************
  85. C
  86.       IMPLICIT INTEGER(I-N)
  87.       INTEGER IRN,IGN,IRETN,IGAMN,ICHPSU,ICHPDI,ICHPVO,INORM
  88.      &     ,MELEMC,MELEMF,MELEFE,ICHRES,ISPG1,ISPG2,NFAC
  89.      &     ,IGEOMF,IGEOMC, NGCEG, NGCED, NGCF, NLCF, NLCF1, NLCEG, NLCED
  90.      &     ,ICEN,NCEN,ICHLIM,MELLIM,NLLIM
  91.       REAL*8 ROG, RUXG, RUYG, UNG, RETG, GAMG, REC, PG, VOLG, DIAMG
  92.      &     , ROD, RUXD, RUYD, UND, RETD, GAMD, PD, VOLD, DIAMD
  93.      &     , CNX, CNY, DCFL
  94.      &     , SURF, SIGMAF
  95.      &     , UNSDT, UNSDTF, UXD, UYD
  96.       CHARACTER*8 TYPE
  97. C
  98. C**** LES INCLUDES
  99. C
  100.  
  101. -INC PPARAM
  102. -INC CCOPTIO
  103. -INC SMCHPOI
  104.       POINTEUR MPOVSU.MPOVAL, MPOVDI.MPOVAL
  105.      &       , MPORES.MPOVAL, MPODTI.MPOVAL, MPOVOL.MPOVAL
  106.      &       , MPRN.MPOVAL, MPGN.MPOVAL, MPRETN.MPOVAL, MPGAMN.MPOVAL
  107.      &       , MPNORM.MPOVAL, MPLIM.MPOVAL
  108. -INC SMELEME
  109. -INC SMLMOTS
  110. -INC SMLENTI
  111.       POINTEUR MLELIM.MLENTI
  112. C
  113. C**** Initialisation des MLENTI des conditions aux limites
  114. C
  115. C
  116.       CALL KRIPAD(MELLIM,MLELIM)
  117. C     SEGINI MLELIM
  118. C
  119. C**** Initialisation des MELEMEs
  120. C
  121. C     'CENTRE', 'FACEL'
  122. C
  123.       IPT2 = MELEFE
  124.       SEGACT IPT2
  125.       NFAC = IPT2.NUM(/2)
  126. C
  127. C**** KRIPAD pour la correspondance global/local de centre
  128. C
  129.       CALL KRIPAD(MELEMC,MLENT1)
  130. C
  131. C**** MLENTI1 a MCORD.XCOORD(/1)/(IDIM+1) elements
  132. C
  133. C     Si i est le numero global d'un noeud de ICEN,
  134. C     MLENT1.LECT(i) contient sa position, i.e.
  135. C
  136. C     I              = numero global du noeud centre
  137. C     MLENT1.LECT(i) = numero local du noeud centre
  138. C
  139. C     MLENT1 déjà activé, i.e.
  140. C
  141. C     SEGINI MLENT1
  142. C
  143. C
  144. C**** KRIPAD pour la correspondance global/local de 'FACE'
  145. C
  146.       CALL KRIPAD(MELEMF,MLENT2)
  147. C     SEGINI MLENT2
  148. C
  149. C
  150. C**** CHPOINTs de la table DOMAINE
  151. C
  152.       CALL LICHT(ICHPSU,MPOVSU,TYPE,IGEOMF)
  153.       CALL LICHT(ICHPDI,MPOVDI,TYPE,IGEOMC)
  154.       CALL LICHT(ICHPVO,MPOVOL,TYPE,IGEOMC)
  155.       CALL LICHT(INORM,MPNORM,TYPE,IGEOMC)
  156. C
  157. C**** LICHT active les MPOVALs en *MOD
  158. C
  159. C     i.e.
  160. C
  161. C     SEGACT MPOVSU*MOD
  162. C     SEGACT MPOVOL*MOD
  163. C     SEGACT MPOVDI*MOD
  164. C     SEGACT MPNORM*MOD
  165. C
  166.       CALL LICHT(ICHRES,MPORES,TYPE,IGEOMC)
  167. C     SEGACT MPORES*MOD
  168. C
  169. C     MPODTI initialisé a zero; MPODTI = 1 / DT
  170. C
  171.       SEGINI, MPODTI=MPORES
  172. C
  173.       CALL LICHT(IRN,MPRN,TYPE,IGEOMC)
  174.       CALL LICHT(IGN,MPGN,TYPE,IGEOMC)
  175.       CALL LICHT(IRETN,MPRETN,TYPE,IGEOMC)
  176.       CALL LICHT(IGAMN,MPGAMN,TYPE,IGEOMC)
  177.       CALL LICHT(ICHLIM,MPLIM,TYPE,MELLIM)
  178. C
  179. C     SEGACT MPRN*MOD
  180. C     SEGACT MPGN*MOD
  181. C     SEGACT MPRETN*MOD
  182. C     SEGACT MPGAMN*MOD
  183. C     SEGACT MPLIM*MOD
  184. C
  185.       IF(IGEOMF .NE. MELEMF)THEN
  186.          WRITE(IOIMP,*) 'Anomalie dedans kfmm1.eso'
  187.          WRITE(IOIMP,*) 'Probleme de SPG'
  188. C        21 2
  189. C        Données incompatibles
  190.          CALL ERREUR(21)
  191.          GOTO 9999
  192.       ENDIF
  193.       IF(IGEOMC .NE. MELEMC)THEN
  194.          WRITE(IOIMP,*) 'Anomalie dedans kfmm1.eso'
  195.          WRITE(IOIMP,*) 'Probleme de SPG'
  196. C        21 2
  197. C        Données incompatibles
  198.          CALL ERREUR(21)
  199.          GOTO 9999
  200.       ENDIF
  201. C
  202. C
  203. C**** BOUCLE SUR FACEL pour le calcul du FLUX
  204. C
  205.       DO NLCF = 1, NFAC
  206. C
  207. C******* NLCF  = numero local du centre de facel
  208. C        NGCF  = numero global du centre de facel
  209. C        NLCF1 = numero local du centre de face
  210. C        NGCEG = numero global du centre ELT "gauche"
  211. C        NLCEG = numero local du centre ELT "gauche"
  212. C        NGCED = numero global du centre ELT "droite"
  213. C        NLCED = numero local du centre ELT "droite"
  214. C
  215.          NGCEG = IPT2.NUM(1,NLCF)
  216.          NGCED = IPT2.NUM(3,NLCF)
  217.          NGCF  = IPT2.NUM(2,NLCF)
  218.          NLCF1 = MLENT2.LECT(NGCF)
  219.          NLCEG = MLENT1.LECT(NGCEG)
  220.          NLCED = MLENT1.LECT(NGCED)
  221. C
  222. C******* NLCF != NLCF1 -> l'auteur (MOI) n'a rien compris.
  223. C
  224.          IF(NLCF .NE. NLCF1)THEN
  225.             WRITE(IOIMP,*) 'Anomalie dedans kfmm1.eso'
  226.             WRITE(IOIMP,*) 'Probleme de SPG'
  227. C           21 2
  228. C           Données incompatibles
  229.             CALL ERREUR(21)
  230.             GOTO 9999
  231.          ENDIF
  232. C
  233.          CNX = MPNORM.VPOCHA(NLCF,1)
  234.          CNY = MPNORM.VPOCHA(NLCF,2)
  235. C
  236. C******* Recuperation des Etats "gauche" et "droite"
  237. C
  238. C
  239.          ROG = MPRN.VPOCHA(NLCEG,1)
  240.          RUXG = MPGN.VPOCHA(NLCEG,1)
  241.          RUYG = MPGN.VPOCHA(NLCEG,2)
  242.          UNG = (RUXG * CNX) + (RUYG * CNY)
  243.          UNG = UNG / ROG
  244.          RETG = MPRETN.VPOCHA(NLCEG,1)
  245.          GAMG = MPGAMN.VPOCHA(NLCEG,1)
  246.          REC = 0.5D0 * ((RUXG * RUXG) + (RUYG*RUYG))
  247.          REC = REC / ROG
  248.          PG = (GAMG - 1.0D0) * (RETG - REC)
  249.          VOLG = MPOVOL.VPOCHA(NLCEG,1)
  250.          DIAMG = MPOVDI.VPOCHA(NLCEG,1)
  251. C
  252.          ROD = MPRN.VPOCHA(NLCED,1)
  253.          RUXD = MPGN.VPOCHA(NLCED,1)
  254.          RUYD = MPGN.VPOCHA(NLCED,2)
  255.          RETD = MPRETN.VPOCHA(NLCED,1)
  256.          GAMD = MPGAMN.VPOCHA(NLCED,1)
  257.          REC = 0.5D0 * ((RUXD * RUXD) + (RUYD*RUYD))
  258.          REC = REC / ROD
  259.          PD = (GAMD - 1.0D0) * (RETD - REC)
  260.          VOLD = MPOVOL.VPOCHA(NLCED,1)
  261.          DIAMD = MPOVDI.VPOCHA(NLCED,1)
  262.          IF(NLCEG .NE. NLCED)THEN
  263.             UND = (RUXD * CNX) + (RUYD * CNY)
  264.             UND = UND / ROD
  265.          ELSE
  266. C           Murs au condition aux limite
  267.             NLLIM=MLELIM.LECT(NGCF)
  268.             IF(NLLIM .EQ.0)THEN
  269. C              Mur
  270.                UND = -1.0D0 * UNG
  271.             ELSE
  272.                ROD = MPLIM.VPOCHA(NLLIM,1)
  273.                UXD = MPLIM.VPOCHA(NLLIM,2)
  274.                UYD = MPLIM.VPOCHA(NLLIM,3)
  275.                PD = MPLIM.VPOCHA(NLLIM,4)
  276.                GAMD = GAMG
  277.                UND = (UXD * CNX) + (UYD * CNY)
  278.                RETD = ((1.0D0/(GAMD - 1.0D0))*PD)+
  279.      &              (0.5D0*ROD*((UXD*UXD)+(UYD*UYD)))
  280.                VOLD = VOLG
  281.             ENDIF
  282.          ENDIF
  283. C
  284.          SURF   = MPOVSU.VPOCHA(NLCF,1)
  285.          SIGMAF = (0.5D0 * (GAMG + GAMD)) * (PG + PD) / (ROG + ROD)
  286.          SIGMAF = SIGMAF ** 0.5D0
  287.          SIGMAF = (0.5D0 * (ABS(UNG) + ABS(UND))) + SIGMAF
  288. C
  289. C******* On a defini (ROg,ROUNg,ROUTg,Pg,(Yg)), (ROd,ROUNd,ROUTd,Pd,(Yd))
  290. C        et on a déjà verifié ROg, ROd, Pg, Pd > 0 et 0<Y_i<1
  291. C
  292.          MPORES.VPOCHA(NLCEG,1) = MPORES.VPOCHA(NLCEG,1) +
  293.      &        (SURF * SIGMAF / (2.0D0 * VOLG))
  294.          IF(NLCED .NE. NLCEG)
  295.      &        MPORES.VPOCHA(NLCED,1) = MPORES.VPOCHA(NLCED,1) +
  296.      &        (SURF * SIGMAF / (2.0D0 * VOLD))
  297. C
  298.          UNSDT=MPODTI.VPOCHA(NLCEG,1)
  299.          UNSDTF = SIGMAF / DIAMG
  300.          IF(UNSDT .LT. UNSDTF) MPODTI.VPOCHA(NLCEG,1)=UNSDTF
  301.          UNSDT=MPODTI.VPOCHA(NLCED,1)
  302.          UNSDTF = SIGMAF / DIAMD
  303.          IF(UNSDT .LT. UNSDTF) MPODTI.VPOCHA(NLCED,1)=UNSDTF
  304.       ENDDO
  305. C
  306.       NCEN=MPODTI.VPOCHA(/1)
  307. C
  308.       DO ICEN=1,NCEN,1
  309.          MPORES.VPOCHA(ICEN,1)=MPORES.VPOCHA(ICEN,1)+
  310.      &        (MPODTI.VPOCHA(ICEN,1) / (0.5D0 * DCFL))
  311.       ENDDO
  312. C
  313.       SEGSUP MPODTI
  314.       SEGDES MPOVSU
  315.       SEGDES MPOVDI
  316.       SEGDES MPOVOL
  317.       SEGDES MPNORM
  318. C
  319.       SEGDES MPORES
  320. C
  321.       SEGDES MPRN
  322.       SEGDES MPRETN
  323.       SEGDES MPGN
  324.       SEGDES MPGAMN
  325. C
  326.       SEGDES IPT2
  327.       SEGSUP MLENT1
  328.       SEGSUP MLENT2
  329. C
  330.       SEGSUP MLELIM
  331.       IF(MPLIM .GT. 0) SEGDES MPLIM
  332. C
  333.  9999 CONTINUE
  334. C
  335.       RETURN
  336.       END
  337. C
  338.  
  339.  
  340.  
  341.  
  342.  
  343.  
  344.  
  345.  
  346.  
  347.  
  348.  
  349.  
  350.  
  351.  
  352.  

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